The effect of climate on population growth in a cold-adapted ungulate at its equatorial range limit

. Climatic changes are affecting the distribution and viability of species worldwide, and the effects may be greatest for heat-sensitive organisms in populations situated near the species ’ equatorial range limit. We studied the population dynamics of a cold-adapted large herbivore, moose ( Alces alces shi-rasi ), in a population located at the extreme southern range limit of the species in Utah, USA, using a long-term dataset of aerial counts conducted between 1958 and 2013. We used a modeling approach that acknowledges the uncertainty in the number of moose counted. To determine how climate in ﬂ uenced this population, we tested four models corresponding to different hypotheses suggested by previous studies of moose population dynamics: (1) High summer rainfall increases population growth, (2) severe winters reduce population growth, (3) high summer and winter temperatures cause heat stress which reduces population growth, and (4) snow conditions favorable to winter ticks reduce population growth. We then ranked these models against two additional hypotheses that tested the combined effects of the best-performing models. The best-supported model included summer rainfall and the number of days with snow cover in late winter, which described the conditions in ﬂ uencing winter tick numbers, a common parasite of moose in the region. Reproductive female ticks drop off their hosts at the end of winter to lay eggs, and fewer ticks survive in years with abundant snow cover. Positive effects of snow and rainfall indicated that moose population growth was higher following summers with more rainfall and late winters with more days of snow cover, the latter because those conditions likely reduced winter tick numbers. In accor-dance with global patterns in which altered precipitation regimes are in ﬂ uencing the dynamics of many species, the top-performing model suggested that both summer and winter precipitation acted together to explain the most variation in moose population growth. Our analysis demonstrates the multiple pathways by which climate and population density can affect the dynamics of temperate species living at their equatorial range limits, including potential parasite-mediated effects.


INTRODUCTION
A substantial body of research indicates that climatic changes are altering the distribution and viability of species across many taxa worldwide (Bellard et al. 2012). Evidence that species in the Northern Hemisphere have shifted their ranges toward higher latitudes and elevations (Chen et al. 2011) indicates that populations at the southern range periphery are at high risk from environmental stress. Anticipating the responses of species to future environmental change is challenging, however, because there is enormous taxonomic variation in species response to changes in climate (Pearce-Higgins 2015), and certain species may be able to compensate for changing environments (Doak and Morris 2010). This is partly because climate can affect organisms through different pathways: either by directly affecting individual physiological processes or indirectly through changes in ecosystem processes . Consequently, the effects arising from either pathway may depend on the life history of the organism (Chen et al. 2011).
Large herbivores inhabiting seasonal environments have forged an intricate relationship with climate over evolutionary time (Coulson et al. 2000, Post and Forchhammer 2002, and empirical evidence indicates both direct and indirect climate effects on their population dynamics. Direct influences of warming temperatures can cause individuals to alter foraging strategies in response to high ambient temperatures (Owen-Smith 1998, Aublet et al. 2009), in some cases at the expense of body mass (van Beest and Milner 2013). Indirect effects of warming can create mismatches between the timing of reproduction and energy balance stemming from phenological advancement of vegetation growth Forchhammer 2008, Moyes et al. 2011). Further complexity occurs when indirect effects of climate influence competition for limited resources and population density acts to modulate the effect of climate (Bonenfant et al. 2009). The influences of climate on large herbivores are therefore nuanced and difficult to predict. Populations within a given species may even experience disparate responses to the same climate pattern (Grøtan et al. 2009). This would especially hold true if the local climates differed between populations with respect to some species-specific optimum.
Moose (Alces alces) are a species that have evolved specific climatic tolerances, favoring cool and wet environments (Kelsall andTelfer 1974, Renecker andHudson 1986). Populations in North America occupy a large region spanning >30°latitude and in which the climate varies dramatically (Telfer 1984). Moose demographic rates often relate positively to precipitation (Thompson 1980, Murray et al. 2006, Monteith et al. 2015 and negatively to warm-season temperatures (Crête and Courtois 1997, Murray et al. 2006, Grøtan et al. 2009) across their distribution. The effects of winter climatic conditions, however, are more inconsistent. Severe winters are limiting to moose in the northern reaches of their distribution (Thompson 1980, Crête and Courtois 1997, Keech et al. 2000, Sivertsen et al. 2012), but a southernedge population in Minnesota has shown just the opposite-that warm winter temperatures negatively affect moose survival (Lenarz et al. 2009; but see Fieberg 2014 andMech et al. 2018). Monteith et al. (2015) showed that warmer annual temperatures decreased recruitment rates of southern moose populations in the western United States. Moose are also susceptible to diverse climate-related diseases and parasites (Murray et al. 2006), and the degree of susceptibility may vary across the species' range.
Moose located near the species' southern range limit offer a unique case study to investigate the effects of climate on peripheral populations of cold-adapted species (Ditmer et al. 2017). Several studies of moose at their southern range limit in the central United States have shown declines in recent decades due in part to climate change (Murray et al. 2006, Lenarz et al. 2009 consistent with global patterns of climate-caused impacts on wild populations (Chen et al. 2011). By contrast, moose in the western United States have expanded their range southward during the last century despite warming temperatures (Darimont et al. 2005, Wolfe et al. 2010). Many of these recently colonized populations experienced phases of rapid growth after exploiting new habitats but have since stabilized or begun to decline (Monteith et al. 2015). Thus, the longterm viability of these recently established populations at the species' southern range limit is not known.
We used aerial count data spanning 56 yr to test hypotheses about the influence of climate and population density on the population growth of a harvested moose population in Utah that occurs at the species' extreme southern range limit. We tested four hypotheses about the effects of climate on moose population growth: Summer rainfall increases population growth, severe winters reduce population growth, high temperatures during both summer and winter induce heat stress and reduce population growth, and late winters with few days of snow cover, that is, conditions favoring winter tick abundances, reduce population growth (Fig. 1, Table 1). We then constructed two additional models composed of the top-ranking individual climate variables to evaluate whether population dynamics were better explained by more than one climate effect. By quantifying the effects of climate and population density on moose in Utah, we hope to shed light on the factors affecting populations of heat-sensitive northern ungulates at their equatorial range limit.

Study area
We analyzed data collected in northern Utah, USA (39.7°to 41.9°N, À109.2°to 112.2°W; Fig. 2). Much of the suitable moose habitat within the study area occurs within the Uinta-Wasatch-Cache and Ashley National Forests. Prominent mountain ranges in the study area included the Wasatch and Uinta Mountains. Habitat was representative of the intermountain west and Rocky Mountains and included riparian communities, shrub-dominated uplands, and coniferous forests (Wolfe et al. 2010). Riparian habitats were generally dominated by willow (Salix spp.) communities, whereas uplands were characterized by Gambel oak (Quercus gambelii), sagebrush (Artemisia spp.), and mountain mahogany (Cercocarpus spp.). Forested habitats were dominated by Douglas fir (Pseudotsuga menziesii), lodgepole pine (Pinus contorta), subalpine fir (Abies lasiocarpa), and quaking aspen (Populus tremuloides). Additional climate and habitat details are provided by Wolfe et al. (2010) and Ruprecht et al. (2016).

Data
In 1957, the Utah Division of Wildlife Resources (UDWR) initiated winter aerial surveys to obtain minimum statewide abundances of various big game species. Moose had only recently colonized the state and were believed to be present in only one management unit (Wolfe et al. 2010). Over time, moose expanded into neighboring management units and presently occupy at least 16 management units. However, many of these units have been colonized only recently, so sufficient count data for the current study were only available from the 10 units in the northeastern part of the state where moose first established (Fig. 2). Surveys were initially conducted with fixed-wing aircraft (1957)(1958)(1959)(1960)(1961)(1962) and subsequently with helicopters . At first, surveys were conducted every year and later switched to an irregular schedule in which only a subset of units was annually surveyed conditional on good survey conditions (i.e., adequate snow cover). Each unit was surveyed on average every third year. All suitable habitat within each management unit was surveyed according to expert opinion of biologists, and routes were flown as consistently as possible across years. Observers counted every moose they sighted and classified moose by sex and age class (juvenile or adult) when possible. We analyzed the total number of moose counted regardless of sex or age because censoring counts of unknown sex or age class would introduce substantial bias. Wolfe et al. (2010) provide more details on survey methodology.

Model specification
Failing to acknowledge uncertainty in the number of animals counted can introduce bias into population growth models with density dependence (Freckleton et al. 2006, Lebreton andGimenez 2013). Further, the periodic nature of the surveys meant that missing data years were common. To accommodate these issues, we used statespace models which in some cases can separate noise arising from biological processes from that which arises due to imperfect detection (Buckland et al. 2004, Dennis et al. 2006, hereafter "process error" and "observation error." State-space models can be implemented using both frequentist and Summer rainfall Higher rainfall in summer tÀ2 increases forage abundance, which in turn increases reproductive success in year tÀ1 . Densitydependent food competition suggests an interaction between rainfall and population size. We expect the interaction term to be negative such that the positive effect of rainfall on population growth will be reduced as moose density increases. Jul-Aug rainfall (tÀ2) + Jul-Aug rainfall (tÀ2) 9 N moose (tÀ1) Thompson (1980), Murray et al. (2006), Monteith et al. (2015) Winter severity Winter mortality increases in years with deep snow due to higher energetic costs of locomotion and limited access to browse. Density-dependent food competition suggests an interaction between snow depth and population size. We expect the interaction term to be negative such that the negative effect of snow depth on population growth increases as moose density increases.
# Days > 0°C winter (tÀ1) + # Days > 32°C summer (tÀ1) Crête and Courtois (1997) Winter ticks Years with poor snow cover favor winter tick survival and reproduction, which increases tick-related moose mortality in subsequent years ( Fig. 1). We expect the number of days with snow cover will increase population growth. Correlation between tick abundance and host density implies an interaction between number of days with snow and moose population size. We expect the interaction to be negative such that the positive effect of snow cover will be reduced as moose density increases # Days snow Feb-Mar (tÀ2) + # Days snow Feb-Mar (tÀ2) 9 N moose (tÀ1) Drew and Samuel (1986), Garner and Wilton (1993), DelGiudice et al. (1997), Samuel (2007) ❖ www.esajournals.org Bayesian approaches, but we opted to use a Bayesian analysis because of the ease in accommodating missing data (K ery and Schaub 2011). We used a discrete-time Gompertz population growth model (Eq. 1) because it has been widely used in state-space models of time series of animal counts and exhibits several desirable qualities. It estimates the strength of density dependence (b N moose ) by incorporating an autoregressive term describing the effect of the population size from the previous year (N tÀ1 ). When b N moose < 0 in our formulation of the Gompertz model, population density decreases population growth; when b N moose > 0, population density increases population growth (i.e., an Allee effect), and when b N moose = 0, density independence is assumed (Dennis et al. 2006).
The parameter b rmax is estimated from the model and is equivalent to the maximum intrinsic rate of increase (Dennis et al. 2006), that is, the maximum growth rate a given species could attain if resources were unlimited (Hone et al. 2010). We denoted the raw counts as Y and the population estimates as N. The model is indexed by the subscript i representing each of the 10 surveyed management units, and the subscript t denoting each of the 56 yr spanning the survey period .
We included a term in Eq. 1 to account for known annual harvests of moose (Colchero et al. 2009, Koons et al. 2015). The parameters b climate measure the effect of climate variables. For shorthand, we refer to Eq. 1 as g(N i,tÀ1 ) and note that this function predicts the mean of the posterior distribution for the number of moose (N i,t ). A lognormal error term (r 2 p ) estimates the process error, or error not explained by the deterministic portion of the Gompertz model (Eq. 2).
Our raw data represented the minimum counts of moose which are almost certainly less than the true population size since it is usually not possible to count every animal on the landscape. Failing to account for this systematic undercounting could bias results. Further, our surveys were conducted using fixed-wing aircraft (1957-1962) and helicopters (1963-2013) and the two types of aircraft likely yielded different rates of detection. We therefore chose to model the observation process as a binomial random variable such that the observed count data Y i,t represented the true, latent population size N i,t after it was modified by detection probability p, and we allowed p to vary according to management unit and aircraft type (i.e., fixed-wing vs helicopter; Eq. 3). The estimates of population size (N i,t ) were rounded to the nearest integer away from zero in order to conform to the support of the binomial distribution.
This study did not have sufficient data to estimate detection probability so we adopted the moose sightability model developed from helicopter survey data by Anderson and Lindzey (1996). That study occurred in western Wyoming between 150 and 250 km from our study area. We believe the proximity and similar habitat justify the use of their sightability model for the years in our study surveyed by helicopter. The best model developed by Anderson and Lindzey (1996) included canopy cover as the single covariate influencing moose detection probability. Our survey data lacked the associated data on canopy cover for each moose detected which precluded the application of their sightability model at the individual level. Instead, we applied the model to the management unit level and, in doing so, made the simplifying assumption that each moose was observed at the average level of canopy cover in that unit.
We used satellite-derived data to estimate canopy cover for each of the 10 management units in our study to predict the detection probability of moose in each unit using the Anderson and Lindzey (1996) model. The resulting predictions of detection probabilities included variance which we sought to acknowledge instead of assuming that sightability was fixed at the predicted values. We therefore used moment matching  to find the parameters of the beta distributions matching the means and variances of the predicted detection probabilities and let this beta random variable serve as the prior distribution for p in the binomial observation model. Formulating the prior distribution in this way directly informs detection probability while also acknowledging uncertainty in the Anderson and Lindzey (1996) model, providing a flexible observation model.
We have no prior information on moose sightability from fixed-wing aircraft in similar habitat types as our study. However, a previous study compared fixed-wing and helicopter searches for moose and found that fixed-wing surveys detected just 56% of the moose counted by helicopter searches in the same area (Gosse et al. 2002). We therefore specified a prior distribution for fixed-wing detection probability to be 0.56 9 the helicopter detection probabilities. Fixed-wing surveys occurred when moose were believed to occupy only one management unit, so p fixed-wing is a single parameter whereas p helicopter is estimated separately for each unit.
We chose vague priors for two parameters for which we had no prior information: process error (r 2 p ) and the effects of climate (b climate ). We used more informative priors for two parameters that we could reliably estimate from ecological theory or could be derived from empirical data: strength of density dependence (b N moose ) and r max in moose (b rmax ). We followed the approach of Koons et al. (2015) in defining a prior distribution for b N moose to exclude values that exceed the mathematical limitations of the effect of population density on growth rate in the Gompertz model (Dennis et al. 2006). Also following Koons et al. (2015), we estimated an informative prior for b rmax (x ¼ 0:304, standard deviation [SD] = 0.08) based on the mean of five previous ❖ www.esajournals.org 6 February 2020 ❖ Volume 11(2) ❖ Article e03058 studies which provided estimates of r max in moose (Bergerud 1981, Keith 1983, Van Ballenberghe 1983, Cederlund and Sand 1991, Sinclair 2003. Informing the parameters with reliable prior information aids identifiability of other model parameters including process and observation errors (Lebreton andGimenez 2013, Koons et al. 2015). Finally, we summed the unit-level population sizes in each year to estimate population-wide abundance (Eq. 4).
A directed acyclic diagram of the model structure is provided as a visual aid in the appendix (Appendix S1: Fig. S1). The full model, including prior distributions, is specified by the following statement in which items in bold represent matrices:

Modeling approach
We first reviewed the literature documenting the effects of climate on North American moose populations, and from this knowledge base, we developed four hypotheses that most plausibly explained how climate may affect moose population growth in our study area (Table 1). We used climate variables from the National Oceanic and Atmospheric Administration's (NOAA) Annual Climatological Summary dataset (National Centers for Environmental Information, NESDIS, NOAA, US Department of Commerce). We calculated each climate variable as the mean of measurements from 7 climate stations within the study area. The winter severity model included data from January to March (winter), and the heat stress model included data from June to August (summer). The summer precipitation model included data from July to August which is when the least rain falls in our study area. The winter tick model included data from February to March which is when adult female ticks drop from their hosts and lay eggs to replenish the tick population for the subsequent year (Fig. 1).
We used a two-step approach to model selection. We first ranked the importance of each climate model individually. Next, we constructed models consisting of combinations of the highest ranked individual climate models to determine whether population dynamics were better explained by multiple mechanisms. To avoid overfitting of models and producing spurious results, we did not compare all subsets of climate variables.
Aerial counts were conducted between December and February according to good survey conditions (i.e., adequate snow cover; Wolfe et al. ❖ www.esajournals.org 2010), but for consistency, we assumed each count was conducted in January of year t; accordingly, any climate variable in the preceding 12 months was defined as a lag of one year because it was a different calendar year. For example, for a count conducted in January 2013, a winter variable with a lag of one year corresponded to climatic conditions between January and March of 2012. We standardized all climate variables to have mean = 0 and SD = 1 to assist with model convergence and interpretation. The coefficient estimates for each climate variable are therefore on the same scale and directly comparable. We also calculated the proportion of the posterior that has the same sign as the mean which we denote as f (Kellner 2018). This quantity represents our confidence that the parameter is indeed positive or negative and not null.

Model implementation
To estimate posterior distributions of the parameters of interest, we conducted Markov chain Monte Carlo (MCMC) simulations in JAGS (v. 4.3.0; Plummer 2012) accessed from program R using the jagsUI package (Kellner 2018). For each model, we ran three chains each consisting of 100,000 iterations with the first 50,000 discarded as burn-in and thinned the sample to retain every 50th simulation. Model convergence was assessed visually using traceplots and by ensuring that the b R value for each parameter of interest was less than 1.1 (Gelman 1996). If models did not reach convergence after 100,000 iterations, we updated the model with additional iterations until convergence was satisfactory. We assessed model goodness-of-fit using posterior predictive checks (Gelman 2004, K ery and Schaub 2011. To do this, we generated hypothetical count data (i.e., Y.new i,t ) from the model and used a chi-square statistic to compare the observed and expected values from the original and new datasets at every MCMC iteration, that is, where Y.est i,t N i,t 9 p i . Calculating the proportion of iterations in which the chi-square statistics arising from the original and hypothetical datasets are more extreme than one another provides a measure of goodness-of-fit; a value of 0.5 indicates perfect fit, and values close to 0 or 1 suggest a lack-of-fit. We then ranked the candidate models using posterior predictive loss (PPL; Gelfand andGhosh 1998, Hooten and, where PPL ¼ P 10 i¼1 P 56 t¼1 y i;t À E y i;t jdata The first term compares the observed data with the posterior predictive mean and provides a measure of model fit, and the second term represents posterior predictive variances and adds a penalty for model complexity (Gelfand and Ghosh 1998, Fieberg et al. 2013. Lower PPL scores indicate better model performance.

RESULTS
Posterior predictive checks indicated adequate fit for all models; that is, values were >0.4 and <0.6. These results imply that each candidate model was capable of generating data consistent with the observed data. Additionally, each model successfully converged such that b R values for each parameter were ≤1.1.
The northern Utah moose population exhibited variable population growth between 1958 and 2013 (Fig. 3). In general, a period of sustained growth was observed as the population grew from low density after the species colonized the area until ca. 1990, at which point the population presumably reached carrying capacity and the trend was subsequently punctuated by several periods of decline. Each of the 10 individual management units also experienced variable population dynamics over the course of the study (Appendix S1: Fig. S2).
In the first round of model selection comparing only the univariate climate hypotheses, the summer rainfall model ranked the highest based on posterior predictive loss scores, followed by the winter tick model, the winter severity model, and finally the heat stress model. We next fit two additional models composed of multiple climate mechanisms identified to be important from the first round and compared them to the single climate hypothesis models. The new models ❖ www.esajournals.org combined summer rainfall + winter ticks, and summer rainfall + winter severity. The summer rainfall + winter ticks model (PPL = 23,488.2; Table 2) ranked the highest out of the entire candidate set, that is, those considering both single and multiple climate mechanisms.
In the top model, the parameter describing the ability of the population to grow from low density, b rmax , was estimated at 0.291 (95% Bayesian credible interval [BCI] = 0.194, 0.391; Table 3), which was slightly lower than the mean of the r max values in previous studies that we provided as a prior. Population growth rate was negatively influenced by moose density in the previous year (b N moose = À0.041, 95% BCI = À0.060, À0.023, f = 1; Table 3), suggesting density dependence was operating in the population.

DISCUSSION
Our results suggest that the lagged effects of both summer and winter precipitation acted in tandem to influence the dynamics of this moose population. Specifically, the combination of summer rainfall and snow conditions influential to winter ticks acted as the primary driver of moose population growth. Our best model implies that moose population growth will be highest 2 yr after winters with more days of snow cover during January-March followed by a summer with abundant precipitation. This result suggests that a parasite-mediated effect may have nearly as strong of an influence on moose population dynamics in Utah than climate factors traditionally thought to influence moose such as winter severity or summer rainfall. Our analysis found little support for the effects of heat stress which are thought to be detrimental to moose across their range but especially at southern latitudes (Murray et al. 2006, Monteith et al. 2015, Ditmer et al. 2017.
Our analysis demonstrates that summer rainfall played an important role in moose population dynamics in our study area. Specifically, total precipitation in July and August increased subsequent population growth. Years with abundant rain during these months likely acted to prolong the growing season and provide moose access to a higher quality and/or quantity of forage for a longer period of time. This likely allowed moose to enter winter in better nutritional condition and buffered them against malnutrition during the time when resources were most limiting. Consistent with density-dependent food limitation, the interaction between summer rainfall and number of moose suggests that population growth is highest at low moose densities and high rainfall, but the positive effect of rainfall diminishes as population size increases (Fig. 4a). Numerous studies have determined summer precipitation to be an important factor influencing moose populations (Thompson 1980, Murray et al. 2006, Monteith et al. 2015, so this was not an unexpected finding. A more surprising result was that the snow conditions influencing winter ticks had nearly as strong of an influence as summer rainfall. However, the detrimental effects of winter ticks on North American moose have long been known. One study in Alberta, Canada, estimated the mean number of winter ticks on a moose was >32,000, though numbers can reach nearly Fig. 4. The relationship between lambda and climate predicted by the top model with harvest held constant at its median value. In both panels, the broken green line represents a low density of moose (2.5% quantile of predicted abundance) and the solid blue line represents a high density of moose (97.5% quantile of predicted abundance). Shaded bands represent the 95% Bayesian credible interval for the predicted effect. Units on the x-axis are represented in terms of the number of standard deviations from the mean of the climate variable. (a) The relationship between summer rainfall 1.5 yr before surveys were conducted, which describes density-dependent resource limitation. (b) The relationship between number of days with snow between February and March two years before surveys were conducted, a proxy for conditions detrimental to winter ticks. 150,000 on a single moose (Samuel and Welch 1991). Engorged ticks can negatively affect moose and other hosts through a number of pathways including behavioral modifications, blood loss, anemia, and possibly transmission of toxins (Samuel and Welch 1991). Moose actively groom themselves in response to irritation arising from engorged ticks, and this behavior may result in less time spent foraging or resting (Samuel and Welch 1991). Grooming is thought to be effective in reducing tick numbers but may have bioenergetic and thermoregulatory consequences if time spent foraging is reduced or if significant hair is lost in the process. Several studies have shown that the cumulative effects of tick infestations can reduce body fat stores and hamper weight gain (Addison and McLaughlin 1988, Addison et al. 1994, Musante et al. 2007). Winter ticks may be the most prominent mortality factor for adult (Bergeron and Pekins 2014) and juvenile (Jones et al. 2017) moose in New Hampshire. Current research on the demography of moose in Utah has also correlated high winter tick loads on adult female moose with lower reproductive success including poorer calf survival among dams with high tick loads (UDWR, unpublished data). The agency has also observed lower tick loads on moose following years with a heavy snowpack.
Few studies have correlated winter tick outbreaks with population crashes in moose despite the serious metabolic detriments of tick infestations to moose on an individual level. A notable exception is the well-studied Isle Royale moose population which has experienced several periods of decline coincident with epizootics of winter ticks (DelGiudice et al. 1997). Other tickrelated die-offs have been noted in Algonquin Provincial Park, Ontario , and Elk Island National Park, Alberta (Samuel 2007). However, no studies to our knowledge have assessed the importance of conditions related to winter ticks on moose population growth.
The positive effect of winter snow cover on subsequent moose population growth is consistent with literature regarding winter tick population dynamics. Numerous studies suggest a direct link between late winter weather conditions and future tick abundances. Specifically, snow cover in late winter reduces survival and reproduction of ticks; when ticks drop off their hosts during this time and land on snow, their chances of surviving and reproducing are Notes: Each estimate of p (detection probability) refers to a different management unit where 1 = North Slope Summit, 2 = North Slope Daggett, 3 = Chalk Creek, 4 = Morgan-Rich, 5 = Cache, 6 = Ogden, 7 = East Canyon, 8 = Kamas, 9 = South Slope, and 10 = Wasatch Mountains. f represents the proportion of the posterior with the same sign as the mean, that is, our confidence that the parameter is indeed positive or negative.
decreased (Drew and Samuel 1986, DelGiudice et al. 1997, Samuel 2007; Fig. 1). This is because ticks have low thermal tolerances and are susceptible to snow or cold temperatures (Samuel 2007). We cannot rule out the possibility of snow cover influencing moose through other pathways than tick-mediated effects, but our results aligned with our a priori hypothesis regarding tick effects. Further, our models describing other plausible climate effects on moose in winter (i.e., the winter severity and heat stress models) did not perform as well.
Reports that tick abundances track moose numbers (Peterson 1955, Samuel and Welch 1991, Samuel 2007, Bergeron and Pekins 2014 motivated our test of the interaction between moose density and winter snow conditions. The model predicted a negative interaction between moose density and number of days with snow which aligned with our expectations that the positive effect of snow cover on population growth rate would be diminished at high moose densities. However, the interaction term was extremely weak and the BCI broadly spanned zero (b days snow 9 N Moose = À0.003, 95% BCI = À0.055, 0.030, f = 0.55). We therefore conclude that the effect of ticks operates largely in a density-independent manner in this population.
The models that tested other hypotheses (i.e., winter severity and heat stress) received less support, as did a model combining the effects of summer rainfall and winter severity. The parameter estimates for snow depth in both models in which it occurred were negative, which indicates a deleterious effect; however, the models still ranked low in the candidate set. We assumed that greater snow depths would negatively influence moose due to increased energetic costs from locomotion, the inability to efficiently access forage, and density-dependent competition for forage due to home range compression as moose sought out refugia with less snow. Although other moose studies have reported such effects (Thompson 1980, Mech et al. 1987, Crête and Courtois 1997, Solberg et al. 1999, Keech et al. 2000, Murray et al. 2006, those concerned higher latitude populations which plausibly experience more severe winters than our study area at the southern extent of the species' distribution. By contrast, we hypothesized that the low latitude position of this population would render moose more vulnerable to heat stress and potentially causing them to alter behavior to seek thermal refugia, reduce feeding rates, and experience amplified physiological processes such as respiration and heart rates (van Beest and Milner 2013, Street et al. 2015, Ditmer et al. 2017. We speculated that the physiological stress induced by warm temperatures should be most evident at populations near the species' equatorial range limit and would therefore be reflected in the dynamics of this population. This was not the case. Not only did this model receive the least amount of support, but the positive sign of the coefficients for both summer and winter temperatures was opposite of what we expected. One potential reason for this finding stems from the lack of detailed historic temperature data in the study area to evaluate the specific temperature thresholds determined to induce heat stress in moose (Renecker andHudson 1986, McCann et al. 2013). Because such specific climate data were unavailable for our entire study period, we used coarser climate summary data that did not explicitly evaluate heat stress thresholds specific to moose. Another possible explanation for the positive effects of temperature we documented is that the moose population generally became larger coincident with temperatures that generally became warmer during the course of the study.
Our results suggest that climate is indeed influencing this moose population and that both summer and winter precipitation are important yet operate through different pathways. Summer rainfall appears to positively influence moose across their geographic range, and this study is no exception. However, the effects of winter precipitation are more varied across the species' range. Consistent with predictions from biogeographic models of species range limits, we found evidence that external factors (i.e., parasites, mediated by climate) affected this peripheral moose population (Caughley et al. 1988). Current climate models predict that cold-season precipitation in this region will increasingly fall as rain, not snow (Safeeq et al. 2016), and this suggests that future conditions may exacerbate the effects of winter ticks on moose. Similarly, summer rainfall in the western United States is projected to decrease in the future (Easterling et al. 2017), and this too could negatively affect moose. Overall, our study illustrates how predictions about the effects of environmental change on wildlife population dynamics may require greater attention to multiple pathways of climate including insight about hostparasite interactions.

ACKNOWLEDGMENTS
We thank the dedicated biologists and pilots who collected the aerial survey data. Funding for the study came from the Utah Division of Wildlife Resources, the Albert W. Franzmann and Distinguished Colleagues Memorial Award, Sportsmen for Fish and Wildlife, Safari Club International, Utah Archery Association, and Rocky Mountain Elk Foundation.