Climate, Deer, Rodents, and Acorns as Determinants of Variation in Lyme-Disease Risk

Risk of human exposure to vector-borne zoonotic pathogens is a function of the abundance and infection prevalence of vectors. We assessed the determinants of Lyme-disease risk (density and Borrelia burgdorferi-infection prevalence of nymphal Ixodes scapularis ticks) over 13 y on several field plots within eastern deciduous forests in the epicenter of US Lyme disease (Dutchess County, New York). We used a model comparison approach to simultaneously test the importance of ambient growing-season temperature, precipitation, two indices of deer (Odocoileus virginianus) abundance, and densities of white-footed mice (Peromyscus leucopus), eastern chipmunks (Tamias striatus), and acorns ( Quercus spp.), in both simple and multiple regression models, in predicting entomological risk. Indices of deer abundance had no predictive power, and precipitation in the current year and temperature in the prior year had only weak effects on entomological risk. The strongest predictors of a current year's risk were the prior year's abundance of mice and chipmunks and abundance of acorns 2 y previously. In no case did inclusion of deer or climate variables improve the predictive power of models based on rodents, acorns, or both. We conclude that interannual variation in entomological risk of exposure to Lyme disease is correlated positively with prior abundance of key hosts for the immature stages of the tick vector and with critical food resources for those hosts.


Introduction
Many emerging and re-emerging infectious diseases of humans are zoonoses transmitted by vectors. Examples include West Nile virus, ehrlichiosis, anaplasmosis, tick-borne encephalitis and Lyme disease. In each case, the vectorusually a mosquito or tick-acquires the pathogen from a vertebrate host during a blood meal taken early in the life cycle and becomes capable of transmitting it to humans during a later blood meal. Risk of human exposure to the disease increases with increasing abundance and infection prevalence of the vectors [1,2]. For virtually all vector-borne zoonoses, disease incidence in humans varies substantially from year to year [3][4][5]. Determining the causes of interannual variation in entomological risk would facilitate the development and deployment of preventative measures, potentially reducing the burden of disease.
Lyme disease is the most frequently reported vector-borne disease in the US [6]. Lyme disease is most prevalent in northeastern and north-central regions where suburban and exurban development encroaches on deciduous forest ecosystems that support the pathogen, vector, and their vertebrate hosts [7]. The etiological agent is a spirochete bacterium, Borrelia burgdorferi, which is transmitted by ticks in the Ixodes ricinus complex. In the eastern and central US, the vector is the blacklegged tick, Ixodes scapularis. I. scapularis is a three-host tick, requiring three blood meals, one each as a larva, nymph, and adult, to fulfill its life cycle. Larval ticks hatch in midsummer, typically uninfected with B. burgdorferi, and begin seeking a host for their initial blood meal. After feeding on a vertebrate host for several days, the larvae drop off the host and molt into the nymphal stage, which undergoes diapause for almost a year before becoming active and seeking a host the following late spring or early summer. Both larvae and nymphs are highly nonspecific in their choice of hosts, parasitizing dozens of species of mammal, bird, and lizard [8]. After feeding to repletion, nymphs drop off the host and molt into the adult stage, which seeks a medium-or large-mammal host in midautumn of the same year. Infection with B. burgdorferi can be acquired from the host during either the larval or nymphal blood meal, and both nymphs and adults are capable of transmitting infection to a vertebrate host, including humans.
Owing to its tiny size, potentially high abundance, and summer feeding, the nymphal stage is most likely to transmit B. burgdorferi to people and hence is responsible for the great majority of Lyme-disease cases [9]. Risk of human exposure to Lyme disease, given entry into habitats where ticks occur (mainly forests; [10]), is a function of the density of infected nymphal ticks, which in turn is the product of the total density of nymphs (DON) and the nymphal infection prevalence (NIP) [11]. Determining the causes of variation in density of infected nymphs (DIN) and its component parts is an important goal with both ecological and epidemiological implications. Human behavioral patterns, such as time spent in forest habitat and protective measures taken against exposure to ticks, also influence Lyme-disease risk but are beyond the scope of this study.
Prior studies of the factors influencing DIN and DON have focused largely on variation in climate and in the abundance and distribution of white-tailed deer (Odocoileus virginianus). Because ticks spend greater than 95% of their lives on the forest floor either digesting the blood meal, undergoing diapause, or seeking a host, exposure to ambient temperature and humidity could be important to survival rates and population dynamics [12][13][14]. In the laboratory, ticks experience high mortality when exposed to low humidity and high temperatures [15]. Consequently, hot and/or dry springs and summers have been postulated to reduce subsequent nymphal tick densities and Lyme-disease risk [16][17][18].
Because adult I. scapularis feed predominantly on whitetailed deer [19], much research has evaluated the impact of variation in abundance of deer on abundance of ticks. When deer are eliminated from some habitats by hunting or fencing, the abundance of ticks typically is strongly reduced [20][21][22]. Studies comparing natural variation in deer abundance with that in tick abundance are less conclusive; some have shown strong associations [23][24][25], whereas others have not [26][27][28].
Less attention has been paid to the potential effects of variation in abundance of hosts for larval ticks in influencing variation in DIN, DON, and NIP. This is perhaps a consequence of the lack of specialization by larvae on any particular host species. However, larval I. scapularis feed abundantly on white-footed mice (Peromyscus leucopus), and this host is the most competent natural reservoir for B. burgdorferi [29,30]. High feeding success on mice combined with high reservoir competence has led some researchers to postulate that Lyme-disease risk will vary with mouse abundance [31,32]. Although eastern chipmunks (Tamias striatus) host many larval ticks and are competent B. burgdorferi reservoirs [33], the impacts of variation in chipmunk abundance on Lyme-disease risk have been even less thoroughly explored (but see [34,35]).
The food resources for tick hosts might also be important to Lyme-disease risk. Oak trees (Quercus spp.) that dominate many forests in the US Lyme-endemic zone are known to produce highly variable acorn crops, a phenomenon known as masting. Acorns comprise a crucial resource for several vertebrate species, including white-footed mice, eastern chipmunks, and white-tailed deer, and can influence population density of the rodents [32,[36][37][38] as well as space use by deer [39]. Jones et al. [40] described a dramatic increase in abundance of larval I. scapularis following experimental simulation of a masting event, but the scale of the experiments was insufficient to assess longer-term impacts on abundance of nymphs. In a follow-up, Ostfeld et al. [41] described a positive correlation between metrics of Lymedisease risk and both prior-year mouse abundance and acorn abundance 2 y previously.
No prior study has assessed simultaneously the effects of variation in temperature, precipitation, deer, mice, chipmunks, and acorns, on variation in entomological risk of exposure to Lyme disease ( Figure 1). Assessments of subsets of these putative causal variables are characterized by relatively short time series, which have limited power to assess the influence of each factor separate from the others. Here we use long-term monitoring of these parameters combined with model comparison approaches to address the causes of variable risk in an area of high Lyme-disease incidence.

Results
The principal entomological Lyme-disease risk factor, DIN, varied by an order of magnitude among years, ranging from 1.07 infected nymphs 3 100 m À2 in 1997 (averaged across all six plots) to 10.00 3 100 m À2 in 1996. This variation was due primarily to variation in DON, which varied 6-fold among years (3.59 to 21.07 3 100 m À2 ). In contrast, NIP varied less than 2-fold among years, from 0.24 in 2005 to 0.45 in 1999.

Models for Total DON
Effect of growing degree days in the previous year (GDD tÀ1 ) on DON in the current year was positive but very weak in both the magnitude (estimated slope of the regression) and the strength of evidence for the effect (Table 1). Total growing season precipitation in the current year (PPT t ) also weakly influenced DON, but the effect was Gaussian with a peak at 223 mm. A model in which these two weather terms are combined multiplicatively was stronger than either of the univariate models, but explained only 15% of the variance ( Table 1). None of the other climate variables, deer variables, or prior year's density of larvae (DOL tÀ1 ) had any effect on DON; i.e., Akaike's information criterion corrected for small sample size (AIC corr ) values were much higher than those of the means model. Abundance of acorns tÀ2 , mice tÀ1 , chipmunks tÀ1 , and rodents tÀ1 (sum of mice and chipmunks) all independently influenced DON, and in all cases linear models were superior to exponential or power functions. The best univariate model for DON was a simple linear model of chipmunks tÀ1 , which explained 40% of the variance (Table 2; Figure 2). In no case was a multiple regression model superior to the chipmunk model ( Table 2).

Models for NIP
Among all univariate models, NIP responded only to the density of acorns in year tÀ2 (Figure 3), and this relationship explained only 16% of the variance in NIP. Models incorporating the climate variables, deer variables, DOL, mice, chipmunks, and total rodents performed no better than the means model (Tables 3 and 4). Because none of the independent variables other than acorns tÀ2 produced an improvement over the means model, there was no justification for testing multiple regression models.

Models for DIN
Similar to the results for DON, the effect of GDD tÀ1 on DIN in year t was positive but very weak in both the magnitude and the strength of evidence for the effect (Table 5). Total PPT t also weakly influenced DIN, and again the effect was Gaussian with a peak at 213 mm of rainfall. A model in which these two weather terms are combined multiplicatively was stronger than either of the univariate models, but explained only 11% of the variance ( Table 5). None of the other weather variables, deer variables, nor DOL tÀ1 had any effect on DIN, with AIC corr values higher than those of the means model. Abundance of acorns tÀ2 , mice tÀ1 , chipmunks tÀ1 , and rodents tÀ1 all independently influenced DIN. For acorns tÀ2 , mice tÀ1 , and rodents tÀ1 , the best models were nonlinear, but the difference in AIC corr (DAIC corr ) between the best model and corresponding linear model was always less than 1 (i.e., they had equivalent levels of support in the data), and the shapes of all nonlinear models were very close to linear (Table  5). Consequently, we tested linear combinations of these predictor variables in multiple regressions.
A model combining mice tÀ1 and acorns tÀ2 multiplicatively was the best model of all multiple regressions combining predictor variables that were supported by univariate analyses, and it explained 57% of the variance in DIN (Table  6; Figure 4A). The model with multiplicative effects of chipmunks tÀ1 and acorns tÀ2 was nearly as good, with DAIC corr ¼ 1.05, and a slightly higher R 2 (0.61) ( Figure 4B).

Host Responses to Acorns
Densities of both mice and chipmunks responded strongly to the prior year's acorn abundance, and in both cases the relationship was best described by saturating power functions ( Figure 5A and 5B). As a result of their similar responses to acorn abundance, the abundance of mice and chipmunks was strongly correlated among years (r ¼ 0.62). A linear model of the mouse-chipmunk correlation was superior to any of the nonlinear models, with a slope of 0.25 indicating that mice are consistently about four times as abundant as chipmunks ( Figure 5C). As a consequence of acorns influencing rodent abundances, and rodents influencing nymphal abundances, we observed a general pattern in which DON tracked rodents, and rodents tracked acorns, each effect displaying a 1-y lag ( Figure 6).  Given their relatively long generation times and low reproductive potential, it is unrealistic to expect deer population abundance to track annual variation in acorn production. However, the potential exists for deer to be attracted from nonoak to oak habitats by the presence of abundant acorns [39,42]. Because deer are an important host for adult ticks, and adult ticks lay eggs where they drop off their hosts, we expected that DOL in any given year would reflect deer space use the prior fall and should be positively correlated with acorn abundance. We found a weak, linear relationship between deer tÀ2 and DOL tÀ1 , with this model being only a slight improvement over the means model (DAIC corr ¼ 0.39; R 2 ¼ 0.04). In addition, no relationship existed between acorns tÀ2 and DOL tÀ1 ; all models, whether linear or nonlinear, were worse than the means model.

Discussion
Climate, deer, and acorns each have been proposed as primary determinants of temporal variation in risk of human exposure to Lyme disease, as measured by abundance and Borrelia-infection prevalence in nymphal Ixodes ticks. Using a model comparison approach and a 13-y dataset, we found weak support for climate variables, no support for deer, and strong support for an effect of acorns, mediated by acorn effects on white-footed mice and eastern chipmunks, which host many larval ticks and are competent reservoirs for B. burgdorferi.

Climate Variables
Of the four climate variables, two influenced DON and DIN, but they did so in unanticipated ways. Both DON and DIN increased linearly, albeit weakly, with increases in the prior year's temperature (GDD tÀ1 ). This result conflicts with the expectation that heat-caused mortality is an important regulator of tick abundance [15,43], but is consistent with the finding that (detrended) incidence of Lyme disease in people is positively correlated with summer temperatures in the prior year [44]. Both DON and DIN were influenced, again weakly, by precipitation in the current but not prior year, but intermediate levels of precipitation favored highest nymphal abundances. Again, this result conflicts with the expectation that tick survival increases linearly with ambient moisture [18]. Low DON and DIN in years of high precipitation could be caused by either flood-induced mortality or that caused by natural enemies (e.g., fungi) that are facilitated by high moisture [45]. As expected, none of the climate variables influenced NIP. NIP is determined by the proportions of larval tick meals taken from the various host species, which vary strongly in their reservoir competence. It seems unlikely that climate variables will influence the choice of hosts by larval ticks or relative abundances of different hosts. It remains possible that climate variables other than the ones we examined influence stage-specific tick survival and abundance and that the weak or absent effects we observed are a result of not including more important variables. Our analytical approach could easily support assessments of other variables for which there is some a priori expectation of an effect.

Deer Variables
The assertion that variable deer abundance is responsible for variable abundance of blacklegged ticks and hence Lymedisease risk has become almost axiomatic [46,47]. This seems   (3) abundance of nymphs in year t is correlated with that of larvae in year tÀ1. Here, we found that more than 3-fold variation in our indices of deer abundance did not affect  subsequent nymph abundance. This observation supports the assertion that, once deer abundance exceeds a low threshold value, further increases in deer density have little if any effect on nymphal densities [48]. We did observe a weak effect of deer on subsequent larval abundance, but larval abundance in any given year had no impact on next year's abundance of nymphs. This clear but surprising result indicates a decoupling of stage-specific abundances and suggests instead that abundance of nymphs depends on larval feeding opportunities the prior year (see below). The lack of demographic forcing from larval to nymphal stage each year also suggests that the effect of larval host abundance in any given year should penetrate only to the next year and not beyond.

Acorn and Rodent Variables
Previous research [32,40] supported an effect of acorns in year tÀ2 on nymphs in year t via two pathways, one in which abundant acorns tÀ2 boosted larvae tÀ1 by enhancing deer tÀ2 , and the other in which acorns tÀ2 boosted rodents tÀ1 , which in turn elevated nymphs t (Figure 1). Our results strongly support the second, rodent-driven pathway ( Figure 6) and refute the  first, deer-driven one. Although univariate models with GDD tÀ1 , PPT t , acorns tÀ2 , mice tÀ1 , and total rodents tÀ1 all had support, the best model had chipmunks in the previous year as the sole predictor of DON. A stronger role for chipmunks than for mice was somewhat surprising given the lower population densities of chipmunks ( Figure 4) and only modestly higher average larval burdens on chipmunks [33] at our sites. However, less intense grooming of larval ticks by chipmunks [49] could facilitate larval survival to the nymphal stage. Elsewhere in the North American range of Lyme disease, chipmunks have been postulated to play a critical role in the enzootic cycle [34,35,50].
In contrast to previous results based on a shorter time series [41] and to expectations based on the high reservoir competence of mice and chipmunks [33], neither mice tÀ1 nor chipmunks tÀ1 influenced NIP. Instead, we found that the univariate model of acorns tÀ2 was the only one supported by the data. Empirically based models [51] indicate that NIP varies strongly with variable species composition in the community of hosts for larval ticks. Consequently, an impact of acorns but not of single host species would be expected if acorns alter the species composition in the host community. In addition to mice and chipmunks, acorns are likely to influence abundance and space use by gray squirrels (Sciurus carolinensis), raccoons (Procyon lotor), and turkeys (Meleagris gallopavo). The former two species are incompetent reservoirs [33], and the latter is unlikely to be a competent reservoir [52]. Therefore, by influencing a suite of hosts that both stimulate and depress NIP, acorns might be expected to explain more of the variation in NIP than would any single host species.
Because DIN is the product of DON and NIP, one might expect that the best explanatory model would be more complex than those for its component parts. Indeed, we found that the model of multiplicative effects of mice tÀ1 (the second best predictor of DON) and acorns tÀ2 (the best predictor of NIP) had the most support. The model with multiplicative effects of chipmunks tÀ1 (the best predictor of DON) and acorns tÀ2 was almost as strongly supported.
Tick abundance and Lyme-disease risk can be high in habitats with few or no oaks [21]; therefore, we do not expect that acorn abundance will be a universal predictor of risk [44]. However, small rodent hosts are almost always involved in risk as both permissive tick hosts and competent Borrelia reservoirs. Population densities of mice and chipmunks vary dramatically among years in both oak-dominated and nonoak-dominated forests [53], and we expect that the increase in Lyme-disease risk that accompanies high rodent densities should be widespread, no matter the causes of rodent fluctuations. Acorns provide a convenient leading indicator of rodent abundance, and seeds of other tree species [53] or predators [54,55] might act similarly.
Previous studies of the determinants of variable Lymedisease risk or incidence have tended to focus on one or a small number of potential independent variables, and statistically significant effects of both climate [16][17][18] and deer [23][24][25] have been described. When we included candidate climate variables, deer indices, and larval tick density together with densities of mice, chipmunks, and acorns, our model comparison methods never selected models incorporating climate, deer, or larvae and always selected models with either rodents, acorns, or both, as explanatory variables. Effects of variable climate or deer abundance on Lyme-disease risk might achieve statistical significance without being biologically important if statistical models fail to assess more potent independent variables. It remains possible that variable climate and deer abundance affect large-scale spatial variation in Lyme-disease risk even if their impacts on temporal variability are weak. Our longterm studies of Lyme-disease risk at the epicenter of the epidemic in North America strongly implicate a role for population density of rodent hosts and their food resources. We suggest that masting indices could be used to alert the public when years of high Lyme-disease risk are anticipated.

Materials and Methods
Study sites. Field studies were conducted on the property of the Institute of Ecosystem Studies (IES) in Dutchess County, southeastern New York (lat 41 8 509N, long 73 8 459W), in the center of the northeastern US endemic zone for Lyme disease. Dutchess County has had among the highest incidence rates of Lyme disease in the US Note that, typically, chipmunk density tracks acorn density with a 1-y lag, and DON tracks chipmunk density also with a 1-y lag. DOI: 10.1371/journal.pbio.0040145.g006 during the past 10 y [56]. IES forests are typical of the eastern deciduous forests of New York and New England, dominated by oaks (Quercus rubra and Quercus prinus) in the overstory (57%-70% oak relative basal area; [42]), with oak and sugar maple (Acer saccharum) seedlings, maple-leaved viburnum (Viburnum acerifolium), witch hazel (Hamamelis virginiana), and ironwood (Ostrya virginiana) common in the understory. Two 2.25-ha plots (150 m 3 150 m) were established in 1991, and four more were added in 1995 to comprise three pairs of plots with ca. 150 m separating members of a pair and more than 700 m separating pairs. Acorn sampling. Acorn abundance was monitored on each of the six plots by placing circular baskets under the canopies of mature oaks distributed throughout the plot. The original two plots had twenty 0.5-m 2 baskets, and the remaining four had twenty-five 1.0-m 2 baskets. Seed baskets were supported by monofilament line attached to nylon stakes and were resistant to seed predators. Intact, mature acorns were counted monthly during autumn of each year, and the total number of acorns from all baskets within a plot was divided by the total basket area to derive an estimate of annual acorn production on each plot. Full-season acorn data were collected from 1993 through 2004 on the original two plots and from 1999 through 2004 on the remaining four plots. Acorn density in year tÀ2 was used as an independent variable potentially affecting Lyme-disease risk factors.
Small-mammal sampling. Each year from 1991 (the two original plots) or from 1995 (remaining plots) through 2005 we have monitored abundance of small mammals at IES using capturemark-recapture techniques. On each plot we established an 11 3 11 point grid of Sherman live traps, with 15 m between trap stations and two traps per station, for a total of 242 traps per grid. Trapping was conducted for 2-3 consecutive days every 3-4 wk, generally from May to November of each year. Traps were baited with crimped oats (with sunflower seeds and cotton batting added during cold weather), set at 1600 hours and checked between 0800 hours and about 1200 hours the following morning. This schedule allowed us to capture both diurnal (chipmunks) and nocturnal (mice) small mammals. These two species comprised more than 90% of captures. Small mammals were marked with individually numbered metal ear tags and released after handling at the point of capture. Data on age, sex, body mass, ectoparasite burden, and trap station were recorded on each capture. Protocols for animal handling were approved annually by an Institutional Animal Care and Use Committee.
We estimated population densities of white-footed mice and eastern chipmunks by inputting data from all trap sessions in a year into the Jolly-Seber open population model in program POPAN5 [57]. We selected the Jolly-Seber model that incorporates individual heterogeneity in capture probability. Because we were interested in assessing the importance of rodent abundance on nymphal tick abundance the following year, we focused on estimating rodent densities in midsummer, which is the time of peak activity of larval I. scapularis at our sites [10]. In order to create a standard metric when actual trapping sessions varied in time, we estimated mouse and chipmunk abundance for each grid and year on August 15 by linear interpolation between rodent abundance estimates for trap sessions immediately before and after August 15. Rodent densities in year tÀ1 were used as independent variables potentially affecting Lymedisease risk factors.
Deer abundance estimates. Deer abundance was estimated at two levels, one at the scale of the IES property and the other at the scale of individual plots. For estimating property-wide deer abundance annually, we used population estimates from the IES limited-access bow-hunting program, which has run continuously from 1987 to the present [58]. Between seven and 11 hunters per year hunted an average of 40 h each (range, 29-54 h) between mid-October and mid-November. All hunters were IES staff members or volunteers working with the staff wildlife biologist. Each was given exclusive access to one of 34 discrete hunting areas averaging 22 ha (range, 9.3-35.7 ha). Virtually all hunting was from tree stands. Each day, hunters reported the number of hours hunted and the number of deer sighted, and these data were converted to deer observed per hour hunted. As a validation of this method, we asked whether deer observed per hour by bow hunters were correlated with deer counts from annual autumn spotlighting surveys conducted from 1987 to 2000 (details in [58]) and found that the two census methods were highly correlated (r ¼ 0.70, df ¼ 11, p ¼ 0.01; [58]).
To estimate deer distribution on a smaller scale, we used deer browse surveys conducted each spring ) at a number of sites (range, 38-50) distributed throughout the IES property. Commonly browsed tree species or genera [59,60] have served as an index to trends in browsing rates. These species/genera include red maple (Acer rubrum), sugar maple (A. saccharum), serviceberry (Amelanchier arborea), black birch (Betula lenta), black cherry (Prunus serotina), oaks (Quercus spp.), hickories (Carya spp.), and ashes (Fraxinus spp.). Browsed woody stems detected in spring reflect the distribution of deer foraging activity during the previous autumn and winter. Data were collected along transects of unrestricted width spaced either 200 m or 300 m apart. Starting points were associated with corners of a 100-m grid overlay of the property. A major compass direction was randomly drawn each year for transect direction. As index species were sighted, all buds below 2 m were counted and examined for browsing by deer. At least 100 buds in total were examined for each site. The percentage of available stems browsed was calculated for each site and for each index species. We selected browse survey plots closest to each of the trapping grids (less than 100 m distant) to estimate deer activity specific to each of the plots. Deer browsing intensity in year tÀ2 was used as an independent variable potentially affecting Lyme-disease risk factors.
Climate variables. An almost infinite number of climate variables (temperature, precipitation, minimum, maximum, variance, mean, specific to months or seasons, etc.) potentially could influence tick survival and densities of nymphs. Consequently, the probability of uncovering spurious relationships between climate and tick abundance is quite high if many explanatory variables are included without clear a priori justification [61]. To avoid this problem, we selected climate variables that have been reported to influence either changes in the abundance of immature blacklegged ticks or in human Lymedisease incidence [16,17]. The climate variables selected represent temperature and precipitation conditions in either the current year (year t) or the prior year (year tÀ1). Values for year t reflect potential impacts of climate on survival of the current year's nymphal stage, whereas those for year tÀ1 represent possible effects of climate on survival of the prior year's larval stage. Following Jones and Kitron [16], for temperature data, we included only those months from the beginning of the growing season through the end of the activity season for the appropriate life stage (July for nymphs, September for larvae; [16,41]). Specifically, we used cumulative growing degree days (GDD) for March 1-June 30 of year t, and cumulative GDD for March 1-September 30 of year tÀ1. For precipitation, we included only those months from the onset of potential soil moisture limitation (May) to the end of the appropriate activity period [16]; i.e., PPT (in mm) for May 1-June 30 of year t; and PPT for May 1-September 30 of year tÀ1 [16]. All climate data came from the IES environmental monitoring station (http://www.ecostudies.org/emp_purp.html), located less than 1.5 km from our field plots Tick and Borrelia sampling. Estimates of the abundance and infection prevalence of nymphal ticks comprised the response variables of interest. In addition, estimates of larval abundance in year tÀ1 were used as an independent variable. We monitored the abundance of larval and nymphal ticks in each plot and year by dragging 1-m 2 white corduroy drag cloths [62] along 450 m of transects approximately every 3 wk from April through November. Drag cloths were examined and all ticks counted and removed every 30 m. Frequent sampling in the early 1990s revealed that peak hostseeking activity for larvae occurred in mid-to late August, and for nymphs in mid-to late June, and we timed our annual sampling to coincide with these peaks. For each plot and year, we estimated larval and nymphal abundance as the peak density (ticks per 100 m 2 ). Peak densities were highly correlated with cumulative seasonal densities (correlation coefficients typically greater than 0.80; R. Ostfeld, unpublished data) on each plot.
Infection of individual ticks with B. burgdorferi was determined using direct immunofluorescence assay (DIA). Ticks were washed once in 70% ethanol and twice in deionized water and ground in phosphate-buffered saline (PBS). Three 5-ml aliquots of tick suspension were placed in separate wells in a multiwell slide, airdried, and fixed in cold acetone for 10 min. Fluorescein rabbit anti-B. burgdorferi conjugate was incubated in wells at 37 8C for 45 min, after which slides were washed in PBS, dried, and mounted with a coverslip. Slides were examined systematically to categorize each tick as either infected or uninfected. On average, 378 nymphs (range, 146-660) were examined for infection each year. Our estimates of tick infection based on DIA have been verified by virtually identical infection prevalence estimates (within 1.6 percentage points) for IES ticks based on PCR and reverse-line blotting [63,64].
Statistical methods. Our goal was to evaluate the strength of evidence for effects of a series of plausible independent variables on temporal variation in entomological risk of human exposure to Lyme disease ( Figure 1). Risk was measured by annual estimates at individual sample locations for the DON, NIP, and their product, the DIN. For each of the three response variables (DON, NIP, DIN) we compared the strength of evidence for a series of alternative, competing regression models using AIC corr . We created both linear and nonlinear (exponential, power function, Gaussian) models, as appropriate. We first compared evidence for each of the 11 independent variables separately by comparing the AIC corr of their regression models to the AIC corr value of an intercept-only (i.e., mean) model. The 11 independent variables included two temperature variables, two precipitation variables, two indices of deer abundance, abundance of larval ticks, acorn abundance, and three measures of small mammal abundance. We checked for colinearity among the independent variables using correlation coefficients. Because there were missing values for some variables in some years and plots, the univariate models were compared against mean models estimated with the same subset of (nonmissing) observations. We then tested series of increasingly complex models by combining sets of independent variables for which there was evidence (as measured by AIC) of univariate effects. Our choices for the forms of the multiple regression models were guided by the dictates of parsimony: while our dataset represents an enormous sampling effort over a 13-y period, the actual sample sizes (a given plot in a given year) were still relatively small, ranging from 42 to 58 observations for the various models.
We used simulated annealing (a global optimization algorithm) to find the maximum likelihood estimates for the parameters of each model. The observations were assumed to be normally distributed with a homogeneous variance. Examination of the residuals indicates that this assumption was appropriate for all three of the response variables. The specific annealing algorithm we used (based on Goffe et al. [65]) allows a bounded search in cases where a range of parameter values is either mathematically inconsistent or biologically unreasonable (i.e., negative densities or variance estimates). The sample locations were widely enough distributed that we are confident in assuming that the error terms in our models are spatially independent, and there was no evidence of temporal autocorrelation in the residuals. The simulated annealing and all of the statistical analyses were done using R (http://www.r-project.org).