# Nested coevolutionary networks shape the ecological relationships of ticks, hosts, and the Lyme disease bacteria of the Borrelia burgdorferi (s.l.) complex

Agustin Estrada-Pena, Hein Sprong, Alejandro Cabezas Cruz, Jose de La Fuente, Ana Ramo, Elena Claudia Coipan

<figure><img src="image_1.png" /></figure>

## To cite this version:

Agustin Estrada-Pena, Hein Sprong, Alejandro Cabezas Cruz, Jose de La Fuente, Ana Ramo, et al. Nested coevolutionary networks shape the ecological relationships of ticks, hosts, and the Lyme disease bacteria of the Borrelia burgdorferi (s.l.) complex. Parasites & Vectors, 2016, 9 (1), 10.1186/s13071-016-1803-z. doi: 10.1186/s13071-016-1803-z

HAL Id: hal-02634961

https://hal.inrae.fr/hal-02634961

Submitted on 27 May 2020

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_2.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
RESEARCH

Open Access

<figure><img src="image_3.png" /></figure>

CrossMark

# Nested coevolutionary networks shape the ecological relationships of ticks, hosts, and the Lyme disease bacteria of the *Borrelia burgdorferi* (s.l.) complex

Agustín Estrada-Peña<fnref n="1" />*, Hein Sprong², Alejandro Cabezas-Cruz³, José de la Fuente⁴,⁵, Ana Ramo<fnref n="1" /> and Elena Claudia Coipan²

## Abstract

**Background:** The bacteria of the *Borrelia burgdorferi* (s.l.) (BBG) complex constitute a group of tick-transmitted pathogens that are linked to many vertebrate and tick species. The ecological relationships between the pathogens, the ticks and the vertebrate carriers have not been analysed. The aim of this study was to quantitatively analyse these interactions by creating a network based on a large dataset of associations. Specifically, we examined the relative positions of partners in the network, the phylogenetic diversity of the tick's hosts and its impact on BBG circulation. The secondary aim was to evaluate the segregation of BBG strains in different vectors and reservoirs.

**Results:** BBG circulates through a nested recursive network of ticks and vertebrates that delineate closed clusters. Each cluster contains generalist ticks with high values of centrality as well as specialist ticks that originate nested sub-networks and that link secondary vertebrates to the cluster. These results highlighted the importance of host phylogenetic diversity for ticks in the circulation of BBG, as this diversity was correlated with high centrality values for the ticks. The ticks and BBG species in each cluster were not significantly associated with specific branches of the phylogeny of host genera (R² = 0.156, P = 0.784 for BBG; R² = 0.299, P = 0.699 for ticks). A few host genera had higher centrality values and thus higher importance for BBG circulation. However, the combined contribution of hosts with low centrality values could maintain active BBG foci. The results suggested that ticks do not share strains of BBG, which were highly segregated among sympatric species of ticks.

**Conclusions:** We conclude that BBG circulation is supported by a highly redundant network. This network includes ticks with high centrality values and high host phylogenetic diversity as well as ticks with low centrality values. This promotes ecological sub-networks and reflects the high resilience of BBG circulation. The functional redundancy in BBG circulation reduces disturbances due to the removal of vertebrates as it allows ticks to fill other biotic niches.

**Keywords:** *Borrelia burgdorferi* (s.l.), Ticks, Vertebrates, Ecological networks, Centrality, Phylogenetic diversity

<footnote>* Correspondence: aestrada@unizar.es</footnote>

<footnote n="1">¹Department of Animal Health, Faculty of Veterinary Medicine, University of Zaragoza, Miguel Servet, 177, Zaragoza, Spain</footnote>

<footnote>Full list of author information is available at the end of the article</footnote>
## Background

The bacteria of the complex *Borrelia burgdorferi* (*sensu lato*) (BBG) are a group of related tick-transmitted pathogens. They are distributed all over the Holarctic biogeographical region [1], with a recent record from South America [2]. A vast range of vertebrates (lizards, birds and mammals) acts as reservoirs for these pathogens. Vertebrates also act as hosts of the ticks that can circulate the pathogens among animals and transmit them to humans [3]. This group of pathogens has been revised considerably in the last decade after the discovery of new strains by phylogenetic, ecological and clinical studies revealed important distinctions among the strains [1, 4, 5]. There is major interest in understanding the ecological processes that lead to infection by a particular species, in the factors that influence the interactions of the reservoirs, their vectors and the prevailing climate [6–8] and in the molecular features that delineate the geographical distribution of the pathogens and their spread [9–11]. There have also been advances in determining the phylogeny of the group using molecular markers [12] and multilocus sequence typing of housekeeping genes [9]. These advances have furthered our knowledge of the complex processes underlying the distribution and the speciation of the pathogens. However, there has not been a major analysis of the ecological relationships of the species of *Borrelia*, the ticks in which they were recorded and the vertebrates that host the bacteria or feed the ticks.

Interactions among species shape biodiversity as communities of organisms that interact and co-evolve across time and space [13–16]. Recent studies have stressed the need to capture the co-evolutionary mechanisms of large sets of interacting species, since a local view or a focus on small groups of species does not encompass their context [17, 18]. The investigation of interacting organisms has produced a corpus, termed the network analysis, that describes the relationships and that stresses the critical need to understand both the persistence and the co-evolution of large assemblages of species in which diversity is the rule [17, 19]. Network theory can describe ecological communities and the distribution of species among the communities and can provide tools to quantify their interactions [20]. Notably, the complexity of a set of interacting species in a network approaches the complexity of the ecological background when there are multiple species in a community [21]. In this context, the interactions between parasites and vertebrates can be considered to be co-evolved structures in networks, a view that supersedes the classic interpretation of diffuse multispecific interactions [22]. The biological interactions of vertebrates, parasites and transmitted pathogens can thus be depicted as a web in which the nodes represent the organisms and the edges or links represent their interactions [23]. In this approach, the properties of the system emerge

from the properties of the organisms, allowing the system to be investigated at both the organismal and the network levels [24]. Accordingly, the extent and the co-occurrence properties of pathogens, ticks and vertebrates can reveal ecological properties that impact the persistence of these complex networks. Generally speaking, networks are characterised by distinct network subgroups (clusters or modules) that are composed of organisms that interact preferentially among themselves rather than with other nodes in the network. In host-parasite networks, closely related host species comprise clusters, in part because they have similar characteristics that determine their compatibility with parasites [25].

The nature of the associations is largely unexplored in tick-transmitted pathogen systems [23, 26, 27]. In this study, we hypothesised that the structure of the network of BBG and its vertebrate hosts and vectors can be investigated using methods derived from network theory. This approach would show the ecological associations and relationships that influence the circulation of the pathogen in active foci. We examined this hypothesis using modularity, which is a network property that is crucial to the ecology of hosts and parasites [23, 25, 28], as a way to investigate the organisation of a large set of records of vertebrates, ticks and BBG pathogens. We further hypothesised that centrality indices can be used to quantify the roles of the partners in this epidemiological network. We demonstrated that the structure of the networks in which BBG circulates is related to the relationships between ticks and hosts and that these relationships can be described based on their centrality indices.

## Methods

### Literature search

Data on the pair-wise systematic associations of ticks in vertebrates, BBG in vertebrates and BBG in ticks were compiled from a literature review of journals that are searchable in the Thomson Reuters, Scopus and PubMed databases. A “record” was a pair-wise BBG-tick-vertebrate combination at a single geographical site. If the same pair-wise combination was observed at the same site several times (for example, in seasonal collections), it was included only once. However, the same combination of partners at different sites was included once for each site in order to estimate how common that association was, and then the number of reports was used for weighting (see below). The literature search extended from January 1990 to December 2014. It was evident in the initial stages of the bibliographical search that the use of a long series of keywords (i.e. a series based on the names of tick species and vertebrate hosts) would miss relevant studies. We therefore decided to perform a deliberately less stringent query that only included the generic names “Ixodes” OR “Borrelia” and then to critically evaluate the abstracts in
order to manually select the papers that dealt only with ecological information. Papers addressing the clinical presentation, treatment or diagnosis of infections in humans were only included in the dataset if adequate information about the association of a species of tick and BBG was reported. We removed records involving livestock because they are considered accidental hosts, and those records would generate spurious information that would distort the relationships of the community of tick-borne pathogens [23]. There are no compilations of records of BBG detected on ticks while feeding because the bacterial DNA cannot be reliably attributed to either the vertebrate or the vector. It is important to stress that the network structure requires bivalent interactions (i.e. the pathogen must be associated with either the tick or the vertebrate). These associations ignore the effects of the vertebrate, which may transmit the pathogen to co-feeding ticks [29, 30].

We also removed records involving species in the recurrent fever group, with some exceptions that were used as controls for the computations (see below), as well as records involving the ticks and vertebrates associated with this group. We updated the scientific names of BBG species to include the most recent and accepted ones [1-8]. Because of the lack of availability of the original publication at the time of the compilation, data on *Borrelia chilensis* and *Ixodes stilesi* [2] were not included in the analysis; we included only the records for *I. pararicinus* [31]. The few data on *B. ruski* were included, although this species has been shown to be *B. afzelii* with a highly variable IGS locus [11]. The final database contained 10,972 pairwise interactions among species of vertebrates, BBG and ticks. This is available as Additional file 1.

## Building and exploring the network of interactions

To address the ecological aspects of the relationships among BBG-ticks-vertebrates, we developed a network of biotic interactions. The network depicts a community in which each node represents a species and the link between two nodes represents a relationship. In this particular network, nodes represent “carriers” (ticks or vertebrates) that are linked to “cargo” (ticks or pathogens, since vertebrates may carry either). Thus, the network is directed: each link or edge connects one organism to another in one direction. Host-parasite data are sensitive to the sampling effort [32]. Consequently, the computation of further indices is strongly influenced by the sampling and reporting intensity. To ensure that our findings were robust, we increased the weight of the least sampled species and decreased the weight of the most sampled species. Specifically, we regressed the weight of each edge according to the number of citations of the least sampled species (vertebrate, tick, pathogen)

in each edge. This regression was highly significant (linear regression: $R^2$ = 0.56, $F_{(2,2)}$ = 34.99, $P$ < 0.0001). Next, we additively rescaled the residuals to be greater than zero, since edges cannot have negative weights [32]. The residuals reflect the number of links relative to the sampling effort assuming that the measure of the sampling effort should be based on the less studied species. We replaced the original weights of the edges with the rescaled residuals and then estimated the indices of centrality [32].

Several metrics have been used to investigate the relationships of nodes in networks. The *weighted degree* (WD) is a simple measure of the number of edges that leave or arrive at a given node, and it provides an estimation of the number of nodes connected to every single node. However, the WD does not reflect the importance of the nodes in the context of the network. A more informative measure can be obtained by extending the definition of WD to include the node strength [33], and we use this measure here to indicate the *circulation capacity* (CC) of the pathogen. The CC is a quantitative extension of the species degree, which is the number of interactions per species in a network. The circulation capacity of a species of tick or vertebrate is defined as the sum of dependences (measured as a function of the number and the weight of each edge) of the BBG species that rely on these partners. That is, the CC is a measure of the importance of the vertebrate or the tick relative to the pathogen.

Centrality measures in ecological networks indicate the presence of high-ranking nodes in the network that have significantly higher-than-average connectivity and/or links that stretch far beyond their local network neighbourhoods. Identifying the most central nodes in the system is an important step in network characterization. To quantitatively account for the role of more central nodes, *betweenness centrality* (BNC) was defined as the number of shortest paths between pairs of nodes that pass through a given node. Therefore, central nodes are part of the shortest paths within the network, but peripheral nodes are not. To be in a central node, a carrier had to be infected by many cargo species that infect many other carriers in the network. The vertebrates with the greatest centrality are super-spreaders [34]. It is intuitive to consider the alternative definition of centrality by looking at the *PageRank* (PR), an index of centrality that recursively assigns a universal rank to nodes based on the importance of the other nodes to which it is linked. The *weighted clustering coefficient* (WCC) is a measure of the degree to which nodes in a graph tend to cluster together. The WCC measures the local group cohesiveness and is defined for any node as the fraction of connected neighbours. Thus, the WCC expresses the statistical level of cohesiveness by measuring the global density of interconnected nodes in
the network. Network computations were carried out using *igraph* [35] for R [36], and the ForceAtlas2 algorithm was used to display the network [37].

Because of the high complexity of the network and to provide coherence to the analyses, we restricted further computations to the genera of vertebrates. The purpose of reducing from species to genera is to obtain an overall description of the network, since it was impossible to obtain adequate information of the phylogenetic relationships among the species of hosts (see below). We acknowledge this may introduce some bias in the particular role of a host in supporting the circulation of BBG, but do not change the general indexes of the network [23]. The resulting network has 358 nodes with 288 genera of vertebrates, 50 species of ticks, 20 species of BBG [including *B. burgdorferi* (s.s.)] and 1,026 links. The nodes of the networks tended to create tightly knit groups that were characterised by a relatively high density of links, which is greater than the average probability of a tie that is randomly established between two nodes [38]. We calculated the modularity of the network using the Louvain algorithm [39], obtaining also a measure of its nestedness [17]. It has been demonstrated that networks of parasites and hosts tend to be highly nested, and clusters of these networks include both generalist (parasites with many links) and specialist (parasites with few links, shared with the generalists) [23, 26]. In this case we wanted to know if the nested structure already detected for a large set of ticks and pathogens [23] is already retained in the particular case of BBG. Other than records of the BBG group and its associated ticks and vertebrates, we also kept the records of three groups of ticks, vertebrates and pathogens that were not directly related with BBG. One of these groups of records is *Ixodes lividus*, which is a tick parasite on birds of the genera *Riparia* and *Delichon*, a system in which no BBG species have been found. The second group includes *B. persica* and *Ornithodoros tho-lozani*. This pathogen does not belong to BBG. The third group includes *B. miyamotoi*, which does not belong to BBG but shares a subset of reservoirs with the complex of target pathogens. As a control of the reliability of calculations, we checked whether these two subsystems were separated from the main clusters and routes of circulation of BBG.

## Calculating the phylogenetic relationships of vertebrates

We aimed to evaluate the relationships between the centrality indices of the ticks in the network, the ticks' ability to transmit species of BBG and the phylogenetic diversity of the exploited vertebrates. To obtain the phylogenetic relationships of the vertebrates, we queried the Open Tree of Life (http://www.opentreeoflife.org/; accessed 19 Sept 2016). We found information for 240

genera of vertebrates, which is summarised in the phylogenetic tree available as Additional file 2. The tree is a representation of the phylogenetic relationships among the genera of vertebrates without information about the times of divergence.

Phylogenetic diversity was calculated using *Faith's phylogenetic diversity* (PD) index [40] as the total branch length spanned by the tree, including all of the species of vertebrates in a “community”. The term community is used here to mean all of the genera of vertebrates that are included in one of the network clusters and that are exploited by a species of tick or on which BBG has been recorded. We also calculated the *mean pairwise distance* (MPD) and the *mean nearest taxon distance* (MNTD) as described previously [41]. Null models were generated that randomised the tips of the phylogeny to calculate the significance of the phylogenetic association among ticks, BBG and vertebrates [41]. These calculations are intended to link the indices of the network to the vertebrate's composition of each cluster.

## Multilocus sequence typing and ecological relationships

Multilocus sequence typing (MLST) is a powerful molecular biology technique that is based on using the sequences of 8 chromosomal housekeeping genes to obtain information on population structure and the evolutionary relationships of BBG [9]. We used MLST to evaluate the significance of the associations among strains of BBG and species of ticks or genera of vertebrates. For the MLST analysis, we used the sequences in the *Borrelia* MLST databases (http://pubmlst.org/borrelia; accessed 19 Sept 2016). The data were updated with 613 additional sequences that were obtained from a variety of sources in Europe [42]. The sequences of individual housekeeping genes were concatenated, and those that differed by one or more nucleotides were assigned allele numbers. We assembled a total of 1,832 single associations between BBG, 11 species of ticks and 21 species of vertebrates for a total of 704 unique MLST sequences. We removed all of the strains that were isolated from humans, which constituted more than the half of the available data, since our focus was on wild vertebrates and ticks. All sequences were aligned using the ClustalW MEGA 4 algorithm [43], and phylogenetic trees were constructed using the Kimura two-parameter model in MEGA. We focused on the segregation of MLST sequences and the species of tick or vertebrate in which they have been detected. This analysis should be considered a partial analysis because the available data are not as complete as the data for the relationships of ticks, vertebrates and BBG. All MLST data and data on the ticks or vertebrates in which the pathogens have been reported was translated into a network in which each unique MLST sequence is associated with ticks or vertebrates. We
calculated the centrality indices, Faith's PD, MPD and MNTD as described previously [40, 41]. These are all measures of the phylogenetic diversity of the BBG as calculated for each species of tick or vertebrate, and the aim was to look for a higher-than-expected phylogenetic clustering of strains of BBG with specific vertebrate or tick taxa. Null models were generated by randomizing the tips of the phylogeny to calculate the significance of the phylogenetic association.

## Results

### BBG circulates in a redundant and nested network

The network that supports the circulation of BBG includes 20 species of pathogens of the complex, 18 species of tick vectors and 32 other species of ticks in which the pathogens have not yet been reported but

which share vertebrate hosts with those in which the bacteria have been reported. The network also has 82 genera of vertebrate reservoirs of BBG and 206 genera of vertebrates on which the pathogen has not yet been detected but which have been reported as hosts of one or several species of ticks. The network includes 687 species of vertebrates, with 1,026 unique pairs of associations between genera of vertebrates, ticks and BBG. The network is thus highly redundant in terms of vertebrates: several species of ticks feed on the same genus of vertebrate, and the same species of tick may feed on several genera of vertebrates.

A schematic view of the network clusters is shown in Fig. 1. Additional files 3 and 4 show the complete network at the level of genera of vertebrates and at the level of species of vertebrates, respectively. The network has

<figure><img src="image_4.png" /></figure>
10 clusters. Each cluster is a group of species that are more frequently reported to interact with each other than with other species. The clusters are categorical and are numbered consecutively. The nestedness of the complete network is very high, with a score of 5.2 in a possible range of 0–100, where 0 is the maximum nest- edness. The clustering algorithms show the differences in the two groups of species that were included as in- ternal controls, referred to in Fig. 1 as clusters 8 and 10. These are not connected to the other clusters because they contain either non-BBG species of *Borrelia* or ticks that are not connected to the spectrum of vertebrates of the BBG network. The remaining 8 clusters contain a variable number of ticks, vertebrates and pathogens. Clusters 1, 2 and 6 have the highest number of ticks, species of BBG and vertebrates and contain the Palaearctic taxa and *I. scapularis*. The most prominent articulation points of clusters 1 and 2 have at least five genera of verte- brates. Cluster 5 is also notable in the network as it con- tains all of the New World taxa (except for *I. scapularis*) that circulate on 59 genera of vertebrates and 8 species of BBG, including *B. japonica*. Other clusters are only per-ipherally connected to the main heart of the network. Cluster 3 contains species of ticks that parasitise birds (*I. frontalis* and *I. arboricola*) without the species of BBG associated with that ecological cluster. These ticks seem to be secondary and perhaps opportunistic trans- mitters of BBG. Cluster 7 contains the specialist tick *I. baergi*, in which BBG has been detected; cluster 9 con- tains species of carnivores and other mammals, the ar- ticulation points of which are vertebrates that are shared with ticks of other clusters and with which no BBG species are ecologically related. Cluster 4 is of spe- cial interest because it is strongly associated with both

*B. sinica* and *B. yangtzensis*, both of which were re- ported in the tick *I. granulatus* and in several verte- brates. Clustering algorithms indicate that this cluster only contains species that have been collected in the Palaearctic region. Additional files 5, 6 and 7 include all the numeric data related to ticks, BBG and vertebrates.

BBG circulates through ticks with high values of centrality and phylogenetically diverse vertebrates

We wanted to determine whether the ticks in the differ- ent clusters had different properties in the context of the network. Figure 2 shows that the ticks in the different clusters were distributed randomly in terms of the range of BNC and PR values and that the Faith's PD of the hosts of the ticks was correlated with both the BNC and the PR values. Specifically, the ticks with the highest Faith's PD values tended to occupy the top positions in the network in terms of BNC and PR values (PD-BNC: R² = 0.713, F*(2,2) = 114.57, P < 0.0001; PD-PR: R² = 0.818, F*(2,2) = 207.95, P < 0.0001). *Ixodes ricinus* and *I. persulca- tus* had the highest centrality (>90) and Faith's PD (18.13 and 5.81, respectively) values, with other prominent spe- cies in the transmission cycles of BBG (like *I. scapularis* and *I. pacificus*) showing relatively high Faith's PD (3.08 and 1.51, respectively) values but low BNC (17 and 35, respectively) and PR (3 and 8, respectively). Only two species of ticks with lower BNC (< 60) and PR (< 1) values have been reported as being infected with BBG (*I. columnae* and *I. sinensis*). All the other reported vectors of BBG had high BNC and PR values. We calculated the CC of BBG by ticks and plotted it against their Faith's PD and centrality indices. The result was an almost lin- ear relationship (R² = 0.799, F*(1,2) = 27.65, P < 0.0001) (Fig. 3). However, this correlation was reduced if Faith's

<figure><img src="image_5.png" /><figcaption>**Fig. 2** The distribution of the betweenness centrality and PageRank values of the tick species in the epidemiological *Borrelia burgdorferi* (BBG) network. This plot shows the number of clusters in which the ticks are allocated (indicated by the colours) and the phylogenetic diversity of the vertebrates of each tick species as measured using Faith's phylogenetic diversity index (PD; indicated by the circle size). **b** shows an expanded view of the species in the square at the bottom left in **a**</figcaption></figure>
<figure><img src="image_6.png" /></figure>

PD was replaced by the crude number of genera of hosts ($R^2$ = 0.211, $F_{(1,2)}$ = 1.97, $P$ = 0.349).

Statistical analysis of the Faith's PD of vertebrates showed that the clusters were not linked to specific branches of the phylogenetic tree of vertebrates. Rather, each cluster included groups of vertebrates with high Faith's PD. There were no statistically significant differences in the associations of the ticks or BBG in each cluster with specific branches of the phylogeny of the genera of vertebrates ($R^2$ = 0.156, $P$ = 0.784 for BBG; $R^2$ = 0.299, $P$ = 0.699 for ticks; Additional file 8). Not all species of ticks could be analysed either because reliable phylogenetic data for the vertebrates was lacking or because the ticks were restricted to only one genus of vertebrates. The following tick species were tightly linked to specific branches of the phylogeny of vertebrates: *I. angustus*, *I. arboricola*, *I. cookei*, *I. frontalis* and *I. trianguliceps*. Notably, all of the species of ticks that were restricted to a set of a few vertebrates had low importance in terms of the circulation of BBG, but all of them occurred in clusters in which generalist ticks were predominant. The same calculations were performed using the BBG species and

the phylogeny of the vertebrates on which they were reported. The results showed higher phylogenetic diversity of vertebrates for BBG than for ticks, and only *B. carolinensis* and *B. yangtzensis* were linked to a specific branch group of vertebrate genera. Additional file 8 shows the relationships of the organisms in clusters 1, 2, 5 and 6 with the branches of the vertebrate phylogenetic tree.

A few genera of rodents and medium-sized birds had some of the highest centrality values of the network. They also had the greatest CC values for BBG. The genera of hosts with the highest centrality index values were *Turdus*, *Apodemus* and *Microtus* (Fig. 4). However, the relationship between the CC and the position of the vertebrate in the range of centrality values was not completely linear ($R^2$ from multiple regression = 0.633, $F_{(1,2)}$ = 11.98, $P$ < 0.0001) since some genera of hosts, like *Crocidura* and *Spermophilus*, also had relatively high PR and BNC values and made low contributions to the circulation of BBG. A few hosts with PR and BNC values in the low range might also make important contributions to the circulation of BBG (Fig. 5), such as small, ground-feeding birds (*Erithacus*, *Anthus*)
<figure><img src="image_7.png" /><figcaption>**Fig. 4** The distribution of the betweenness centrality and PageRank values of the genera of hosts in the *Borrelia burgdorferi* (BBG) epidemiological network. This plot shows the circulation capacity of each host genus for BBG (indicated by the size and colour of the circle). An expanded view of the rectangle that contains low values of centrality is shown in Fig. 5</figcaption></figure>

and *Sylvia*), a lizard (*Lacerta*) and two large ungulates (*Odocoileus* and *Cervus*). Data from regressions performed using centrality indices and CC suggested that the PR is the variable with the most influence. The ecological implication is that vertebrates with the highest values of centrality for ticks are also those that support the highest circulation of BBG.

## Ticks do not share strains of BBG

We further evaluated the links between strains of BBG (as defined by MLST analysis) and the taxa of ticks and vertebrates. This test could not be performed for several species of carriers because only one strain was available; therefore, phylogenetic diversity could not be calculated. The results must thus be considered a proof of concept that is still preliminary and incomplete due to the bias in recording and typing strains from the most commonly surveyed ticks and vertebrates. Table 1 shows the phylogenetic diversity of the strains of BBG in the different carriers. Faith’s PD was very low for *I. granulatus* (0.07483), *I. pacificus* (1.517) and *I. pavlosvkyi* (1.991), it was medium for *I. scapularis* (3.0) and very high for *I.*

*ricinus* (18.13). Negative MPD values and values of MPD-P near 0.01 indicate strong phylogenetic clustering of the MLST strains in the carrier species. Although the results should be interpreted cautiously due to the small sample size, we found a strong phylogenetic link between the BBG strains and the tested species of *Ixodes* (0.01 for the species with more than one strain reported) but not *I. ricinus*. Figure 6 shows the network of MLST strains that were linked to unique species of ticks or vertebrates without overlapping BBG strains in the ticks. Additional files 9, 10, 11, 12 and 13 report the MLST-derived phylogenetic tree of the BBG strains and show the links of the ticks to specific portions of the BBG tree. Additional file 14 complements Fig. 4 by labelling the BBG species that were linked to the nodes of ticks or vertebrates. *Ixodes ricinus* and *I. persulcatus* (Palaearctic species) shared only two strains, MLST numbers 86 and 244; however, *I. persulcatus* and *I. pavlovskyi* shared several strains, and these strains were similar to those found in *I. turdus* and *I. granulatus*, which are partially sympatric species. *Ixodes pacificus* (western Nearctic) strains are more similar to strains found in eastern Palaearctic
<figure><img src="image_8.png" /><figcaption>**Fig. 5** The distribution of the betweenness centrality and PageRank values of the genera of hosts in the *Borrelia burgdorferi* (BBG) epidemiological network. This plot shows the circulation capacity of each host genus for BBG (indicated by the circle size and colour). This is an expanded view of the rectangle in Fig. 4</figcaption></figure>

ticks. Fewer strains of BBG have been detected in *I. scapularis* (eastern Nearctic), and they are highly repeated. The strains of BBG in *I. scapularis* and *I. pacificus* (Nearctic) were very different in the context of the network.

## Discussion

Here we present the largest BBG epidemiological network that has been analysed to date. This network showed that the pathogens circulate in closed clusters in a highly nested manner and that they are supported by a few species of ticks that exploit a range of phylogenetically diverse hosts. These ticks are restricted to closed clusters that closely resemble geographical components. We found that the vertebrates that support BBG circulation are those that are most commonly used as hosts by generalist species of ticks. Specialist ticks interact with a few species of vertebrates to generate sub-networks that link other vertebrates that would otherwise remain unconnected to the main web structure.

The constructed network of pathogens, their vertebrate hosts and their vectors formed a structure of redundant interconnections. Some species of BBG linked a

few otherwise isolated clusters, suggesting invasion events and further spread; these have been described previously [10]. However, undetermined BBG genospecies were included in the network in order to use data obtained before modern molecular techniques became available, and this could result in bias in terms of the structure of the clusters. A combination of centrality indices was used to quantify the relative importance of each tick in the network, with centrality being positively correlated with the phylogenetic diversity of the species of ticks but not with its host richness. Vertebrates that had the highest centrality with ticks also contributed more to the circulation of the pathogens. Circulation of the pathogens is thus supported by a wide range of vertebrates and by the generalist ticks that they share. The framework revealed the co-occurrences and interactions of specialised and generalised species that are key to maintaining functional diversity [44], promoting community stability [45] and ensuring network persistence [17].

The existence and underlying causes of clusters in ecological networks are still a matter of debate. In parasite networks, the hosts of the most specialised species of
<table><caption>**Table 1** The phylogenetic diversity of *Borrelia burgdorferi* (s.l.) strains recorded in vertebrates or ticks based on their characterization by multilocus sequence typing (MLST). Negative MPD values and values of MPD-P near 0.01 indicate strong phylogenetic clustering of the MLST strains in the tick or vertebrate; positive MPD values and higher MPD-P values indicate weak phylogenetic clustering of the MLST strains in the tick or vertebrate</caption><thead><tr><th>Species</th><th>PD</th><th>n</th><th>MPD</th><th>MPD-P</th></tr></thead><tbody><tr><td>Apodemus agrarius</td><td>na</td><td>1</td><td>na</td><td>na</td></tr><tr><td>Apodemus speciosus</td><td>2.86735</td><td>8</td><td>-2.10119</td><td>0.05</td></tr><tr><td>Apodemus sylvaticus</td><td>na</td><td>1</td><td>na</td><td>na</td></tr><tr><td>Apodemus uralensis</td><td>0.13776</td><td>4</td><td>-3.64226</td><td>0.01</td></tr><tr><td>Crocitara watasei</td><td>0.02041</td><td>3</td><td>-3.08071</td><td>0.01</td></tr><tr><td>Dipodomys californicus</td><td>na</td><td>1</td><td>na</td><td>na</td></tr><tr><td>Eothenomys smithii</td><td>na</td><td>1</td><td>na</td><td>na</td></tr><tr><td>Microtus pennsylvanicus</td><td>0.0102</td><td>2</td><td>-2.40207</td><td>0.01</td></tr><tr><td>Mus caroli</td><td>0.06633</td><td>4</td><td>-4.52009</td><td>0.01</td></tr><tr><td>Myodes glareolus</td><td>1.70408</td><td>2</td><td>0.57406</td><td>0.72</td></tr><tr><td>Myodes rufocanus</td><td>0.11224</td><td>3</td><td>-4.04499</td><td>0.01</td></tr><tr><td>Neotoma floridana</td><td>na</td><td>1</td><td>na</td><td>na</td></tr><tr><td>Neotoma fuscipes</td><td>0.0068</td><td>2</td><td>-2.0568</td><td>0.01</td></tr><tr><td>Niviventer fulvescens</td><td>na</td><td>1</td><td>na</td><td>na</td></tr><tr><td>Peromyscus gossypinus</td><td>0.0068</td><td>2</td><td>-2.19929</td><td>0.025</td></tr><tr><td>Rattus norvegicus</td><td>na</td><td>1</td><td>na</td><td>na</td></tr><tr><td>Rattus rattus</td><td>0.04762</td><td>2</td><td>-2.06397</td><td>0.02</td></tr><tr><td>Suncus murinus</td><td>0.01361</td><td>2</td><td>-2.44465</td><td>0.02</td></tr><tr><td>Thryothorus ludovicianus</td><td>na</td><td>1</td><td>na</td><td>na</td></tr><tr><td>Haemaphysalis longicornis</td><td>na</td><td>1</td><td>na</td><td>na</td></tr><tr><td>Ixodes granulatus</td><td>0.07483</td><td>6</td><td>-5.10505</td><td>0.01</td></tr><tr><td>Ixodes minor</td><td>na</td><td>1</td><td>na</td><td>na</td></tr><tr><td>Ixodes ovatus</td><td>na</td><td>1</td><td>na</td><td>na</td></tr><tr><td>Ixodes pacificus</td><td>1.51701</td><td>20</td><td>-9.68056</td><td>0.01</td></tr><tr><td>Ixodes pavlovskyi</td><td>1.9915</td><td>9</td><td>-5.07941</td><td>0.01</td></tr><tr><td>Ixodes persulcatus</td><td>5.81122</td><td>116</td><td>-6.09055</td><td>0.01</td></tr><tr><td>Ixodes ricinus</td><td>18.13265</td><td>307</td><td>4.61072</td><td>1</td></tr><tr><td>Ixodes scapularis</td><td>3.0085</td><td>104</td><td>-25.75207</td><td>0.01</td></tr><tr><td>Ixodes spinipalpis</td><td>na</td><td>1</td><td>na</td><td>na</td></tr><tr><td>Ixodes stilesi</td><td>na</td><td>1</td><td>na</td><td>na</td></tr><tr><td>Ixodes turdus</td><td>na</td><td>1</td><td>na</td><td>na</td></tr></tbody></table>

Abbreviations: MPD mean pairwise distance, MPD-P significance of MPD index, na not available, PD Faith's phylogenetic diversity, n number of haplotypes recorded

tick tend to be a subset of the hosts of the ticks that are slightly more generalised; in turn, these hosts are a subset of the hosts of ticks that are even more generalised, and so on [46]. We tested the efficiency of the cluster-detecting algorithm by retaining the records of (i) species of ticks that shared vertebrate hosts with other ticks in which BBG has been detected; (ii) ticks and

vertebrates in which the pathogen has never been reported; and (iii) species of *Borrelia* that are not BBG. In every case, the algorithm assigned the partners to the correct clusters. Specifically, the algorithm placed *B. miyamotoi* (which is not a member of BBG) in a subset of vertebrates that are shared by BBG species, as reported previously [47]; it separated *B. persica* from the circulation of the main network; and it found that *I. lividus* was not involved in BBG circulation. These findings supported other results derived from the network.

The findings from the BBG ecological network agreed with those of other studies of the phylogenetic structure of the populations of other species in BBG. For example, it has been reported that *B. bavariensis* comprises two populations, namely a genetically diverse Asian population and a genetically homogeneous European population. The observed differences in genetic heterogeneity cannot be attributed to different geographic scales of sample origin [48]. For example, other BBG species with a low prevalence or a focal distribution in Europe, such as *B. spielmanii* or *B. lusitaniae*, show genetic heterogeneity according to MLST analysis, while *B. bavariensis* does not [48]. The very low diversity and short branches of the phylogeny of the European population suggest a recent founding event that is probably linked to the shift of a population of *B. bavariensis* to the European vector. *I. ricinus* [48]. Our network results support the idea that *B. bavariensis* is not native to Europe, since it circulates within a cluster of 99 genera of vertebrates and 15 species of ticks that are distributed in the eastern Palaearctic. *Ixodes ricinus*, a Western Palaearctic species, is not a “native” tick in that cluster.

This study has some limitations. The network approach ignores the abilities of the vertebrates to act as reservoirs, a feature that has been evaluated elsewhere using both regional [46, 49, 50] and larger sets of empirical data [11, 51–53]. This approach also ignores transmission inhibition in vertebrates that are infected by several strains [54]. Each vertebrate has different reservoir capacities and thus infects ticks at different rates and contributes in varying ways to the circulation of the pathogen. Since our approach was based on the genera of vertebrates, and since there is a scarcity of experimental data about the reservoir capacity of vertebrates, we could not address these issues in the current study. Another limitation of the study is the phylogenetic construction of the genera of vertebrates. Because of the large number of vertebrate species, the inclusion of taxa from both the Palaearctic and Nearctic regions and the unavailability of data for some species, a complete phylogenetic approach was to use a mega-phylogeny [55, 56]. We eliminated a potential source of errors that commonly arise from literature searches by weighting the number of records [23, 32].
<figure><img src="image_9.png" /></figure>

The ecological interpretation of our findings is that groups of allopatric species of ticks delineate clusters of interacting organisms. These clusters are not phylogenetically linked to the genera of vertebrates, indicating that they are not derived purely from the relationships of ticks with specific branches of vertebrate evolution. The results revealed direct relationships between the CC of the pathogen, the phylogenetic diversity of the vertebrates that feed ticks (but not the raw richness of the host taxa) and the ticks' centrality values. The results also showed the relative

importance of some vertebrates for BBG circulation. Notably, the prevalence of BBG in ticks is due to a complex relationship based on the tick densities on available hosts, on the capacity of the hosts to feed the ticks and on the capacity of the hosts to transmit BBG to the feeding ticks [57]. Thus, the tick burdens of vertebrates are positively correlated with the infection prevalence of BBG, suggesting that the majority of ticks are fed by only a few species of hosts. Our results are in line with previous reports linking the life-cycle traits of the hosts with their tick burden
and to the infection prevalence of BBG in hosts [52, 58]. Considering our simplification, which included the use of the genera of hosts rather than the species of hosts, the results from the network still overlapped and agreed with previous quantitative reports on the topic. In particular, they highlighted the importance of medium-sized birds, some rodents and Insectivora in supporting active BBG foci. It is important to stress that the same basic ecological structure is found using different methods, including quantitative reviews of the BBG-vertebrate system [57], the life-cycle traits of the hosts [52] and our network, which was constructed using crude association data on ecological relationships. Vertebrates that were identified as being the most important for the circulation of ticks were those with the highest importance for pathogen circulation, an intuitively obvious finding that now has a mathematically tractable framework. It must to be stressed that the circulation of BBG is supported by either the reservoir ability of the hosts and/or the feeding of large number of ticks. This may be the case of some large ungulates, which are hosts of large numbers of ticks but are recognised as inefficient reservoirs of the pathogen.

We speculate, however, that the importance of under-represented hosts (e.g. those with low centrality) has been neglected. Alone, the individual contributions of hosts in the low range of centrality would probably not be able to support active BBG foci in the absence of the most central vertebrates, but their collective contributions, even with low values of centrality, would support the successful circulation of BBG in biodiverse foci. The network could be revised to integrate previous quantitative approaches into the ecological background in order to further study the influence of vertebrates and ticks in the complex relationships of BBG. Further research on the topic is needed, including the individual contributions of species of hosts (i.e. according to body size and number of feeding ticks) although this seems to be a very local feature, that may change even from year to year according to the relative faunal composition of available hosts and the tick-host relationships observed at a given site.

Our results suggest that the BBG network depends on the persistence of the selective forces required for the coevolution of these assemblages. We postulate that two forces act in combination to lead to the observed architecture of the BBG network: one, functional redundancy (i.e. the many vertebrates and ticks involved in pathogen circulation), with varying degrees of relative importance; and two, the co-evolutionary overlap of the environmental niches of ticks and vertebrates. We postulate that interactions between ticks and vertebrates were initially derived from their complementary ecological requirements, thereby increasing functional redundancy and allowing the pathogens to circulate. If the environmental dimensions of the interacting species expanded, or if

new species converged towards such niches due to interactions, then other species became part of the evolving network, leading to the observed sub-networks. We hypothesise that sub-networks are the result of secondary niches between ticks and hosts that led to subdivisions of the central core of the network.

We updated the publicly available MLST dataset and produced 704 unique sequences with the aim of checking the links to ticks and/or vertebrates. These results are preliminary due to bias stemming from the reporting of the pathogens in only a few species of ticks and the scarcity of data from vertebrates. However, it is important to notice that (i) *I. ricinus* (western Palaearctic) and *I. scapularis* (eastern Nearctic) do not share any BBG strains; (ii) *I. persulcatus*, *I. pavlovskyi* and *I. granulatus*, species that partially overlap in their Palaearctic distribution range, share several BBG strains; (iii) *I. pacificus* (western Nearctic) share more strains of BBG with ticks of the eastern Palaearctic; and (iv) the strains of BBG in *I. scapularis* are highly repeated. Phylogenetic diversity calculations corroborate the higher diversity of the BBG strains associated with the most surveyed Palaearctic ticks, *I. ricinus* and *I. persulcatus*. The results show many of the same strains in *I. scapularis*, which are slightly related to the *I. ricinus* and *I. persulcatus* strains, and this could be interpreted as the result of a recent BBG invasive event. It is unknown whether this is the result of the geographical segregation of BBG or whether it is an adaptation of BBG strains to the carriers, since the strains were more similar in ticks that are partially sympatric. It has been suggested that BBG adapts to the hosts they find most frequently, resulting in the observed MLST patterns [57]. However, we hypothesise that the ticks could be an important part of the evolutionary pressure exerted on the housekeeping genes of the pathogens. This would be in agreement with previous reports of the tick-borne pathogens *Anaplasma marginale* [59] and *A. phagocytophilum* [26, 60, 61].

## Conclusions

The network constructed in this study, which showed the relationships of ticks, vertebrates and BBG, supports the hypothesis that the pathogens circulate through a nested and recursive network in almost closed clusters. In most of the clusters in which BBG circulates, there were one or more species of generalist ticks high values of centrality within the network that feed on a large range of phylogenetically diverse vertebrates. The specialist ticks in each cluster had low values of centrality and originated nested sub-networks that linked secondary vertebrates to the main nodes of the cluster in an expanding and branched web. Furthermore, vertebrates with higher positions in the range of centrality values tended to be more important in the circulation of BBG because they are the most
prominent tick feeders. We postulate that two forces act in concert to produce the observed architecture of the BBG network: functional redundancy, i.e. the many vertebrates and ticks involved in pathogen circulation, and the co-evolutionary overlap of the environmental niches of ticks and vertebrates. Interactions between ticks and vertebrates were initially derived from their complementary ecological requirements, which increased the functional redundancy and allowed the pathogens to circulate. We hypothesise that sub-networks represent secondary niches between ticks and hosts that lead to subdivisions of the central core of the network due to expansion of the environmental dimensions of the interacting species or due to convergence of new species towards such niches. This functional redundancy may buffer the effects of the disturbance caused by the removal of vertebrates (i.e. the ticks are able to fill remaining biotic niches), thereby forming a highly resilient network. Although the analysis is incomplete, the results we obtained highlight the lack of sharing of BBG strains among the ticks. Further studies focusing on modelling the impact of removing vertebrate phylogenetic diversity on the adaptability of the system are necessary in order to identify a threshold at which the network would collapse.

## Additional files

**Additional file 1:** Complete set of data related to the associations of species of ticks, vertebrates and BBG as obtained by a literature search. (CSV 580 kb)

**Additional file 2:** Phylogenetic tree of the genera of vertebrates in the epidemiological network of *B. burgdorferi* (BBG). The tree was generated using the available data from a mega-phylogeny and was computed and delivered by the Open Tree of Life (http://www.opentreeoflife.com/; accessed 19 Sept 2016). (PDF 50 kb)

**Additional file 3:** The epidemiological network of *B. burgdorferi* (BBG) as visualised using the ForceAtlas2 algorithm at the level of the genera of vertebrates. The clusters are randomly coloured. Each cluster comprises a group of species that are more frequently reported to interact with each other, than with other species. The circles (nodes) represent organisms (species of pathogens and ticks or genera of vertebrates), and the lines with the same colour as a cluster represent interactions among the organisms of the cluster. The size of each node is proportional to its betweenness centrality, and the size of the label is proportional to the PageRank of the node. The width of each link is proportional to the weighted number of interactions between each pair of nodes. Although the network is directed, links lack arrows for improving visualization. (PDF 83 kb)

**Additional file 4:** The epidemiological network of *B. burgdorferi* (BBG) as visualised using the ForceAtlas2 algorithm at the level of the species of vertebrates. The clusters are coloured randomly. Each cluster comprises a group of species that are more frequently reported to interact with each other, either as a parasite or as a transmitter, than with other species. The circles (nodes) represent organisms (species of pathogens and ticks or species of vertebrates), and the lines with the same colour interactions among organisms of the cluster. The size of the node is proportional to its betweenness centrality, and the size of the label is proportional to the PageRank of the node. The width of each link is proportional to the weighted number of interactions between each pair of nodes. Although the network is directed, links lack arrows for improving visualization. (PDF 139 kb)

**Additional file 5:** Complete information for the ticks in the *B. burgdorferi* (BBG) epidemiological network. The table includes values of centrality (betweenness centrality, eigenvector centrality and PageRank), the clustering coefficient, the cluster in which the ticks circulate, the phylogenetic diversity (PD) of the hosts, the richness of the genera of hosts (GR), the circulation capacity (CC) for *B. burgdorferi*, the mean pairwise distances of hosts in the phylogenetic tree of the genera of hosts (MPD) and the significance of the MPD index (MPD-P). Negative MPD values and values of MPD-P near 0.01 indicate strong phylogenetic clustering of the hosts of the tick; positive MPD values and higher MPD-P values indicate weak phylogenetic clustering of the hosts of the tick. (CSV 15 kb)

**Additional file 6:** Complete information for the *Borrelia* species in the *B. burgdorferi* (BBG) epidemiological network. The table includes values of centrality (betweenness centrality, eigenvector centrality and PageRank), the clustering coefficient, the cluster in which the ticks circulate, the phylogenetic diversity (PD) of the hosts, the richness of the genera of hosts (GR), the mean pairwise distances of hosts in the phylogenetic tree of the genera of hosts (MPD) and the significance of the MPD index (MPD-P). Negative MPD values and values of MPD-P near 0.01 indicate strong phylogenetic clustering of the hosts of the tick; positive MPD values and higher MPD-P values indicate weak phylogenetic clustering of the hosts of the tick. (CSV 7 kb)

**Additional file 7:** Complete information for the genera of vertebrates (hosts and/or reservoirs) in the *B. burgdorferi* epidemiological network. The table includes values of centrality (betweenness centrality, eigenvector centrality and PageRank), the clustering coefficient, and the cluster in which the organisms circulate. (CSV 82 kb)

**Additional file 8:** The records of ticks and species of *B. burgdorferi* (BBG) in the phylogenetic tree of the genera of vertebrates (see Additional file 2) for the significant clusters detected in the epidemiological network. Each blue circle represents a record of either a tick or a species of BBG on the host. A: cluster 1; B: cluster 2; C: cluster 5; D: cluster 6. (PDF 46 kb)

**Additional file 9:** The phylogenetic tree for the concatenated multilocus sequence typing (MLST) sequences of *B. burgdorferi* (BBG) as recorded on the most represented ticks in the MLST dataset. The tip of each branch is labelled with the MLST classification. Each blue dot indicates that it was found in the reference tick. The size of the dot is proportional to the number of times that the specific MLST type was recorded on that tick. (PDF 272 kb)

**Additional file 10:** The phylogenetic tree for the concatenated multilocus sequence typing (MLST) sequences of *B. burgdorferi* (BBG) as recorded on the most represented ticks in the MLST dataset. The tip of each branch is labelled with the MLST classification. Each blue dot indicates that it was found in the reference tick. The size of the dot is proportional to the number of times that the specific MLST type was recorded on that tick. (PDF 272 kb)

**Additional file 11:** The phylogenetic tree for the concatenated multilocus sequence typing (MLST) sequences of *B. burgdorferi* (BBG) as recorded on the most represented ticks in the MLST dataset. The tip of each branch is labelled with the MLST classification. Each blue dot indicates that it was found in the reference tick. The size of the dot is proportional to the number of times that the specific MLST type was recorded on that tick. (PDF 272 kb)

**Additional file 12:** The phylogenetic tree for the concatenated multilocus sequence typing (MLST) sequences of *B. burgdorferi* (BBG) as recorded on the most represented ticks in the MLST dataset. The tip of each branch is labelled with the MLST classification. Each blue dot indicates that it was found in the reference tick. The size of the dot is proportional to the number of times that the specific MLST type was recorded on that tick. (PDF 272 kb)

**Additional file 13:** The phylogenetic tree for the concatenated multilocus sequence typing (MLST) sequences of *B. burgdorferi* (BBG) as recorded on the most represented ticks in the MLST dataset. The tip of each branch is labelled with the MLST classification. Each blue dot indicates that it was found in the reference tick. The size of the dot is proportional to the number of times that the specific MLST type was recorded on that tick. (PDF 272 kb)
**Additional file 14:** The network of strains of *B. burgdorferi* (BBG), ticks and vertebrates as characterised by multilocus sequence typing (MLST) analysis using the ForceAtlas2 algorithm to determine the clusters and retaining only the strains recorded from questing ticks or hosts. This figure complements the data in Fig. 4. Each number represents a MLST strain that is linked to either the tick or the vertebrate where it was recorded. The width of the link is proportional to the weighted number of times that the link was recorded. The label size is proportional to the centrality of the node as measured using PageRank. The name of the species of BBG is also included in the link. The colour of the link represents the species of BBG. (PDF 4212 kb)

**Abbreviations**

BBG: *Borrelia burgdorferi* (s.l.) group; BNC: Betweenness centrality; CC: Circulation capacity; MLST: Multilocus sequence typing; CCL: Chemokine ligand; PD: Faith's phylogenetic diversity; PR: PageRank; WCC: Weighted clustering coefficient; WD: Weighted degree

**Acknowledgements**

Not applicable.

**Funding**

Parts of this work were performed under the umbrella of the COST project (TD1303). The funding agency did not have any role in the design of the study, in the collection, analysis and interpretation of data, or in the writing of the manuscript.

**Availability of data and material**

The datasets that support the conclusions of this article are included in the article and its Additional files.

**Authors' contributions**

AEP designed the study, supervised the literature search and record collection and performed the network analyses. ACC performed the calculations used to construct the phylogenetic trees and to determine phylogenetic diversity. AR performed the literature search. ECC and HS obtained and sequenced the strains of *Borrelia* and performed the additional MLST sequence analyses for the pathogens included in this study. AEP and JF analysed the data and wrote the paper. All authors read and approved the final manuscript.

**Competing interests**

The authors declare that they have no competing interests.

**Consent for publication**

Not applicable.

**Ethics approval and consent to participate**

Not applicable.

**Author details**

¹Department of Animal Health, Faculty of Veterinary Medicine, University of Zaragoza, Miguel Servet, 177, Zaragoza, Spain. ²National Institute for Public Health and the Environment (RIVM), Antonie van Leeuwenhoeklaan 9, 3721 MA Bilthoven, The Netherlands. ³Center for Infection and Immunity of Lille (CILL), INSERM U1019 – CNRS UMR 8204, Université Lille Nord de France, Institut Pasteur de Lille, Lille, France. ⁴SaBio, Instituto de Investigación de Recursos Cinegéticos IREC-CSIC-UCLM-JCCM, Ronda de Toledo s/n, 13005 Ciudad Real, Spain. ⁵Department of Veterinary Pathobiology, Center for Veterinary Health Sciences, Oklahoma State University, Stillwater, OK 74078, USA.

Received: 23 June 2016 Accepted: 14 September 2016

Published online: 23 September 2016

**References**

1. Margos G, Vollmer SA, Ogden NH, Fish D. Population genetics, taxonomy, phylogeny and evolution of *Borrelia burgdorferi* sensu lato. *Infect Gen Evol.* 2011;11:1545–63.

2. Ivanova LB, Tomova A, González-Acuña D, Murúa R, Moreno CX, Hernández C, et al. *Borrelia chilensis*, a new member of the *Borrelia burgdorferi* sensu lato complex that extends the range of this genospecies in the Southern Hemisphere. *Environ Microbiol.* 2014;16:1069–80.

3. Kurtenbach K, De Michelis S, Etti S, Schäfer SM, Sewell HS, Brade V, et al. Host association of *Borrelia burgdorferi* sensu lato—the key role of host complement. *Trends Microbiol.* 2002;10:74–9.

4. Wang G, Van Dam AP, Le Fleche A, Postic D, Peter O, Baranton G, et al. Genetic and phenotypic analysis of *Borrelia valaisiana* sp. nov. (*Borrelia* genomic groups VS116 and M19). *Int J Syst Evol Microbiol.* 1997;47:926–32.

5. Pritt BS, Mead PS, Johnson DKH, Neitzel DF, Respicio-Kingry LB, Davis JP, et al. Identification of a novel pathogenic *Borrelia* species causing Lyme borreliosis with unusually high spirochaetaemia: a descriptive study. *Lancet Infect Dis.* 2016;16:556–64.

6. Kurtenbach K, Hanincová K, Tsao JJ, Margos G, Fish D, Ogden NH. Fundamental processes in the evolutionary ecology of Lyme borreliosis. *Nat Rev Microbiol.* 2006;4:660–9.

7. Gatewood AG, Liebman KA, Vourc'h G, Bunikis J, Hamer SA, Cortinas R, et al. Climate and tick seasonality are predictors of *Borrelia burgdorferi* genotype distribution. *Appl Environ Microbiol.* 2009;75:2476–83.

8. Rudenko N, Golovchenko M, Grubhoffer L, Oliver JH. Updates on *Borrelia burgdorferi* sensu lato complex with respect to public health. *Ticks Tick-borne Dis.* 2011;2:123–8.

9. Margos G, Gatewood AG, Aanensen DM, Hanincová K, Terekhova D, Vollmer SA, et al. MLST of housekeeping genes captures geographic population structure and suggests a European origin of *Borrelia burgdorferi*. *PNAS.* 2008; 105:8730–5.

10. Castillo-Ramírez S, Fingerle V, Jungnick S, Straubinger RK, Krebs S, Blum H, et al. Trans-Atlantic exchanges have shaped the population structure of the Lyme disease agent *Borrelia burgdorferi* sensu stricto. *Sci Rep.* 2015;6:22794.

11. Coipan EC, Fonville M, Tijssé-Klasen E, van der Giessen JW, Takken W, Sprong H, et al. Geodemographic analysis of *Borrelia burgdorferi* sensu lato using the 5S-23S rDNA spacer region. *Infect Genet Evol.* 2013;17:216–22.

12. Wang G, Van Dam AP, Schwartz I, Dankert J. Molecular typing of *Borrelia burgdorferi* sensu lato: taxonomic, epidemiological, and clinical implications. *Clin Microbiol Rev.* 1999;12:633–53.

13. Ehrlich PR, Raven PH. Butterflies and plants: a study in coevolution. *Evolution.* 1964;18:586-608.

14. Thompson JN. The geographic Mosaic of coevolution. *Chicago: Univ. of Chicago Press*; 2005.

15. Bascompte J, Jordano P, Olesen JM. Asymmetric coevolutionary networks facilitate biodiversity maintenance. *Science.* 2006;312:431–3.

16. Thompson JN. The coevolutionary process. *Chicago: University of Chicago Press*; 1994.

17. Bascompte J, Jordano P, Melián CJ, Olesen JM. The nested assembly of plant-animal mutualistic networks. *Proc Natl Acad Sci USA.* 2003;100:9383–7.

18. Jordano P, Bascompte J, Olesen JM. Invariant properties in coevolutionary networks of plant-animal interactions. *Ecol Lett.* 2003;6:69–81.

19. Vázquez DP, Aizen MA. Asymmetric specialization: a pervasive feature of plant-pollinator interactions. *Ecology.* 2004;85:1251–7.

20. Streicker DG, Fenton A, Pedersen AB. Differential sources of host species heterogeneity influence the transmission and control of multihost parasites. *Ecol Lett.* 2013;16:975–84.

21. Poulin R. Network analysis shining light on parasite ecology and diversity. *Trends Parasitol.* 2010;26:492–8.

22. Rigaud T, Perrot-Minnot M-J, Brown MJF. Parasite and host assemblages: embracing the reality will improve our knowledge of parasite transmission and virulence. *Proc R Soc London B Biol Sci.* 2010;277:3693–702.

23. Estrada-Peña A, de La Fuente J, Ostfeld RS, Cabezas-Cruz A. Interactions between tick and transmitted pathogens evolved to minimise competition through nested and coherent networks. *Sci Rep.* 2015;5:10361.

24. Neumann MEJ. Networks: an introduction. New York: Oxford University Press; 2010.

25. Krasnov BR, Fortuna MA, Mouillot D, Khokhlova IS, Shenbrot GI, Poulin R. Phylogenetic signal in module composition and species connectivity in compartmentalized host-parasite networks. *Am Nat.* 2012;179:501–11.

26. de la Fuente J, Estrada-Peña A, Cabezas-Cruz A, Brey R. Flying ticks: anciently evolved associations that constitute a risk of infectious disease spread. *Parasit Vectors.* 2015;8:1–11.

27. de la Fuente J, Estrada-Peña A, Cabezas-Cruz A, Kocan KM. *Anaplasma phagocytophilum* uses common strategies for infection of ticks and vertebrate hosts. *Trends Microbiol.* 2015;24:173–80.
28. Pilosof S, Fortuna MA, Cosson J-F, Galan M, Kittipong C, Ribas A, et al. Host-parasite network structure is associated with community-level immunogenetic diversity. *Nat Commun.* 2014;5:1–9.

29. Ogden NH, Nuttall PA, Randolph SE. Natural Lyme disease cycles maintained via sheep by co-feeding ticks. *Parasitology*. 1997;115:591–9.

30. Piesman J, Happ CM. The efficacy of co-feeding as a means of maintaining *Borrelia burgdorferi*: a North American model system. *J Vector Ecol.* 2001;26:216–20.

31. Nava S, Barbieri AM, Maya L, Colina R, Mangold AJ, Labruna MB, et al. *Borrelia* infection in *Ixodes pararicinus* ticks (Acari: Ixodidae) from northwestern Argentina. *Acta Trop.* 2014;139:1–4.

32. Gómez JM, Nunn CL, Verdú M. Centrality in primate-parasite networks reveals the potential for the transmission of emerging infectious diseases to humans. *Proc Natl Acad Sci U S A.* 2013;110:7738–41.

33. Barrat A, Barthelemy M, Pastor-Satorras R, Vespignani A. The architecture of complex weighted networks. *Proc Natl Acad Sci USA.* 2004;101:3747–52.

34. Craft ME, Caillaud D. Network models: an underutilized tool in wildlife epidemiology? *Interdiscip Perspect Infect Dis.* 2011;2011:676949.

35. Csardi G, Nepusz T. The igraph software package for complex network research. *Int Journal Complex Systems* 2006;1695. http://igraph.org/. Accessed 19 Sept 2016.

36. R Core Team. R: A language and environment for statistical computing. Vienna: R Foundation for Statistical Computing; 2014. http://www.R-project.org/. Accessed 19 Sept 2016.

37. Jacomy M, Heymann S, Venturini T, Bastian M. ForceAtlas2, a continuous graph layout algorithm for handy network visualization. *Gephi's* documentation manual, available at http://wiki.gephi.org/index.php/Main_Page. Accessed Dec 2012.

38. Kourtellis N, Alahakoon T, Simha R, Iamnitchi A, Tripathi R. Identifying high betweenness centrality nodes in large social networks. *Soc Netw Anal Min.* 2013;3:899–914.

39. Blondel VD, Guillaume JL, Lambiotte R, Lefebvre E. Fast unfolding of communities in large networks. *J Stat Mech: Theory and Experiment*. 2008. Online at stacks.iop.org/JSTAT/2008/P10008. Accessed Dec 2014.

40. Faith DP. Conservation evaluation and phylogenetic diversity. *Biol Conserv.* 1992;61:1–10.

41. Webb CO, Ackerly DD, Mcpeak MA, Donoghue MJ. Phylogenies and community ecology. *Ann Rev Ecol Syst.* 2002;33:475–505.

42. Coipan EC, Jahfari S, Fonville M, Oei GA, Spanjaard L, Takumi K, Hovius JW, et al. Imbalanced presence of *Borrelia burgdorferi* s.l. multilocus sequence types in clinical manifestations of Lyme borreliosis. *Infect Genet Evol.* 2016; 42:66–76.

43. Tamura K, Dudley J, Nei M, Kumar S. MEGA4: molecular evolutionary genetics analysis (MEGA) software version 4.0. Mol Biol Evol. 2007;24:1596–9.

44. Devictor V, Clavel J, Julliard R, Lavergne S, Mouillot D, Thuiller W, et al. Defining and measuring ecological specialization. *J Appl Ecol.* 2010;47:15–25.

45. Hassell M, May R. Generalist and specialist natural enemies in insect predator-prey interactions. *J Anim Ecol.* 1986;55:923–40.

46. Bastolla U, Fortuna MA, Pascual-García A, Ferrera A, Luque B, Bascompte J. The architecture of mutualistic networks minimizes competition and increases biodiversity. *Nature*. 2009;458:1018–20.

47. Wagemakers A, Staarink PJ, Sprong H, Hovius JW. *Borrelia miyamotoi*: a widespread tick-borne relapsing fever spirochete. *Trends Parasitol.* 2015;31:260–9.

48. Gatzmann F, Metzler D, Krebs S, Blum H, Sing A, Takano A, et al. NGS population genetics analyses reveal divergent evolution of a Lyme borreliosis agent in Europe and Asia. *Ticks Tick-borne Dis.* 2015;6:344–51.

49. Dubska L, Literak I, Kocianova E, Taragelova V, Sychra O. Differential role of passerine birds in distribution of *Borrelia* spirochetes, based on data from ticks collected from birds during the postbreeding migration period in Central Europe. *Appl Environ Microbiol.* 2009;75:596–602.

50. Heylen D, Tijssé E, Fonville M, Matthysen E, Sprong H. Transmission dynamics of *Borrelia burgdorferi* s.l in a bird tick community. *Environ Microbiol.* 2013;15:663–73.

51. Shaw MT, Keesing F, McGrail R, Ostfeld RS. Factors influencing the distribution of larval blacklegged ticks on rodent hosts. *Am J Trop Med Hyg.* 2003;68:447–52.

52. Barbour AG, Bunikis J, Fish D, Hanincová K. Association between body size and reservoir competence of mammals bearing *Borrelia burgdorferi* at an endemic site in the northeastern United States. *Parasit Vectors.* 2015;8:299.

53. van Duijvendijk G, Sprong H, Takken W. Multi-trophic interactions driving the transmission cycle of *Borrelia afzelii* between *Ixodes ricinus* and rodents: a review. *Parasit Vectors.* 2015;8:1.

54. Devevey G, Dang T, Graves CJ, Murray S, Brisson D. First arrived takes all: inhibitory priority effects dominate competition between co-infecting *Borrelia burgdorferi* strains. *BMC Microbiol.* 2015;15:1.

55. Roquet C, Thuiller W, Lavergne S. Building megaphylogenies for macroecology: taking up the challenge. *Ecography*. 2013;36:13–26.

56. Smith SA, Beaulieu JM, Donoghue MJ. Mega-phylogeny approach for comparative biology: an alternative to supertree and supermatrix approaches. *BMC Evol Biol.* 2009;9:37–49.

57. Hofmeester TR, Coipan EC, van Wieren SE, Prins HHT, Takken W, Sprong H. Few vertebrate species dominate the *Borrelia burgdorferi* s.l. life cycle. *Environ Res Lett.* 2016;11:1–16.

58. Huang ZYX, de Boer WF, van Langevelde F, Olson V, Blackburn TM, Prins HHT. Species' life-history traits explain interspecific variation in reservoir competence: a possible mechanism underlying the dilution effect. *PLoS One*. 2013;8:e54341.

59. Estrada-Peña A, Naranjo V, Acevedo-Whitehouse K, Mangold AJ, Kocan KM, de la Fuente J. Phylogeographic analysis reveals association of tick-borne pathogen, *Anaplasma marginale*, MSP1a sequences with ecological traits affecting tick vector performance. *BMC Biol.* 2009;7:1.

60. Rar VA, Epikhina TI, Yakimenko VV, Malkova MG, Tancev AK, Bondarenko EI, et al. Genetic variability of *Anaplasma phagocytophilum* in ticks and voles from *Ixodes persulcatus/ Ixodes trianguliceps* sympatric areas from Western Siberia, Russia. *Ticks Tick-borne Dis.* 2014;5:854–63.

61. Jahfari S, Coipan EC, Fonville M, Van Leeuwen AD, Hengeveld P, Heylen D, et al. Circulation of four *Anaplasma phagocytophilum* ecotypes in Europe. *Parasit Vectors.* 2014;7:365.

Submit your next manuscript to BioMed Central and we will help you at every step:

* We accept pre-submission inquiries

* Our selector tool helps you to find the most relevant journal

* We provide round the clock customer support

* Convenient online submission

* Thorough peer review

* Inclusion in PubMed and all major indexing services

* Maximum visibility for your research

Submit your manuscript at
www.biomedcentral.com/submit

BioMed Central