Cryptic vector divergence masks vector-specific patterns of infection: an example from the marine cycle of Lyme borreliosis.

Vector organisms are implicated in the transmission of close to a third of all infectious diseases. In many cases, multiple vectors (species or populations) can participate in transmission but may contribute differently to disease ecology and evolution. The presence of cryptic vector populations can be particularly problematic as differences in infection can be difficult to evaluate and may lead to erroneous evolutionary and epidemiological inferences. Here, we combine site-occupancy modeling and molecular assays to evaluate patterns of infection in the marine cycle of Lyme borreliosis, involving colonial seabirds, the tick Ixodes uriae, and bacteria of the Borrelia burgdorferi s.l. complex. In this cycle, the tick vector consists of multiple, cryptic (phenotypically undistinguishable but genetically distinct) host races that are frequently found in sympatry. Our results show that bacterial detection varies strongly among tick races leading to vector-specific biases if raw counts are used to calculate Borrelia prevalence. These differences are largely explained by differences in infection intensity among tick races. After accounting for detection probabilities, we found that overall prevalence in this system is higher than previously suspected and that certain vector-host combinations likely contribute more than others to the local dynamics and large-scale dispersal of Borrelia spirochetes. These results highlight the importance of evaluating vector population structure and accounting for detection probability when trying to understand the evolutionary ecology of vector-borne diseases.

Introduction

Vectors are key components in the transmission of many important infectious diseases, but the spatial and temporal dynamics of many vector-borne pathogens still remain poorly understood. This is particularly true for vector-borne diseases with spatially extensive distributions and where multiple vector species participate in transmission (e.g. malaria, Lyme disease, dengue, etc.). In these systems, different vector types can vary in their contribution to local infection dynamics. Indeed, the relative role of each vector will depend on its capacity to acquire, maintain, and transmit the pathogen (i.e. vector competence; Lane 1994; Edman 2002). Likewise, the distribution, diversity, and abundance of vector species will strongly affect disease risk (Kurtenbach et al. 2006). For example, at the within-species scale, the subdivision of the vector into discrete populations (i.e. population structure) can affect the temporal and spatial scale at which transmission occurs (McCoy 2008). Although rarely tested, several cases of cryptic fine scale population structure in vector organisms have now been demonstrated (e.g. (Ravel et al. 2002; McCoy et al. 2003; Ravel et al. 2007; Kempf et al. 2009a,b), as has population-specific variation in vector competence (e.g. Lambrechts et al. 2005; Joy et al. 2008). In some cases, this structure can involve different sympatric vector populations that exploit different local hosts (i.e. host races; McCoy et al. 2005; Whiteman et al. 2006; Johannesen et al. 2008). As vector population structure can influence disease dynamics and evolution, patterns of infection need to be considered explicitly in such cases.

Lyme disease or Lyme borreliosis (LB), caused by spirochetes of the Borrelia burgdorferi species complex, is the most common vector-borne disease in the northern hemisphere and has major medical and economic impacts (Stanek and Strle 2003; Steere et al. 2004). Lyme disease spirochetes are maintained in complex natural transmission cycles including numerous tick vector and reservoir host species. Previous work suggested that the generalist nature of the common tick vectors (Ixodes spp.) acts as a homogenizing mechanism linking the different ecological niches of the bacteria (Kurtenbach et al. 2006). However, recent evidence shows that host specialization and cryptic divergence (i.e. hidden genetic structure) may frequently occur in ticks (McCoy et al. 2001, 2005; Kempf et al. 2009b; De Meeûs et al. unpublished data). Whereas much work has focused on the reservoir hosts of LB bacteria (e.g. LoGiudice et al. 2003, 2008; Brisson et al. 2008), less is known about tick population structure and associated dynamics, or the importance of this structure for the evolutionary ecology and epidemiology of Lyme disease (Qiu et al.2002; De Meeûs et al. 2004).

The marine cycle of LB was discovered more than a decade ago when Olsen et al. (1993) demonstrated the circulation of one of the pathogenic species of the LB complex, Borrelia garinii, among seabird reservoirs via the tick Ixodes uriae. Although current evidence suggests a role for seabirds in the global circulation of Lyme disease bacteria (Olsen et al. 1995; Olsen 2003, 2007; Duneau et al. 2008; Staszewski et al. 2008), the importance of the marine cycle in Lyme disease epidemiology has yet to be explicitly considered (e.g. Kurtenbach et al. 2006). As is still the case for the terrestrial tick vectors of Lyme disease bacteria (Ixodes ricinus species complex), I. uriae was traditionally considered as a generalist ectoparasite exploiting numerous species of colonial seabirds (Rothschild and Clay 1961; Wilson 1970; Guiguen 1988). However, genetic work has revealed that this tick consists of locally distinct seabird species-specific host races, which are often found in sympatry (McCoy et al. 2001, 2005). That is, McCoy et al. (2001, 2005) have demonstrated higher genetic differentiation between ticks from different sympatric seabird host species than between ticks from allopatric populations of the same seabird species. Furthermore, host-associated divergence in I. uriae seems to be a recurrent and spatially dynamic process, such that tick races have evolved recently and independently over different geographic regions (Kempf et al. 2009a). This system therefore provides a pertinent biological model to examine differences among vector races in pathogen transmission and the consequences of race evolution for the evolutionary ecology of Lyme disease bacteria.

Here, we examine the implications of host-associated genetic structure in this vector system by examining differences in local pathogen infection among different sympatric host races. Evaluating prevalence is often challenging due, in part, to difficulties in diagnosing infection, i.e. imperfect detection. Moreover, pathogen detection probabilities may be population-specific leading to vector-specific biases in estimated parameters and thus to erroneous epidemiological inferences. We therefore employ recently developed site-occupancy models (Mackenzie et al. 2002), to explicitly test for differences in pathogen detection probabilities among different vector races and to estimate vector-specific prevalence taking into account these detection probabilities. We evaluate the potential source of imperfect detection by combining statistical modeling with a quantitative PCR assay. Finally, we discuss these results in terms of their importance for the disease dynamics of this pathogen and for other vector-borne systems.

Materials and methods

Tick sampling

The tick sample included 58 adult ticks of unknown infection status collected from three different host races (ticks associated with Black-legged kittiwakes Rissa tridactyla, n = 25; Atlantic puffins Fratercula arctica, n = 15; and Common murres Uria aalge, n = 18) in four geographic locations [three large seabird colonies in Iceland: Skrudur (64°85′40″N, 13°83′80″W), Grimsey (66°83′30″N, 18°80′00″W), and Breidafjordur (65°82′30″N, 22°85′40″W), and one in Norway: Hornøya (70°82′20″N, 31°81′00″E)] (Table S1). In addition, we included 37 adult ticks that were ‘known positives’ from prior serological tests and that were sampled from the same colonies and host species (kittiwake ticks = 7, puffin ticks = 16, murre ticks = 14; Staszewski et al. 2008). These 37 ‘known positives’ are assumed to have the same detection probabilities as the other 58 ticks, and thus were used to improve the estimation of detection probabilities, but not in estimating prevalence. All ticks were sampled from as many individual host birds as possible (average 1.3 ± 0.1 ticks per bird) and were stored using the same procedures. We assumed that these ticks were representative samples from their populations and that their health-status did not affect their probability of being sampled (i.e. infected ticks were as likely to be collected as uninfected ticks).

PCR detection of Borrelia

Borrelia prevalence was estimated by a standard nested PCR procedure that involves amplifying a portion of the flagellin B (flaB) gene using DNA extractions of individual ticks (Johnson et al. 1992; Clark 2004). A positive infection was revealed by the presence of an amplification product on an agarose gel. All amplified products were subsequently sent for direct sequencing (see Duneau et al. 2008 for protocols).

Quantification of Borrelia infection

We applied a highly sensitive target-specific qPCR procedure to quantify relative differences in infection among the ticks collected from the three different seabird hosts. The quantitative real-time PCR assay utilized the LightCycler platform and SYBR Green I detection system (Applied Biosystems, Courtaboeuf, France). The tick sample included the same set of 95 tick extracts analyzed by nested PCR. Given that ticks were not homogeneously engorged and DNA concentrations differed among tick extracts, DNA templates were adjusted to a stock concentration of 40–50 ng/μL. For this analysis, we designed forward and reverse primers to amplify a conserved 214-pb region within the flagellin gene (qFlaB-F 5′-GGAATGCAACCTGCAAAAAT-3′ and qFlaB-R 5′-GGCTGTTGAGCTCCTTCTTG-3′). Based on the genome sequence of B. burgdorferi type strain B31, it was assumed that only a single copy of flaB was present per spirochete (Fraser et al. 1997). Reactions were performed in a 10 μL volume containing 1× FastStart Taq DNA polymerase mixture (Applied Biosystems), 4 mm MgCl2, 1 μm concentrations of each primer, and 1 μL of tick DNA (stock dilution of 40–50 ng/μL). The plasmid pB31/41-9 (kindly provided by Dr R. Wallich, Department of Applied Immunology, Heidelberg, Germany), containing the flagellin gene of B. burgdorferi B31 (X16833), was used as a quantification standard. A standard curve was generated with serial dilutions from 5 to 500 000 copies of the reference plasmid pB31/41-9. Two negative controls were included each time. The amplification program included the initial activation step of the FastStart Taq polymerase at 95°C for 10 min and 50 cycles of denaturation at 95°C for 10 s, annealing at 60°C for 15 s, and extension at 72°C for 30 s. Fluorescence was measured at the end of each extension step. After amplification, a melting curve was acquired by heating the product at 20°C/s to 95°C, cooling it at 20°C/s to 60°C, keeping it at 60°C for 20 s, and then slowly reheating it at 0.1°C/s to 95°C. Fluorescence was measured through the slow reheating phase.

The relative number of spirochetes in each PCR mixture was calculated by comparing the crossing points of the samples with those of the standards using the Lightcycler software. Specificity of amplification products were tested by comparing the melting curves and the mean Tm values of the samples with those of different reference strains [B. garinii 20047 (GenBank accession no. D82846), B. burgdorferi B31 (X15661), Borrelia lusitaniae PotiB1 (DQ111035), and Borrelia afzelii VS461 (L29237)]. Post-PCR analyses were also performed using conventional gel electrophoresis followed by direct sequencing of all amplification products (Cogenics, Meylan, France).

Reliable and reproducible quantification of infection was achieved over five orders of magnitude, from 5 to 500 000 plasmid copies. Each sample was run up to four times. Fifty-four of the 95 adult ticks tested positive by qPCR and all runs were highly repeatable (Pearson R2 = 0.98, P < 0.001). We averaged the calculated number of spirochetes from each run for each sample and used this number as an index of the quantity of bacteria in each individual tick (i.e. infection intensity: Margolis et al. 1982).

Site-occupancy models: ‘naive’ versus ‘estimated’ prevalence

To determine the influence of imperfect detection, i.e. false negatives, on PCR prevalence estimates, we applied a site-occupancy modeling framework to the dataset (Mackenzie et al. 2002, 2006). Although more frequently applied in ecological modeling (e.g. Bailey et al. 2004; Altwegg et al. 2008; Weller 2008; Winchell and Doherty 2008), site-occupancy, and related, models are rarely incorporated into parasitological or epidemiological studies (Jennelle et al. 2007; Thompson 2007; Conn and Cooch 2009). Site-occupancy studies rely on a re-sampling of sites that are likely to harbor the species of interest to estimate the proportion of sites occupied, where it is explicitly acknowledged that species may often be missed during a survey even when they are present. Detection probabilities are estimated from the repeated visits. In our study, each tick is a ‘site’ that the pathogen may occupy and re-sampling entails testing an individual vector several times for the presence of infection. Patterns of presence and apparent absence over multiple tests are used to estimate the probability of detection (p) of the pathogen or, in other words, the sensitivity of the test. Essentially the apparent (or ‘naive’) prevalence is corrected for the detection probability and an unbiased estimate of the number of vector individuals occupied by the pathogen (ψ), i.e. estimated prevalence, is obtained.

Using the samples described above, we repeated nested PCRs five times on all 95 individuals, one PCR plate per day for 5 days using the same tick DNA extractions, PCR reagents, and thermocycler. From this data, we constructed histories for each sampled tick where each tick was analogous to a single site and each nested PCR test for the presence of Borrelia to a site visit in a standard site-occupancy analysis. For instance, a tick with a history (h) of 10101 indicates that the bacterium was detected on the first, third, and fifth PCR, but not on the second and fourth PCR. Thus, this tick was occupied by the bacterium (ψ), but the probability of detection was not equal to 1. A probability statement associated with this history would be

If no bacterium was detected in the five PCRs, then the history would be 00000. In this case, either the tick could be occupied by the bacterium, but it was not detected, or the tick could be unoccupied by the bacterium (i.e. uninfected). Thus, the probability statement associated with this history is

Once the probability statement for each of the s observed histories is formed, the model likelihood for the observed data can be constructed, such that

and parameters can be estimated using maximum likelihood methods. Parameters can also be modeled as a function of covariates using a logit link. For example, occupancy can be modeled as a function of an index of bacteria quantity such that logit(ψ) = β0 + β1(quantity), where β0 is the y-intercept and β1 is the slope parameter associated with the quantity of bacteria. In our model set below, we use notation such that this structure would be designated as ψquantity. If we consider a parameter to be constant (equivalent to having only a y-intercept term in a model), then we use the label ‘constant’. Because we had a sample of known positive ticks, we were able to share information on the detection probability from this sample with the unknown prevalence sample by fixing prevalence for these ticks to 1 and modeling the detection parameters in common across the two groups of ticks. We constructed a number of models (see Model set section) to test hypotheses and explain variation in the detection and prevalence parameters.

Model set

We modeled prevalence (ψ) as constant across all tick races (constant) or as differing for each seabird-specific tick race (tick race). We modeled the detection probability (p) as a constant (constant), as a function of tick race (tick race), and as a function of the bacteria quantity index using both linear (quantity) and quadratic (quantity2) relationships. We also included an effect of engorgement status (engorge) for detection in each tick. Ixodes uriae feed only once per life stage, thus engorgement status was considered as a binomial variable designating whether the tick had had a blood meal when collected (0: nonengorged, 1: engorged). We expected that detection would be higher in engorged ticks as bacteria start to multiply in the midgut of the tick vector once the blood meal is engaged, thus increasing the overall spirochete density in the tick (Piesman et al. 1990). We ran all combinations of additive effects and constructed a total of 24 alternative models.

Model selection and parameter estimation were performed using Program mark (White and Burnham 1999; http://welcome.warnercnr.colostate.edu/~gwhite/mark/mark.htm). We calculated a value (i.e. overdispersion factor) using the median procedure available in Program mark. The can help correct for overdispersion in the data and improve model goodness-of-fit. We used Akaike's Information Criterion, adjusted for and with a small sample size correction (QAICc) to rank models (Burnham and Anderson 2002). We used QAICc weights to help interpret the strength of evidence for each model. These model weights sum to 1. We also computed cumulative QAICc weights (0 ≤ ∑ w ≤ 1), or importance values, for each variable. Cumulative weights evaluate the strength of evidence for each variable and are a sum of the model weights associated with each model in which the variable appears. We also present model-averaged prevalence estimates across the entire model set to include model uncertainty in our estimates (Burnham and Anderson 2002; Anderson 2008b).

Results

Borrelia detection and infection intensity

Of the 95 I. uriae ticks (both of known and unknown infection status) subjected to single-nested PCR and qPCR analyses, 29 (30.5%) ticks were positive and 39 (41.1%) negative by both methods. The remaining 27 (28.4%) ticks were positive only by qPCR. The sequence analysis of Borrelia isolates revealed the majority of strains corresponded to the B. garinii genospecies. However, some atypical strains were detected and grouped most closely with B. lusitaniae in phylogenetic analyses (Duneau et al. 2008).

The index of the number of spirochetes in the tick DNA extracts used in PCR varied from 0 to a maximum value of 1957. Among infected ticks, this index showed a median of 2.97 and a geometric mean of 5.01 ± 0.53 (no. bacteria ± SE). These are relative index values of the total number of spirochetes per tick because quantitative analyses were performed on standardized tick DNA concentrations (40–50 ng of DNA). These values cannot be easily extrapolated to the entire tick because spirochetes are not homogeneously distributed in the tick's body and tick DNA concentrations in the extracts change with the engorgement status. However, as there was no strong association between host type and tick engorgement status (Spearman's rho = 0.02, N = 91, P = 0.82), we assume this underestimation would not affect our overall results (i.e. testing differences among host-associated tick races in infection prevalence and intensity). Overall, the relative infection intensity was greater in puffin ticks (148.89 ± 90.57) than in murre (52.47 ± 45.74) and kittiwake ticks (13.50 ± 6.83) (mean no. bacteria ± SE) (Fig. 1).

Figure 1

Intensity of infection of Lyme borreliosis bacteria (mean number of spirochetes per infected sample) in three tick races associated with puffins, murres, and kittiwakes, respectively. The error bars are standard errors (SE).

Borrelia detection probabilities

We estimated to be 1.7 and adjusted our model set accordingly. Model-selection results (Table 1) showed that Borrelia detection was primarily a function of tick race (cumulative QAICc wt = 1.0 (βmurre = −2.44, 95% CI: −3.71,−1.16; βkittiwake = −1.56, 95% CI: 0.01, 0.41) and bacteria quantity (quantity2; cumulative QAICc weight = 1.0, Fig. 2) and, to a lesser extent, engorgement status (cumulative QAICc weight = 0.31). The best model that contained an engorgement effect, predicted a positive relationship (β = 0.46, 95% CI: −0.53, 1.46), but the 95% confidence interval widely overlapped zero and was probably a ‘pretending variable’ (Anderson 2008a). Such variables rank highly because they are only one parameter different from other high-ranking models and by definition AIC will only drop by a maximum of 2 with an addition of a single parameter. The bacteria quantity index was also positive (, 95% CI: 0.01, 0.04) and differed by tick race (Fig. 2). The puffin tick race had the highest detection probability and the strongest relationship with bacteria quantity (Fig. 2A), whereas the murre tick race had the lowest (Fig. 2C).

Table 1

Model selection results for models estimating prevalence (ψ) of Borrelia spp. in Ixodes uriae host races

Model	ΔQAICc	QAICc weights	Model likelihood	Number of parameters	Deviance
	0.00	0.48	1.00	8	168.68
	1.59	0.21	0.45	9	167.83
	1.65	0.21	0.44	6	175.05
	3.14	0.10	0.21	7	174.21
	13.16	0.00	0.00	5	188.84

We modeled Borrelia occupancy of ticks as a constant regardless of which seabird species they fed on (ψconstant) and as a function of tick race (ψtickrace). We modeled detection as a constant (pconstant), a function of tick race (ptickrace), a function of the relative quantity of bacteria in a linear (pquantity), and quadratic () fashion, as well as a function of engorgement status (pengorge). We ran all combinations of additive effects. Only the top five models are presented. ΔQAICc is the QAICc scaled so the lowest value is zero (for ease of interpretation). QAICc weights are the model weights. The model likelihood is the QAICc weight of the model of interest divided by the QAICc weight of the best model. This value is the strength of evidence of one model relative to other models in the set. Number of parameters and deviance [difference in the −2log(likelihood) of the model of interest and the −2log(likelihood) of the saturated model] are also given.

Figure 2

Detection probability of Borrelia spp. as a function of the quantity of bacteria in tick races associated with three seabird species (A: puffin; B: kittiwake; C: murre). The dotted lines are standard errors (SE).

Modeled Borrelia prevalence

Our model-selection analysis suggested a difference in prevalence (corrected for detection) according to tick race (Table 1), as indicated by our top two models (cumulative QAICc weight = 0.69). However, the third ranked model provided some support for no difference among tick races (QAICc weight = 0.21). We therefore calculated tick race-specific prevalence estimates by model-averaging estimates over our entire model set (Table 2). This accounted for model-selection uncertainty (Burnham and Anderson 2002). Average prevalence across all tick races was 0.40 ± 0.14, with the highest prevalence in the puffin tick race and lowest in the kittiwake tick race (Table 2).

Table 2

Comparisons of naïve and corrected prevalence for three host races of Ixodes uriae (±1 SE)

Tick race	Naïve prevalence	Corrected prevalence (Ψ)	Bias
Puffin	0.53 (0.13)	0.57 (0.15)	0.04
Kittiwake	0.12 (0.07)	0.22 (0.14)	0.10
Murre	0.11 (0.08)	0.41 (0.32)	0.30
All races	0.22 (0.06)	0.40 (0.14)	0.18

The naïve prevalence was calculated from the raw data as the number of ticks observed with Borrelia spp. divided by the total number of ticks tested and the SE was based on a binomial distribution. The corrected prevalence takes into account the detection probability and is the model-averaged estimate across the entire model set, and thus accounts for both uncertainty in the estimation process, the detection process, and the model set. The corrected estimate for all tick races combined was from the highest ranking model in which prevalence was considered constant across tick races. Bias was calculated as the corrected prevalence minus the naïve prevalence. Naïve and corrected prevalence are based on the unknown tick sample only (see Materials and methods section).

Accounting for detection was important in this dataset as a naive overall estimate would have been 0.22 if all five PCR tests were considered (Table 2) and would have ranged from 0.10 to 0.16 if only a single-nested PCR test was considered. Race-specific prevalence estimates changed rank order of importance after accounting for detection (Table 2), further underlining the importance of accounting for nondetection in parameter estimates.

Discussion

Imperfect pathogen detection

Site-occupancy analyses showed that conventional PCR methods can strongly underestimate the presence of LB Borrelia in vector ticks and that this underestimation varies among host-specific tick groups. Several factors could contribute to the low performance, or sensitivity, of the diagnostic test (Banoo et al. 2008). For example, the intrinsic variability of the nested PCR approach, i.e. reproducibility, may lead to the imperfect detection of LB spirochetes (Altwegg 1995; Valkinas et al. 2006). However, in the present study, experimental conditions were carefully standardized to minimize this effect and such variability should not lead to race-associated differences in detection. More likely, the sensitivity of the detection technique to the relative concentrations of bacterial and vector DNA within the extract played a role (Cogswell et al. 1996; Banoo et al. 2008). Indeed, imperfect detection is often associated with poorly adapted pathogen strains or species as they tend to be present in low numbers within a given host (McKenzie et al. 2003). The qPCR approach we applied here showed that LB detection was related to the quantity of bacteria within an individual vector. The effect of infection intensity on amplification success was also supported by experimental tests using variable combinations of pure bacterial and host DNA; as host DNA concentrations increase relative to that of the bacteria, pathogen detection is considerably reduced (results not shown). However, despite the importance of infection intensity for detection, this factor could not entirely explain vector-specific differences. In particular, detection probability was the lowest in murre ticks, even though these tick vectors did not show the lowest average infection intensity. Detection probabilities could therefore also be linked to the presence of inhibitors that originate from the tick itself and that vary among races (Cogswell et al. 1996; Schwartz et al. 1997). Likewise, there may be differences among Borrelia strains found in the different races that alter detection probabilities. In the terrestrial cycle, previous evidence shows differences in host specificity and selectivity of Borrelia genospecies (Kurtenbach et al. 2002). In our system, a few atypical strains (non-B. garinii) have been previously isolated from I. uriae ticks (B. lusitaniae and B. burgdorferii-like isolates). These atypical strains tend to show lower detection probabilities compared with B. garinii, but the prevalence is too low to make inferences regarding patterns of infection among vector groups. In addition, no support for vector-related or host-related strain specificity in this system has yet been found (Duneau et al. 2008).

In most current studies, researchers employ PCR procedures (simple or nested) to determine the infective status of LB vector ticks (e.g. Ferquel et al. 2006; Halos et al. 2006; Kipp et al. 2006; Brisson et al. 2008). Despite the improvement that these techniques represent compared with immunofluoresence methods (Kahl et al. 1998), our results show that they can still have low sensitivity and result in a significant underestimate of pathogen presence. qPCR procedures are becoming more common and can improve detection sensitivity (as seen here). However, given the restricted conditions for the use of such methods, the targeted pathogen range can be large and the amplified sequence cannot always be used for isolate identification. Therefore, the method of choice will depend on the questions to be addressed and the system under consideration. Regardless of the chosen detection method, no method is likely to be perfect. In this way, the inclusion of site-occupancy modeling can be an efficient means to examine the accuracy of naïve estimates of pathogen prevalence and for determining potential sources of detection bias. We therefore advocate the inclusion of such models for all pathogen survey techniques.

Vector-specific patterns of infection

Model estimates indicated that, after taking into account vector-specific detection probabilities, Borrelia prevalence varied among the different host-associated tick races. Puffin ticks showed the highest prevalence and kittiwake ticks the lowest. As the current study was specifically designed to address the implications of cryptic vector structure on pathogen infection, we were not able to include spatial structure in our models. However, infected ticks of different host species came from all four geographic locations and preliminary analyses indicated no obvious effect of geography on infection intensity (results not shown). Future studies, with a more complete geographical coverage, will be designed to address spatial variation in detection and prevalence in this system.

Several alternative hypotheses could explain vector-specific patterns of infection. For example, differences in prevalence among the different tick groups could be due to differences in the life stage of the sampled ticks (Wang et al. 2003) or to their engorgement state (i.e. engorged or unengorged; Piesman et al. 1990). In the present study, all ticks were adults that were directly sampled from their specific hosts. Ixodes uriae host races are very faithful to their seabird host and therefore blood meals at the larval and nymphal stages have most likely been collected from the same host type (McCoy et al. 2001, 2005; Kempf et al. 2009a). Similarly, there were no statistical differences among seabird-tick races in their engorgement status (see Materials and methods section). Thus, inter-vectorial differences in infection are more likely associated with intrinsic differences in either the reservoir capacity of the different seabird host species or the relative competence of the different vector groups.

Successful transmission depends on the capacity of the vector to acquire, maintain, and transmit the pathogen (i.e. vector competence; Lane 1994; Edman 2002). Differences among tick races in Borrelia infection intensity and prevalence may therefore indicate differences in the ability of the bacteria to survive and replicate in the three vector groups. Similar patterns among different vector species have been found for a variety of zoonoses, including for example mosquito vectors of Plasmodium vivax (Adak et al. 1999), fleas transmitting Yersinia pestis (Lorange et al. 2005), or tick vectors of Anaplasma phagocytophilum (Teglas and Foley 2006). However, differences in infection among vector groups could also be linked to the seabird host species, rather than to the vector itself. Other studies support the idea that host susceptibility to infection and reservoir potential may be more important than variability in vector competence for the epidemiology of the pathogenic agent (Teglas and Foley 2006; Eisen et al. 2008). For example, under the ‘dilution effect’ model, the diversity and composition of the host community can have profound implications for disease epidemiology (LoGiudice et al. 2003). Based on this model, different host species show different reservoir capacities and vary in their level of infestation by the vector such that the presence of certain host species in the community limits the transmission potential of the pathogen and thus lowers disease risk. The differences that we find in Borrelia infection may thus be related to differences in the reservoir abilities of the various sympatric seabird host species (Ostfeld et al. 2006; Brunner et al. 2008). For instance, higher bacterial intensities in the tick could result from a higher infection dose received from the infecting host. Alternatively, host-associated blood factors could differentially favor the maintenance and replication of Lyme disease bacteria within the tick (Eisen et al. 2008). In particular, seabird-specific immune factors (i.e. bacteriolytic or inhibitory factors) or the seabird blood itself may affect bacterial fitness in the vector. Significant inter-host heterogeneity in infection can also occur within host populations due to individual differences in pathogen exposure and susceptibility (e.g. Staszewski et al. 2007) and is thus an important factor to take into account in sampling designs. Here, although adult ticks were sampled from numerous host individuals within each population, our dataset remains limited to fully explore this type of heterogeneity.

Vector-specific infection and disease dynamics

Regardless of whether vector or host-associated factors are responsible for the race-specific patterns we find in detection and infection, cryptic vector population structure can have important epidemiological implications. In particular, certain vector–host combinations are going to be more important than others in pathogen transmission. For example, among the three host–vector combinations considered here, our results indicate that a particular effort should be made to survey the presence and structure of LB spirochetes associated with Atlantic puffins. The potential role of this species in Lyme disease epidemiology is further amplified by the fact that puffins seem to move at larger spatial scales than either kittiwakes or common murres, and therefore may disperse vectors and disease agents more widely (McCoy et al. 2003). Pathogen strains carried by this host–vector combination should therefore be distributed at larger spatial scales, should be less impacted by the effects of local drift, and may be exchanged more readily with strains from the terrestrial cycle. Large-scale dispersal may also favor pathogen emergence by generating novel diversity via the horizontal transfer of plasmid-borne genes (Qiu et al. 2004). These factors, in combination with other important variables that may affect pathogen transmission rates such as host community composition, host density, and vector basic reproductive number (R0) (LoGiudice et al. 2003, 2008; Kurtenbach et al. 2006; Ostfeld et al. 2006; Hartemink et al. 2008), will likely have a significant impact on LB dynamics over space and time and should therefore be considered more generally in epidemiological models of this disease.

Concluding remarks

The differences we find in pathogen infection and detection among different vector groups are relevant to many vector-borne zoonoses and should be considered more explicitly in other systems of medical and economic importance. Such differences can lead to vector-specific patterns of transmission which, in turn, can affect the ecological and evolutionary dynamics of the system, e.g. selection on the pathogen, and ultimately, the epidemiological patterns of the associated disease. Likewise, prevalence underestimation and biases in detection may be common in many parasitological and epidemiological studies. Efforts to standardize pathogen detection methods are difficult under certain sampling conditions (e.g. Banoo et al. 2008), and thus estimating detection directly is advocated. The site-occupancy, and related, modeling methods continue to be developed and can be applied in many novel ways to host–parasite studies and disease epidemiology (Mackenzie et al. 2006, 2009; Kendall 2008; Conn and Cooch 2009). Such analyses can provide more robust parameter estimates and can reveal new avenues to explore for our understanding of host–pathogen systems.