Experimental Evidence for Opposing Effects of High Deer Density on Tick-Borne Disease Prevalence and Hazard

doi:10.21203/rs.3.rs-190503/v1

Download PDF

Research Article

Experimental Evidence for Opposing Effects of High Deer Density on Tick-Borne Disease Prevalence and Hazard

https://doi.org/10.21203/rs.3.rs-190503/v1

This work is licensed under a CC BY 4.0 License

Version 1

posted

You are reading this latest preprint version

Identifying the mechanisms driving disease risk is challenging for multi-host pathogens, such as Borrelia. burgdorferi s.l., the tick-borne bacteria causing Lyme disease. Deer are tick reproduction hosts but do not transmit B. burgdorferi s.l., resulting in potentially opposing effects on transmission. Here, we use a deer exclosure experiment to test three hypotheses for how high deer density shapes B. burgdorferi s.l. prevalence in ticks: (H1) high transmission on rodents due to higher tick densities; alternatively, (H2) low B. burgdorferi s.l. prevalence because more ticks feed on deer rather than transmission-competent rodents (dilution effect); (H3) ecological cascades, whereby lower vegetation decreases rodent abundance thus reducing transmission. Although we found support for all three mechanisms, prevalence was reduced almost 3-fold in high deer density plots compared to exclosures, suggesting that the dilution (H2) and cascade (H3) mechanisms outweighed the increased opportunities for transmission (H1). High deer density led to lower vegetation and fewer rodents, providing evidence for an ecological cascade. However, Lyme disease hazard (density of infected ticks) was increased 5-fold at high deer densities due to an 18-fold rise in tick density. This demonstrates that reproduction hosts like deer can drive up vector-borne disease hazard at high densities, despite simultaneously reducing pathogen prevalence.

Agroecology

Terrestrial Ecology

Animal Science

Borrelia burgdorferi s.l.

dilution effect

Lyme disease

ecological cascades

Ixodes ricinus

Ecological trophic cascades occur when a change in population or activity of a keystone species precipitates a cascade of alterations through trophic levels in an ecosystem (Paine 1980). For example, in coastal kelp forests in the Northeast Pacific, sea otters (Enhydra lutris) prey on herbivorous sea urchins (Strongylocentrotus spp.), which releases grazing pressure on kelp and has further consequences for benthic plant and animal communities (Carter et al. 2007; Watson and Estes 2011).

While there are many examples for trophic cascades in ecology, tests of their extension to infectious disease risk have received less attention. This concept could be particularly significant for vector-borne diseases, which rely upon vectors that carry and transmit pathogens from one host to another, as trophic effects could manifest from changes in both vector and reservoir host populations. For instance, one study showed that the introduction of Burmese pythons (Python bivittatus) resulted in a decline in some mammal species but not in those acting as transmission hosts for the Everglades virus. This affected the feeding behaviour of mosquitoes and increased blood meals taken from transmission hosts (Hoyer et al. 2017). Other researchers found a positive association between coyote (Canis latrans) density and Lyme disease incidence because coyotes have negative effects on red fox (Vulpes vulpes) populations, which caused an increase in small mammals that transmit Lyme disease pathogens (Levi et al. 2012). While these studies show how a change in predators can affect vector-borne pathogen risk through predation of transmission hosts, no study to our knowledge has tested for trophic cascades effects on pathogens arising from changes in large herbivore populations that alter habitat (e.g. vegetation structure) for transmission hosts. It is already established that large herbivores alter habitat (Buesching et al. 2011), with cascading effects on small mammal communities (Flowerdew and Ellwood 2001; van Wieren and Bakker 2008), which are transmission hosts for many pathogens, including tick-borne pathogens (e.g. Tick-borne encephalitis virus and the bacteria causing Lyme disease). This potential mechanism of vector-borne disease risk seems likely and requires testing since it has not been studied before.

Lyme disease, the most prevalent vector-borne disease in humans in the northern hemisphere, is an emerging disease in Europe and North America (Stanek et al. 2012). It is caused by the Borrelia burgdorferi sensu lato complex of bacteria (hereafter referred to as B. burgdorferi s.l.) and transmitted by Ixodid ticks. In Europe, the primary vector is Ixodes ricinus (Stanek et al. 2012), a generalist three-host tick that feed on a wide range of vertebrate species (Hofmeester et al. 2016). There are several genospecies of the bacteria within the B. burgdorferi s.l. complex, each with different transmission host associations: in Northern Europe, rodents transmit B. afzelii (Hanincová et al. 2003) while birds are associated with B. valaisiana and B. garinii (Heylen et al. 2014). B. burgdorferi sensu stricto, the fourth most abundant genospecies, can be transmitted by both mammals and birds. Transmission is therefore predicted to increase as more immature ticks feed on infected hosts such as rodents (Daniels and Fish 1995; Perkins et al. 2006).

Deer, on the other hand, do not transmit B. burgdorferi s.l. (Bolzoni et al. 2012; Pacilly et al. 2014). However, deer are the most important hosts for feeding adult female I. ricinus ticks prior to egg laying and, as such, high deer densities are often associated with high tick population densities (Gray et al. 1992; Daniels et al. 1993; Gilbert et al. 2012; Pacilly et al. 2014; Mysterud et al. 2016). We can therefore predict that high deer densities may increase tick burdens on transmission hosts such as rodents, thereby increasing pathogen transmission (Daniels and Fish 1995). Deer can also feed immature tick life stages (Hofmeester et al. 2016) so we could predict that, at high deer densities, a higher proportion of immature ticks feed on deer instead of transmission hosts, thus lowering B. burgdorferi s.l. prevalence in the tick population through a dilution effect (Norman et al. 1999; Ostfeld and Keesing 2000). An alternative mechanism could also result in lower NIP from high deer densities since deer grazing can lead to sparser vegetation cover (Buesching et al. 2011), with cascading negative effects on the density of transmission hosts such as rodents (Flowerdew and Ellwood 2001; van Wieren and Bakker 2008) and birds (Allombert et al. 2005). While some studies have demonstrated an association between deer abundance and Lyme disease hazard (Vourc’h et al. 2016; Takumi et al. 2019), the mechanisms driving this association are equivocal, and a particular gap in knowledge is the potential role of deer in regulating NIP and Lyme disease hazard through their cascading effects on vegetation and therefore transmission hosts such as rodents. Therefore, a key novelty of this study is to test for each component of the ecological cascade and to understand how deer at high density might affect NIP. Such cascading effects on vegetation and rodents might be predicted to occur even at small spatial scales as plant communities, herbivore grazing pressure, rodent activity and tick distribution are all spatially highly heterogeneous (Saïd and Servanty 2005; Jones et al. 2011; James et al. 2014; Vourc’h et al. 2016). Consequently, rodents are expected to spend more time in patches with favourable vegetation structure and less grazing pressure. These heterogeneities provide the opportunity to test cascading effects experimentally, at small spatial scales most amenable to manipulation. Here, we use a replicated deer exclosure experiment to test three hypotheses for the effects of high deer density on NIP:

H1-Increase in transmission potential

high deer density will produce a high density of questing nymphs (DON) and thus a higher tick burden on transmission hosts (Daniels and Fish 1995) and higher NIP in questing ticks than areas without deer.

H2-Dilution effect

lower prevalence in areas of high deer density due to a higher proportion of larval ticks feeding on deer instead of transmission hosts such as rodents.

H3-Ecological cascade

intense grazing pressure from high deer densities will result in shorter vegetation and, therefore, fewer rodents (Flowerdew and Ellwood 2001; van Wieren and Bakker 2008; Buesching et al. 2011), and lower prevalence. Rodent activity might also be affected by high deer density directly through disturbance.

These different mechanisms are not mutually exclusive and the effect of high deer density on NIP will depend on the relative strength of the first (transmission potential) hypothesis compared to the other two (dilution and ecological cascade hypotheses) (Fig. 1). A further key aim was to test the impact of high deer density on Lyme disease hazard, which is defined as the density of infected nymphs (DIN) in the environment and is the product of NIP and DON. It is difficult to predict the effect of deer on Lyme disease hazard because it will depend on the relative strengths of the effect of deer on DON and the mechanisms driving NIP (Fig. 1).

The use of experimental deer exclosures is particularly suitable for testing the ecological cascades hypothesis as it maximises habitat impacts while avoiding noise being introduced from heterogeneities in land use, habitat, topography and climate that typify landscape-scale surveys. In addition, the exclusion of deer has high applied relevance since deer fencing is a common land management tool and is increasingly being used in mitigating the impacts of ticks (Gilbert et al. 2012).

2.1. Experimental design

The experiment took place at Glensaugh research farm in Aberdeenshire, Scotland (56.914217N, -5.532070E), in a 15 ha enclosure of upland moorland where five red deer (Cervus elaphus) stags were kept (32.5. km^− 2). Within this moorland enclosure, four replicates of fenced-unfenced pairs of plots of 0.23 ha were set up in 2004/2005. The four fenced plots (hereafter referred to as deer exclusion plots) excluded deer while the four unfenced plots (hereafter referred to as high deer density plots) were accessible to the deer and subjected to high grazing pressure. As part of a different study, each plot was divided into five subareas each of approximately 15mx15m to create five habitats: high density birch, low density birch, a single birch in the centre of the plot, high density pine, and a control which was not planted and consisted of Calluna vulgaris-dominated heathland. Tree saplings were 9 years old by the time this study started in 2013.

2.2. Quantifying questing I. ricinus nymph density (DON)

Questing nymphs were surveyed in 2013, 2014, 2018 and 2019 using a standard dragging method (Falco and Fish 1992) which consists of dragging a 1m x 1m square of woollen blanket material over the ground vegetation for 10m linear transects. In 2013 and 2014, 10 transects were surveyed in each of the heather, high density birch and pine habitats for all plots (n = 720 transects/year) while in 2018 and 2019, six transects were surveyed in each of the same three habitats (n = 432 transects/year). Tick surveys were done in May, July and September between 0900 hours and 1900 hours. Air temperature and relative humidity were recorded for each transect to be taken account of in statistical models as these factors affect tick activity (Gilbert 2010). Ground vegetation height was recorded at the beginning (0m), middle (5m) and end (10m) of each transects using a sward stick. Questing nymphs were counted, collected and stored at -20°C for later pathogen analysis.

2.3. Quantifying rodent activity and tick burdens

To test the effect of vegetation on rodents and the effect of rodents on NIP (H2- dilution hypothesis and H3 – ecological cascade hypothesis) as well as to quantify tick burdens on rodents (H1- increase in transmission hypothesis), rodent activity was estimated in July 2017 and 2018 using a live-trapping method. In 2017, four non-selective Sherman live traps (16 x 5 x 6.5 cm. Sherman Inc., Tallahassee, Florida) were baited with oats and placed in each of the five habitats within each plot with 10m spacing, resulting in a total of 20 traps per plot (i.e. 200 trap nights (TN) per deer treatment). Traps were set for three consecutive nights in two of the paired fenced/unfenced plots (replicates 1 and 2) and for two consecutive nights in the other two paired fenced/unfenced plots (replicates 3 and 4). In 2018, four traps were installed for four consecutive nights in the three habitats in which ticks were surveyed (heather, high density birch and pine habitats) in each plot, resulting in 12 traps per plot (i.e. 192 trap nights (TN) per high deer treatment). Traps were activated after 1600 hours and checked every morning before 1000 hours. For all captures, species, sex, weight and approximate age (juveniles or adult) of each individual were recorded. Ticks attached to rodents were counted and collected from around the head and ears and stored in 100% ethanol. All captured rodents were released at the capture site. For analysis, we took the proportion of traps with a capture as our rodent activity index. Both tick burden and rodent captures data were needed to test the dilution effect (H2) which requires a comparison between high deer density plots and deer exclosures in the relative proportion of the larval tick population that feeds on rodents versus deer.

2.4. Measurement of B. burgdorferi s.l. prevalence (NIP)

As part of testing all three hypotheses of the effects of high deer density on NIP, and to estimate Lyme disease hazard (DIN), DNA was individually extracted from questing nymphs using an ammonia extraction (Gern et al. 2010). B. burgdorferi s.l. was detected using two methods for this study: nested PCR and qPCR. Ticks collected in 2013 and 2014 were tested using a nested PCR targeting the 5S-23S intergenic spacer region using the protocol described by (Rijpkema et al. 1995). Following an issue arising from the nested PCR protocol, samples collected in 2018 and 2019 were tested using a qPCR method. A qPCR protocol on fragments of OspA genes (Heylen et al. 2013) was optimized using the IQ™ Supermix (Bio-Rad Laboratories, Hercules, USA) in a Stratagene Mx3005P thermal cycler (Agilent, Santa Clara, US). Each reaction contained IQ™ Supermix, two primers at 200nM (B-OspA_modF: AATATTTATTGGGAATAGGTCTAA and B-OspA_borAS:-CTTTGTCTTTTTCTTTRCTTACAAG), the probe (B-OspA_mod-probe:-FAM-AAGCAAAATGTTAGCAGCCTTGA-BHQ-1™) at 100nM and 3µL of DNA. One positive and one negative control were added for every plate.

We confirmed the correspondence between the two PCR protocols by using 61 known positive samples (by nested PCR, mix of genospecies) and 344 known negative samples. In all cases, results from the qPCR matched those of the nested PCR. Positive samples from the qPCR protocol were subjected to the nested PCR protocol to identify B. burgdorferi s.l. genospecies. All PCR products from positive samples were Sanger sequenced to identify the genospecies of B. burgdorferi s.l. For analysis, we used the proportion of questing nymphs infected (NIP).

2.1.1. Statistical analyses

All statistical analyses were performed in R version 3.5.1 (R-core-team, 2013). For Generalized Linear Mixed Effects Models (GLMMs) and Generalized Linear Models (GLMs), we tested for potential collinearity between explanatory variables by calculating variance inflation factors (VIFs) and variables for which the VIF was above 4 were discarded from the model (Zuur et al. 2009). Model selection was done using the dredge function from the MuMIn package (Barton 2019) based on the corrected Akaike Information Criterion (AICc) (Brewer et al. 2016). When appropriate, we conducted post hoc Tukey tests to assess pairwise comparisons between levels of categorical variables.

The effects of high deer density on nymphal infection prevalence (NIP) withB. burgdorferis.l.

Our three hypotheses are all concerned with potential mechanisms through which high deer density might affect NIP. To test for this effect of high deer density, we used a binomial GLMM with a logit link and NIP for B. burgdorferi s.l. (3 prevalence estimates per habitat and per year) as the response variable. The full model included deer treatment (deer exclusion or high deer density) as our main predictor, as well as month (May, July or September), habitat (high density pine, high density birch or heather control) and year (2013, 2014, 2018 or 2019). Random effects of a combined plot and habitat parameter (e.g. high density birch in fenced replicate 1) and an observation level random effect were included, the latter to account for overdispersion (Harrison 2015).

Hypothesis 1: increase in transmission potential

To test the hypothesis that high deer density results in high questing nymph density, which could result in higher tick burdens on rodents, we investigated the effects of high deer density on DON. We used a hurdle generalized linear mixed effect model (GLMM) with a Poisson distribution. We tested for zero-inflation using the DHARMa package (Hartig 2020) and we chose a hurdle model as the number of zeros were underestimated with a Poisson GLMM (Lewis et al. 2011). The response variable was the number of questing nymphs collected per 10m transect and the full model included deer treatment, month, year, habitat, ground vegetation height, relative humidity, temperature and whether the ground was dry during collection. We also included the interactions between deer treatment and month and between deer treatment and habitat, as the effects of deer on ticks might vary between months and habitat type. A random effect for each plot habitat combination (e.g. high density birch in fenced replicate 1) and an observation level random effect were included to account for overdispersion (Elston et al. 2001; Harrison 2014).

Although we initially planned to compare tick burdens on rodents between deer density treatments, only two rodents were caught in high deer density plots (see H3) so this analysis could not be performed.

Hypothesis 2: dilution effect

The dilution effect hypothesis implies that the relative proportion of the larval tick population feeding on deer vs rodents (i.e. tick burden x relative abundance of each host type) is higher in high deer density plots compared to exclosures. Testing this formally was not possible because only two rodents were captured in high deer density plots (as predicted by H3), which precluded us from fitting a model of the effects of high deer density on tick burden on rodents. We did not count tick burdens on deer. To evaluate this hypothesis, we therefore report descriptive statistics for tick burdens on rodents and relative host abundance.

Hypothesis 3: ecological cascade

To test the first part of the cascade, that high deer density negatively affects ground vegetation height (H3), we used a GLMM with a Gaussian distribution and ground vegetation height as our response variable. The full model included deer treatment, month and year as fixed effects and a random effect of each plot habitat combination.

To test the second part of the cascade, that lower vegetation negatively affects rodent activity, we conducted two analyses, both using a binomial GLMM with a logit-link and our rodent index as the response variable. For the first, the full model included vegetation height, habitat and year as fixed covariates and a random effect of each plot habitat combination. Vegetation height was averaged per habitat type and per year and we used vegetation height from the heather habitat for the single birch and sparse birch habitats, because there was no difference between these three (Gilbert et al. 2012). For the second analysis, testing for a direct impact of deer on rodent activity, the full model included deer treatment, habitat and year as fixed covariates and a random effect of a random effect of each plot habitat combination.

To test the last part of the cascade, that rodent abundance positively affects NIP, we used a binomial GLMM. The response variable was NIP for B. afzelii and the full model included the number of rodents captured the previous year (data for 2017 and 2018), habitat, year and month. We added the combined plot and habitat parameter variable as a random effect and an observation level random effect was included to account for overdispersion (Harrison 2015).

The effects of high deer density on Lyme disease hazard

To test how high deer density affects Lyme disease hazard (DIN), we used a hurdle GLMM with a Poisson distribution, as the number of zeros was underestimated with a Poisson GLMM (Lewis et al. 2011). The response variable was the density of infected nymphs and, as this variable was averaged at the plot level (3 estimates per habitat and per year), we used an offset for the area to transform the values to integers (Zuur et al. 2009). The full model included deer treatment, month, year, habitat, temperature and whether the ground was dry during tick survey as fixed variables and the combined plot and habitat level variable and an observation level random effect were included as random terms.

Over four years of data collection, 12,310 questing nymphs were sampled (Table 1). A random subset of 1,000 nymphs was examined under the microscope for species identification using specific keys (Hillyard 1996) and all were found to be Ixodes ricinus. Thus, we assumed that all questing ticks collected in this study were I. ricinus. We tested all 585 questing nymphs collected from deer exclusion plots and a subset of 1,042 from high deer density plots for B. burgdorferi s.l. Prevalence (NIP) was 0.9% [95%CI: 0.4–1.6] in high deer density plots and 2.7% [95%CI: 1.6–4.4] in deer exclusion plots (Table 1). The dominant B. burgdorferi s.l. genospecies was the rodent-associated B. afzelii (92%, n = 23/25) followed by the bird-associated B. valaisiana (4%, n = 1/25) and B. garinii (4%, n = 1/23). Over two years of trapping, 24 individual bank voles were caught in deer exclusion plots (6.69/100TN, SD = 2.41), while only two were captured in high deer density plots (0.51/100TN, SD = 0.97) and no other rodent species were detected (Table 1).

Table 1

Mean values across all plots and habitat-sub areas, separated by year and deer treatment (high deer density and deer exclusion plots), for the density of questing nymphs (DON), number of bank voles caught per 100 trap nights, prevalence of *B. burgdorferi* s.l. in questing nymphs (NIP) and density of nymphs infected with *B. burgdorferi* s.l. (DIN). Ticks were not collected in 2017. Rodent trapping was conducted in 2017 and 2018 only. SD - standard deviation; 95% CI − 95% confidence interval.
Treatment	Year	DON. 10m^− 2 ± SD	Bank vole abundance /100TN ± SD	NIP (%) [95% CI]	DIN. 1000m^− 2 ± SD
High deer density	2013	9.78 ± 9.99		0.56 [0-3.1]	4.81 ± 14.42
	2014	9.67 ± 9.27		0	0
	2017		0.50 ± 1.12
	2018	25.57 ± 25.84	0.52 ± 0.90	0.24 [0-1.3]	3.21 ± 9.64
	2019	6.21 ± 7.60		2.47 [0.9–4.4]	6.12 ± 9.74
	Average	12.10 ± 15.46	0.51 ± 0.97	0.86 [0.4–1.6]	3.86 ± 10.17
Deer exclusion	2013	0.37 ± 0.84		1.92 [0.2–6.8]	0.96 ± 1.50
	2014	0.54 ± 1.01		4.23 [1.9–7.9]	2.09 ± 3.20
	2017		6.95 ± 2.51
	2018	0.67 ± 1.35	6.25 ± 2.70	0.72 [0-3.9]	0.46 ± 1.39
	2019	0.65 ± 1.11		3.10 [0.9–7.7]	1.96 ± 3.20
	Average	0.53 ± 1.06	6.69 ± 2.41	2.74 [1.6–4.4]	1.41 ± 2.54

The effects of high deer density on nymphal infection prevalence (NIP) with B. burgdorferi s.l.

While the best fit model for NIP only included habitat type, year and month while the second-best model (0.44 increase of AICc compared to the best model) included deer treatment, habitat type and year. Based on this second model, predicted NIP was lower in high deer density plots (1.0%, 95%CI: 0.3–4.2) compared to deer exclusion plots (2.2%, 95%CI: 0.6–7.6) (Fig. 2). NIP also varied across years, habitats and months and the full model results are provided in the Supplementary Information 1.

Hypothesis 1: increase in transmission potential

Effects of high deer density on DON: The final model included deer treatment, month, habitat type, year, relative humidity, temperature and whether the ground was wet. On average, predicted DON was 18 times higher in high deer density plots (DON = 10.8, 95%CI: 0.1–36.2) compared to deer exclusion plots (DON = 0.6, 95%CI: 0-3.5) (Fig. 3). DON varied across months, years and habitats and was influenced by temperature, relative humidity and whether the ground was wet (SI 1).

Effects of high deer density on tick burdens on rodents: Bank voles in deer exclusion plots (n = 24) harboured, on average, 14.2 larvae (SD = 11.66) while bank voles in high deer density plots (n = 2) had 17.5 larvae (SD = 10.61). Low captures in high deer density plots prevented us from statistically comparing tick burdens between treatments further. No nymphs were found attached to a bank vole.

Hypothesis 2: Dilution effect

Comparison between high deer density plots and exclosures on the relative proportion of the larval population feeding on rodents vs deer: In deer exclusion plots, bank voles fed 340.8 larvae (24 individuals feeding 14.2 larvae each on average). In high deer density plots, bank voles fed 35 larvae (2 individuals feeding 17.5 larvae each on average). These data suggest that 9.7 times more larvae fed on bank voles in deer exclusion plots compared to high deer density plots. In deer exclusion plots, no larvae fed on deer (there were no deer in exclosures). High deer density plots contained a deer density of 32.5. km^− 2 but we were not able to measure tick burdens on these deer.

Hypothesis 3: Ecological cascade

Effects of high deer density on ground vegetation: The final model included deer treatment, month, habitat type and year. Ground vegetation was predicted to be 14.4 cm lower in high deer density plots (37.0 cm, 95%CI: 31.6–42.4) compared to deer exclusion plots (51.4 cm, 95%CI: 46.0-56.9) (Fig. 4A). Ground vegetation height was also influenced by month and habitats (SI 1).

Effects of vegetation height on rodents: The final model included deer treatment and year and the predicted number of bank voles captured was positively correlated with ground vegetation height and increased by 1.17 captures per 100 TN for every 1 cm increase in vegetation height (Fig. 4B). Bank vole capture rate also varied across both years (SI 1).

Effects of high deer density on rodent activity: The final model included deer treatment only and the model predicted 13 times more bank voles to be captured in deer exclusion plots (6.6/100TN, 95%CI: 4.5–9.5) compared to high deer density plots (0.5/100TN, 95%CI: 0.1-2.0) (Fig. 4C).

Effects of rodent capture rate on NIP for B. afzelii: While the best model only included year and month, rodents were included as predictors in the second-best model (ΔAICc = 1.9) and predicted NIP increased by 0.01% for every vole caught per 100 TN (Fig. 4D) (SI 1).

The effects of high deer density on Lyme disease hazard (DIN)

The final model included deer treatment, habitat type and temperature. Predicted Lyme disease hazard (DIN) was 5.1 times higher in high deer density plots (5.1/1000m⁻², 95%CI: 0-30.8) compared to deer exclusion plots (1.0/1000m⁻², 95%CI: 0-6.7) (Fig. 5). DIN was also influenced by habitat and temperature (Supplementary materials 1) and the same results were obtained when examining DIN for B. afzelii alone (SI 2).

A replicated paired fencing experiment allowed us to test three hypotheses (increased transmission potential, dilution and ecological cascades) for how high deer density may affect B. burgdorferi s.l. infection prevalence in questing ticks (NIP). Furthermore, we could test how these effects on NIP interact with the role of deer in supporting tick populations and shaping Lyme disease hazard. We found that nymphal infection prevalence was lower in high deer density plots along with evidence pointing to a combination of dilution and trophic cascading effects. Despite the lower NIP, Lyme disease hazard (DIN) was five times higher in high deer density plots compared to deer exclusion plots, due to a strong, positive association between deer and nymph density.

Our transmission potential hypothesis (H1) predicted that a high deer density would result in an increased tick density, causing higher tick burdens (and therefore increased transmission potential) on individual rodents. We did find a strong positive effect of high deer density on DON, highlighting the importance of deer in driving tick populations, as found in previous studies (Gray et al. 1992; Daniels et al. 1993; Gilbert et al. 2012; Pacilly et al. 2014; Mysterud et al. 2016). All else being equal, this would have predicted increased opportunities for pathogen transmission. In contrast however, we found high deer density plots had half the predicted NIP of exclusion plots, suggesting that increased tick densities and transmission were counter-acted by other factors driving NIP down.

Several studies have reported positive correlations between DON and tick burdens on rodents (Daniels and Fish 1995; Hofmeester et al. 2017a) and other hosts (Gilbert et al. 2017). This may also have been the case for rodents in our experiment but, as only two bank voles were captured in high deer density plots (as predicted by the ecological cascade hypothesis), we were not able to test for differences between deer treatments in individual rodent tick burdens. This dramatic reduction in rodent captures in high deer density plots is a key finding and invalidates H1’s assumption of similar rodent densities with deer density. Therefore, even if individual rodents had higher tick burdens in high deer density areas, rodent density was probably too low for maintaining effective pathogen transmission in this environment.

Our observation that B. burgdorferi s.l. prevalence was almost three times higher in deer exclusion plots, despite estimates being associated with considerable uncertainty, is consistent with predictions from both the dilution effect (H2) and ecological cascade (H3) hypotheses. The dilution effect predicts lower NIP at high deer densities due to a lower proportion of the larval tick population feeding on rodents, which are transmission hosts for B. burgdorferi s.l., than on deer that do not transmit the pathogens (Gray et al. 1992; Ostfeld and Keesing 2000; Vourc’h et al. 2016). This is challenging to test because, ideally, it requires density estimates and tick burdens of the key hosts. We could not obtain tick burden measurements on deer and we caught only two voles in high deer density plots, precluding statistical test. However, based on differences in vole abundance between deer treatments, we estimated that voles fed almost 10 times more larvae in deer exclusion plots compared to high deer density plots. This compares to no larvae feeding on deer in exclusion plots (as deer were absent), in contrast to the situation in high deer density plots where it is highly likely that most larvae fed on deer. This is due to the (i) very high deer density (32.5. km^− 2) and (ii) high densities of questing nymphs (a result of high larval survival) in high deer density plots. We therefore suggest that a dilution effect was one of the mechanisms operating on NIP in response to high deer densities, as it has been previously suggested for B. burgdorferi s.l. (Gray et al. 1992; Ostfeld and Keesing 2000; Vourc’h et al. 2016) and predicted for other tick-borne pathogens that deer do not transmit, including tick-borne encephalitis virus (Bolzoni et al. 2012) and Louping ill virus (Gilbert et al. 2001).

While we did not survey other known hosts for larval ticks such as birds or shrews (Sorex spp.) (Klaus et al. 2016; Hofmeester et al. 2017b), previous studies suggest that the effect of high deer density on these groups is likely to be negative as well (Flowerdew and Ellwood 2001; Allombert et al. 2005; Herder et al. 2016). It is therefore likely that the populations of alternative hosts were low in high deer density plots and that the majority of larvae must have been feeding on red deer. In addition, the fact that almost all nymphs that tested positive (92%) were infected with the rodent associated pathogen, B. afzelii, suggests that other hosts (e.g. birds) did not contribute much to B. burgdorferi s.l. transmission in our system.

To the best of our knowledge, this is the first test of the ecological trophic cascade hypothesis of B. burgdorferi s.l. prevalence, predicting that grazing pressure from high deer densities will reduce the vegetation, resulting in fewer rodents and therefore lower NIP. We found support for most of the expected trophic links: high deer density resulted in shorter vegetation, highlighting the effects of deer grazing on vegetation structure (Flowerdew and Ellwood 2001; Buesching et al. 2011). Lower ground vegetation was associated with fewer bank voles in support of denser and higher ground vegetation providing better food, shelter and protection from predators (Flowerdew and Ellwood 2001; Eccard et al. 2008). We captured 13 times fewer bank voles in high deer density plots, consistent with our predictions and with previous work that showed higher densities of bank voles (Buesching et al. 2011) and wood mice (Buesching et al. 2011; Smit et al. 2011) in plots excluding deer. It is possible that direct disturbance from deer (which we did not quantify) could also be contributing in addition to the effect of reduced vegetation cover (Flowerdew and Ellwood 2001). However, we could not confirm the last link for the trophic cascade hypothesis (i.e. a strong link between rodent activity and NIP), as bank vole abundance the previous year was only a weak positive predictor of NIP for B. afzelii (Fig. 4D). The weak association between vole abundance and NIP in this experiment could be due to several factors, including low prevalence of B. burgdorferi s.l. overall, providing insufficient signal for detecting unequivocal statistical effects. Another contributing factor might have been the small spatial scale of our experiment plots (0.23 ha), facilitating likely movements of rodents between fenced and unfenced plots, which were 30–100 m apart.

The small spatial scale of our plots, while a potential issue for confirming a link between rodent abundance and NIP, proved sufficient to confirm strong effects of high deer density with DON, NIP, Lyme disease hazard, vegetation height and rodent activity. While a previous meta-analysis of deer exclusion effects on tick abundance (Perkins et al. 2006) suggested that deer exclusion areas of at least 2.5 ha may be necessary to have an effect, our results show smaller plots sizes to be sufficient for testing impacts of deer on ticks, consistent with other more recent studies (Gilbert et al. 2012 - plots of 0.2–0.25 ha; Mysterud et al. 2016 - plots of 0.04 ha)

Irrespective of the mechanisms driving NIP, the most critical parameter governing public health and policy importance for Lyme disease is the density of infected nymphs (DIN = NIP x DON), which is the key proxy for Lyme disease hazard in the environment. DIN was five times higher in high deer density plots compared to exclosures due to deer having a strong positive effect on DON. Similarly, other studies have found a positive correlation between deer density and Lyme disease hazard (Vourc’h et al. 2016; Takumi et al. 2019) or Lyme disease incidence in humans (Mysterud et al. 2016). Using an experimental system with high deer densities (32.5. km^− 2), we were able to demonstrate that this role of deer as tick reproduction hosts is, at least at high deer densities, more important than their role in lowering NIP through dilution and trophic cascade effects in shaping Lyme disease hazard.

In contrast to our experiment, previous studies have shown Lyme disease hazard to be associated with higher densities of transmission hosts (van Duijvendijk 2016; Takumi et al. 2019). However, such an association requires enough vectors in the environment to transmit the pathogen effectively (Logiudice et al. 2008), whereas in our experimental situation, the plots with high densities of transmission hosts did not have deer, and therefore had few ticks to aid transmission. Thus, we might expect Lyme disease hazard to be highest in an environment supporting both high numbers of transmission hosts and tick reproduction hosts. However, based on our findings and previous research showing that high deer densities have strong negative effects on rodent abundance (Flowerdew and Ellwood 2001; van Wieren and Bakker 2008), such a combination may not commonly occur in nature. This demonstrates the importance of taking a systems approach, including considering potential ecological cascades when predicting disease risk. While our experiment used extreme contrasts in deer densities (zero vs 32.5. km^− 2) it is notable that we were able to detect effects on infection dynamics even on very small spatial scales (0.23 ha plots).

Our results highlight the need for taking a systems approach to such a complex disease system where different host types might affect each other’s densities through habitat modification, or other means such as predation (Levi et al. 2012). Gaps in knowledge that can now be addressed include testing for cascading effects from deer on B. burgdorferi s.l. prevalence at landscape scales, and including a full range of intermediate deer densities, to investigate non-linearities and thresholds of effects of deer and to test whether the mechanisms supported here operate in a non-experimental setting.

In conclusion, we found that high deer density could lower Lyme disease pathogen prevalence by a combination of dilution and trophic cascading effects. Despite this, Lyme disease hazard was five times higher in high deer density plots due to a strong positive effect of deer on tick density. This study is, to our knowledge, the first to test for cascading effects of deer on B. burgdorferi s.l. prevalence via grazing and suppressing rodent abundance. This study highlights the need for a systems approach to understand disease dynamics and risks that could arise from complex ecological interactions between host types and habitat.

Acknowledgements

We would like to thank Donald Barrie, manager at Glensaugh research farm for letting us conduct our field work there and Alison Hester and Richard Hewison for setting up the original fenced plots.

Funding

SG was supported by a PhD studentship joint funded by the University of Glasgow and The Macaulay Development Fund. LG was supported by the Scottish Government’s Rural and Environmental Science and Analytical Services Division.

Conflicts of interest/Competing interests

The authors declare that there is no conflict of interest.

Ethics approval

Rodent trapping was conducted under the Home Office Regulations Project license 70/8543 and under the personal license I2B15B1E9.

Consent to participate

Not applicable.

Consent for publication

Not applicable.

Availability of data and material

The datasets generated during and/or analysed during the current study are available in the [NAME] repository, [PERSISTENT WEB LINK TO DATASETS]

Allombert S, Gaston AJ, Martin JL (2005) A natural experiment on the impact of overabundant deer on songbird populations. Biol Conserv 126:1–13. https://doi.org/10.1016/j.biocon.2005.04.001
Barton K (2019) MuMIn: multi-model infectence. In: R package version 1.43.6 http://mumin.r-forge.r-project.org/
Bolzoni L, Rosà R, Cagnacci F, Rizzoli A (2012) Effect of deer density on tick infestation of rodents and the hazard of tick-borne encephalitis. II: Population and infection models. Int J Parasitol 42:373–381. https://doi.org/10.1016/j.ijpara.2012.02.006
Brewer MJ, Butler A, Cooksley SL (2016) The relative performance of AIC, AICC and BIC in the presence of unobserved heterogeneity. Methods Ecol Evol 7:679–692. https://doi.org/10.1111/2041-210X.12541
Buesching CD, Newman C, Jones JT, Macdonald DW (2011) Testing the effects of deer grazing on two woodland rodents, bank voles and wood mice. Basic Appl Ecol 12:207–214. https://doi.org/10.1016/j.baae.2011.02.007
Carter SK, VanBlaricom GR, Allen BL (2007) Testing the generality of the trophic cascade paradigm for sea otters: A case study with kelp forests in northern Washington, USA. Hydrobiologia 579:233–249. https://doi.org/10.1007/s10750-006-0403-x
Daniels TJ, Fish D (1995) Effect of deer exclusion on the abundance of immature Ixodes scapularis (Acari: Ixodidae) parasitizing small and medium-sized mammals. J Med Entomol 32:5–11. https://doi.org/10.1093/jmedent/32.1.5
Daniels TJ, Fish D, Schwartz I (1993) Reduced abundance of Ixodes scapularis (Acari: Ixodidae) and Lyme disease risk by deer exclusion. J Med Entomol 30:1043–1049. https://doi.org/10.1093/jmedent/30.6.1043
Eccard JA, Pusenius J, Sundell J, et al (2008) Foraging patterns of voles at heterogeneous avian and uniform mustelid predation risk. Oecologia 157:725–734. https://doi.org/10.1007/s00442-008-1100-4
Elston DA, Moss R, Boulinier T, et al (2001) Analysis of aggregation, a worked example: Numbers of ticks on red grouse chicks. Parasitology 122:563–569. https://doi.org/10.1017/S0031182001007740
Falco RC, Fish D (1992) A comparison of methods for sampling the deer tick, Ixodes dammini, in a Lyme disease endemic area. Exp Appl Acarol 14:165–173. https://doi.org/10.1007/BF01219108
Flowerdew JR, Ellwood SA (2001) Impacts of woodland deer on small mammal ecology. Forestry 74:277–287. https://doi.org/10.1093/forestry/74.3.277
Gern L, Douet V, López Z, et al (2010) Diversity of Borrelia genospecies in Ixodes ricinus ticks in a Lyme borreliosis endemic area in Switzerland identified by using new probes for reverse line blotting. Ticks Tick Borne Dis 1:23–29. https://doi.org/10.1016/j.ttbdis.2009.11.001
Gilbert L (2010) Altitudinal patterns of tick and host abundance: a potential role for climate change in regulating tick-borne diseases? Oecologia 162:217–225. https://doi.org/10.1007/s00442-009-1430-x
Gilbert L, Brunker K, Lande U, et al (2017) Agriculture , Ecosystems and Environment Environmental risk factors for Ixodes ricinus ticks and their infestation on lambs in a changing ecosystem : Implications for tick control and the impact of woodland encroachment on tick-borne disease in livestock. Agric Ecosyst Environ 237:265–273. https://doi.org/10.1016/j.agee.2016.12.041
Gilbert L, Maffey GL, Ramsay SL, Hester AJ (2012) The effect of deer management on the abundance of Ixodes ricinus in Scotland. Ecol Appl 22:658–667. https://doi.org/10.1890/11-0458.1
Gilbert L, Norman R, Laurenson KM, et al (2001) Disease persistence and apparent competition in a three-host community: an empirical and analytical study of large-scale, wild populations. J Anim Ecol 70:1053–1061. https://doi.org/10.1046/j.0021-8790.2001.00558.x
Gray JS, Kahl O, Janetzki C, Stein J (1992) Studies on the Ecology of Lyme Disease in a Deer Forest in County Galway, Ireland. J Med Entomol 29:915–920
Hanincová K, Schäfer SM, Etti S, et al (2003) Association of Borrelia afzelii with rodents in Europe. Parasitology 126:11–20. https://doi.org/10.1017/S0031182002002548
Harrison XA (2015) A comparison of observation-level random effect and Beta-Binomial models for modelling overdispersion in Binomial data in ecology & evolution. PeerJ 3:e1114. https://doi.org/10.7717/peerj.1114
Harrison XA (2014) Using observation-level random effects to model overdispersion in count data in ecology and evolution. PeerJ 2:e616. https://doi.org/10.7717/peerj.616
Hartig F (2020) DHARMa: Residual Diagnostics for Hierarchical (Multi-Level / Mixed) Regression Models
Herder DM, Helle S, Niemelä P, et al (2016) Large Herbivore Grazing Limits Small-Mammal Densities in Finnish Lapland Large herbivore grazing limits small-mammal densities in Finnish Lapland. Ann Zool Fennici 53:154–164
Heylen D, Matthysen E, Fonville M, Sprong H (2014) Songbirds as general transmitters but selective amplifiers of Borrelia burgdorferi sensu lato genotypes in Ixodes rinicus ticks. Environ Microbiol 16:2859–2868. https://doi.org/10.1111/1462-2920.12304
Heylen D, Tijsse E, Fonville M, et al (2013) Transmission dynamics of Borrelia burgdorferi s.l. in a bird tick community. Environ Microbiol 15:663–673. https://doi.org/10.1111/1462-2920.12059
Hillyard P (1996) Ticks of North-west Europe. Field Studies Council
Hofmeester TR, Coipan EC, van Wieren SE, et al (2016) Few vertebrate species dominate the Borrelia burgdorferi s.l. life cycle. Environ Res Lett 11:043001. https://doi.org/10.1088/1748-9326/11/4/043001
Hofmeester TR, Jansen PA, Wijnen HJ, et al (2017a) Cascading effects of predator activity on tick-borne disease risk. Proc R Soc B Biol Sci 284:20170453. https://doi.org/10.1098/rspb.2017.0453
Hofmeester TR, Sprong H, Jansen PA, et al (2017b) Deer presence rather than abundance determines the population density of the sheep tick, Ixodes ricinus, in Dutch forests. Parasit Vectors 10:433. https://doi.org/10.1186/s13071-017-2370-7
Hoyer IJ, Blosser EM, Acevedo C, et al (2017) Mammal decline, linked to invasive Burmese python, shifts host use of vector mosquito towards reservoir hosts of a zoonotic disease. Biol Lett 13:20170353. https://doi.org/10.1098/rsbl.2017.0353
James MC, Gilbert L, Bowman AS, Forbes KJ (2014) The Heterogeneity, Distribution, and Environmental Associations of Borrelia burgdorferi Sensu Lato, the Agent of Lyme Borreliosis, in Scotland. Front public Heal 2:129. https://doi.org/10.3389/fpubh.2014.00129
Jones EO, Webb SD, Ruiz-Fons FJ, et al (2011) The effect of landscape heterogeneity and host movement on a tick-borne pathogen. Theor Ecol 4:435–448. https://doi.org/10.1007/s12080-010-0087-8
Klaus C, Gethmann J, Hoffmann B, et al (2016) Tick infestation in birds and prevalence of pathogens in ticks collected from different places in Germany. Parasitol Res 115:2729–2740. https://doi.org/10.1007/s00436-016-5022-5
Levi T, Kilpatrick AM, Mangel M, Wilmers CC (2012) Deer, predators, and the emergence of Lyme disease. Proc Natl Acad Sci 109:10942–10947. https://doi.org/10.1073/pnas.1204536109
Lewis F, Butler A, Gilbert L (2011) A unified approach to model selection using the likelihood ratio test. Methods Ecol Evol 2:155–162. https://doi.org/10.1111/j.2041-210X.2010.00063.x
Logiudice K, Duerr STK, Newhouse MJ, et al (2008) Impact of host community composition on Lyme disease risk. Ecology 89:2841–2849
Mysterud A, Easterday WR, Stigum VM, et al (2016) Contrasting emergence of Lyme disease across ecosystems. Nat Commun 7:11882. https://doi.org/10.1038/ncomms11882
Norman R, Bowers RG, Begon M, Hudson PJ (1999) Persistence of tick-borne virus in the presence of multiple host species: tick reservoirs and parasite mediated competition. J Theor Biol 200:111–118. https://doi.org/10.1006/jtbi.1999.0982
Ostfeld RS, Keesing F (2000) Biodiversity and disease risk: the case of Lyme disease. Conserv Biol 14:722–728. https://doi.org/10.1046/j.1523-1739.2000.99014.x
Pacilly FCA, Benning ME, Jacobs F, et al (2014) Blood feeding on large grazers affects the transmission of Borrelia burgdorferi sensu lato by Ixodes ricinus. Ticks Tick Borne Dis 5:810–817. https://doi.org/10.1016/j.ttbdis.2014.06.004
Paine RT (1980) Food Webs : Linkage , Interaction Strength and Community Infrastructure. J Anim Ecol 49:666–685. https://doi.org/10.2307/4220
Perkins SE, Cattadori IM, Tagliapietra V, et al (2006) Localized deer absence leads to tick amplification. Ecology 87:1981–1986. https://doi.org/10.1890/0012-9658(2006)87[1981:LDALTT]2.0.CO;2
Rijpkema SG, Molkenboer MJ, Schouls LM, et al (1995) Simultaneous detection and genotyping of three genomic groups of Borrelia burgdorferi sensu lato in Dutch Ixodes ricinus ticks by characterization of the amplified intergenic spacer region between 5S and 23S rRNA genes. J Clin Microbiol 33:3091–3095
Saïd S, Servanty S (2005) The influence of landscape structure on female roe deer home-range size. Landsc Ecol 20:1003–1012. https://doi.org/10.1007/s10980-005-7518-8
Smit R, Bokdam J, Ouden JD, et al (2011) Effects of Introduction and Exclusion of Large Herbivores on Small Rodent Communities Effects of introduction communities and exclusion of large herbivores on small rodent. Plant Ecol 155:119–127. https://doi.org/10.1023/A:1013239805915
Stanek G, Wormser GP, Gray J, Strle F (2012) Lyme borreliosis. Lancet 379:461–473. https://doi.org/10.1016/S0140-6736(11)60103-7
Takumi K, Sprong H, Hofmeester TR (2019) Impact of vertebrate communities on Ixodes ricinus-borne disease risk in forest areas. Parasit Vectors 12:1–12. https://doi.org/10.1186/s13071-019-3700-8
van Duijvendijk G (2016) The ecology of Lyme borreliosis risk. Wageningen University
van Wieren SE, Bakker JP (2008) The Impact of Browsing and Grazing Herbivores on Biodiversity. In: Gordon IJ, Prins HHT (eds) The Ecology of Browsing and Grazing. Springer Netherlands, pp 263–292
Vourc’h G, Abrial D, Bord S, et al (2016) Mapping human risk of infection with Borrelia burgdorferi sensu lato, the agent of Lyme borreliosis, in a periurban forest in France. Ticks Tick Borne Dis 7:644–652. https://doi.org/10.1016/j.ttbdis.2016.02.008
Watson J, Estes JA (2011) Stability, resilience, and phase shifts in rocky subtidal communities along the west coast of Vancouver Island, Canada. Ecol Monogr 81:215–239. https://doi.org/10.1890/10-0262.1
Zuur AF, Ieno EN, Walker NJ, et al (2009) Mixed effects models and extensions in ecology with R. Springer Netherlands

GandyetalESM1.docx

Download PDF

Version 1

posted

You are reading this latest preprint version

Experimental Evidence for Opposing Effects of High Deer Density on Tick-Borne Disease Prevalence and Hazard

Status:

Version 1

Abstract

Figures

1. Introduction

2. Materials And Methods

2.1. Experimental design

2.2. Quantifying questing I. ricinus nymph density (DON)

2.3. Quantifying rodent activity and tick burdens

2.4. Measurement of B. burgdorferi s.l. prevalence (NIP)

2.1.1. Statistical analyses

3. Results

Discussion

Declarations

References

Supplementary Files

Status:

Version 1