Behavioral changes and nutritional consequences to elk ( Cervus canadensis ) avoiding perceived risk from human hunters

. The life-and-death stakes of predator – prey encounters justify the high price of many antipredator behaviors. In adopting these behaviors, prey incur substantial non-consumptive costs that can have population-level consequences. Because prey knowledge of risk is imperfect, individuals may even adopt these costly behaviors in the absence of a real threat. For example, rather than only avoid hunters, many species categorically avoid all anthropogenic activity. Although hunting seasons only increase risk for speci ﬁ c individuals (e.g., males), non-target individuals may still perceive human hunters as a source of risk and respond accordingly. Here, we used a large-scale experiment including 89 animal-years of location data from 62 unique individuals over 6 yr to quantify the duration, magnitude, and energetic consequences of changes to movement and resource selection behavior adopted by non-target female elk ( Cervus canadensis ) in response to human hunters during three separate experimental 5-d hunts (elk archery, deer ri ﬂ e ( Odocoileus hemionus or Odocoileus virginianus ), and elk ri ﬂ e). We predicted that elk response to hunters would be brief, but that strong behavioral responses to hunters (e.g., strengthened avoidance of roads and trails) would carry nutritional costs. We measured the duration of hunt-related changes in elk speed using quantile regression, further quanti ﬁ ed the strength of elk behavioral responses to hunters using population-level resource selection functions, and evaluated whether anti-predator resource selection behavior translated to measurable metabolic costs in the form of reduced body fat head-ing into winter. Elk responses to human hunters were stronger in the day than at night and were generally more pronounced during the elk hunts than during deer hunts. During hunts, elk shifted their diurnal behavior to avoid forage and intensi ﬁ ed their avoidance of roads and trails. The combination of these changes in behavior led to a predicted pattern of distribution during the hunt that differed substantially from the distribution prior to the hunt. Lactating females that more strongly avoided roads entered winter in poorer nutritional condition, suggesting that the changes in resource selection we describe carry corresponding nutritional costs that have the potential to impact subsequent population performance.


INTRODUCTION
Nearly all species face the risk of predation, and this omnipresent possibility of death has driven the evolution of a wide variety of costly antipredator behaviors (Lima and Dill 1990). Employing these behaviors generally entails energetic or opportunity costs (or risk effects) that can constitute the single largest expense of foraging (Brown and Kotler 2004). These behaviors are not linked to objective threats but are instead elicited by the perception of risk. How prey respond to predators therefore offers key insights into how animals perceive risk. Even in the absence of real threats, perceived risk can reduce fitness (Zanette et al. 2011) and alter prey distributions and activity (Schmitz et al. 2004). Although the existence of risk effects is well established, there has been comparatively little work testing mechanistic explanations or the scales at which these effects are influential (Ferrari et al. 2009, Beauchamp and Ruxton 2011, Moll et al. 2017. Consequently, we still know little about how most species perceive risk, leaving uncertainty around when we should expect risk effects to influence population performance. The population-level consequences of risk effects can be considered as the product of the duration of anti-predator responses, the magnitude of costs associated with these responses, and the distribution of these responses within the affected population. Crucially, individual responses to predation risk are constrained by the need to meet nutritional thresholds associated with survival and reproduction (Lima and Dill 1990). To this premise, the risk allocation hypothesis (Lima and Bednekoff 1999) adds the (implicit) assumption that prey are familiar with temporal patterns of risk (Ferrari et al. 2009) and makes two predictions. First, because it is costly, anti-predator behavior should be of short duration and closely linked to temporal spikes in risk. Thus, like predation risk itself, prey responses to risk are expected to be temporally dynamic (Lima and Bednekoff 1999). Temporal shifts in prey behavior, for example, are common in systems where predator diel schedules are fixed (Schmitz et al. 2017). Second, the risk allocation hypothesis predicts that individuals in better nutritional condition can afford to be more risk averse, spending energetic reserves or incurring opportunity costs to avoid perceived threats. This creates opportunity for a range of risk avoidance behaviors with varying cost. Furthermore, the distribution of these behaviors is expected to depend on the population's nutritional health. Testing these predictions therefore requires detailed data on nutritional condition, which is seldom available. Unraveling the price of perceived risk depends on the difficult dual task of accounting for temporal patterns of risk and individual variation in behavioral response Bednekoff 1999, Moll et al. 2017).
Prey perception of risk may also be influenced by predator characteristics and landscape features. Schmitz (2005) proposed that sit-andwait and sit-and-pursue predators should generate stronger prey responses because their relatively sedentary strategies mean that cues to their presence should be strong indicators of risk compared to coursing predators, whose presence should be more weakly linked to cues. Alternatively, recent game theory suggests that prey responses should be predicted by predator lethality (i.e., the probability of death given detection; Gal and Casas 2014). The difficulty of separating lethal and risk effects has largely limited tests of these hypotheses to small scales more amenable to experimental manipulation (Schmitz 2005). Practical and ethical concerns often preclude the manipulation of wild predator populations, restricting opportunities to test these hypotheses across larger spatial scales.
The idea that risk effects vary spatially (Berger 2007) is central to the management of many hunted populations. Elk (Cervus canadensis) management, for example, has focused on security areas most often identified on the basis of distance from roads and canopy cover, which, respectively, influence hunter access and elk hiding opportunities (Hillis et al. 1991, Proffitt et al. 2013. Previous attempts to test elk reliance on security areas have been hampered by elk spatial avoidance of human hunters through movement onto private lands or other anthropogenic refuges (Conner et al. 2001, Vieira et al. 2003, Proffitt et al. 2013, Ranglack et al. 2017. Identifying the roles that hunting mode and landscape features play in regulating risk effects is central to anticipating their consequences and understanding how prey responses to predators differ from other types of disturbance (Frid and Dill 2002).
Here, we used a large-scale experiment (sensu Walters and Holling 1990) to quantify the duration, magnitude, and energetic consequences of changes to movement and resource selection behavior adopted by non-target female elk in response to human hunters. For six years, we instituted separate five-day hunting seasons for elk archery, deer rifle (Odocoileus hemionus or Odocoileus virginianus), and elk rifle hunts. Human hunters are ideal study predators because (1) their density is often known (i.e., they can be feasibly censused); (2) hunting seasons and regulations provide detailed knowledge of the temporal pattern of predator activity needed to assess the duration of prey responses; (3) hunters' effects on non-target individuals are exclusively non-consumptive (i.e., true risk to female elk in our study is zero), allowing isolation of risk effects; and (4) multiple hunts (e.g., archery and rifle) with different associated predator behaviors and successful harvest rates allow for tests of competing hypotheses regarding the influence of predator characteristics on risk effects (Schmitz et al. 2004, Moll et al. 2017. Additionally, elk within our study did not have access to private land, limiting strategies to spatially avoid risk to those based on associated landscape features. Finally, detailed information on nutritional condition and reproductive history from elk in our population allowed us to test whether differences in individual behavioral responses affected the energetic reserves with which elk entered winter. We predicted (1) following the risk allocation hypothesis, elk behavioral responses to human hunters would be brief, and closely associated with human activity both within the diel cycle and across hunt seasons; (2) the strength of nontarget elk response to hunters would increase with increasing lethality (i.e., the true risk posed to male elk); and (3) individuals that showed stronger behavioral responses to hunters would enter winter in poorer nutritional condition (see Appendix S1: Table S1 for more detailed predictions). To test these predictions, we used three complementary analyses of elk movement and resource selection. We applied quantile regression to measure the duration of hunt-related changes in the speed of elk movements and population-level resource selection functions (RSFs) to quantify the magnitude of other behavioral changes in response to hunters. Finally, we used linear regression to evaluate whether individuallevel risk aversion (as measured by individuallevel RSFs) translated to measurable metabolic costs in the form of reduced body fat entering winter.

Study area
In 1987, the U.S. Forest Service created a large enclosure to facilitate experimental research in Starkey Experimental Forest and Range (Starkey). The main study area is 7768 ha surrounded by ungulate-proof fencing (Fig. 1). The area is composed of dry forest and grassland vegetation typical of large areas of the interior western United States (Rowland et al. 1997).

Experimental treatments
From 2008 to 2013, three hunts were held annually at Starkey for male elk or deer, in chronological order: elk archery, deer rifle, and elk rifle (Fig. 1). The number of tags issued for each hunt was consistent across seasons and years (25 hunters per hunt period, i.e., 0.32 hunters per km 2 , and all hunts lasting five days). We defined a 20-d period surrounding each hunt, which we divided into four 5-d treatments: (1) a null pre-hunt treatment that preceded the arrival of hunters; (2) a scout treatment during which hunters were present and allowed to move freely throughout the study area, but could not yet harvest animals; (3) a hunt treatment (the hunting season) when hunters were allowed to harvest animals; and (4) a post-hunt treatment in which hunters were absent or inactive (having departed prior to or immediately following the hunt's conclusion; Fig. 1). The temporal distribution of hunter activity varied by year and hunt, with hunters arriving and departing on different days. Nearly all hunters arrived during the scout treatment and stayed through the majority of the hunt treatment (M. Rowland, unpublished data). We therefore interpret the scout and hunt treatments as representing both an increase in intensity and a change in the type of human activity (e.g., hunter density was higher during the hunt treatment compared with the scout treatment). Pre-hunt and post-hunt treatments represented baseline levels of anthropogenic disturbance (i.e., human ❖ www.esajournals.org activity unrelated to hunting). The elk archery post-hunt and the deer rifle pre-hunt treatments were separated by 15 d, while the deer rifle posthunt and elk rifle pre-hunt treatments were separated by 5 d (exact dates in Appendix S1: Table S2). Hunters could legally harvest males from 30 min before sunrise to 30 min after sunset but may have also been active outside of legal Pre-hunt roads were closed when hunters were present. Distance to roads and canopy cover have been used in traditional definitions of "security areas" that provide elk (Cervus canadensis) refuge from human hunters. Our study focuses on elk behavior across three hunt periods that each consisted of four 5-d treatments (B). During these treatments, hunters were, respectively: absent (pre); present, but not hunting (scout); actively hunting (hunt); and absent or inactive (post).
❖ www.esajournals.org 4 September 2019 ❖ Volume 10(9) ❖ Article e02864 harvest hours as they moved to or from preferred hunting locations. During rifle hunts, harvest was prohibited within a 402-m buffer surrounding the perimeter fence; during the archery hunts, this buffer was reduced to 46 m. Data on hunters in our system support differences in behavior between archery and rifle hunters, with rifle hunters exploiting areas in closer proximity to roads (M. Rowland, unpublished data). Based on these differences in behavior and distribution, we expect that rifle hunters provide prey more reliable cues to their presence (e.g., gunshots; similar to sit-and-wait predators) compared to archery hunters, whose broader patterns of space use and hunting style mean their presence should be more weakly linked to cues (similar to coursing predators; Schmitz et al. 2004). Consequently, if elk anti-predator behavior is dictated by the reliability of cues to predator presence, we would predict that elk show similar responses to rifle hunters (regardless of target species) and a comparatively weaker response to elk archery hunters (Appendix S1: Table S1).

Treatment response
Hunts for male elk and deer allowed us to evaluate the non-consumptive effects of predation on female elk. Female elk at Starkey were hunted in the past, but not during the years described here (last antler-less hunt: 2006; Davidson et al. 2012). To assess the response of female elk to each treatment and hunt period, we used telemetry locations and animal condition data collected on adult females 2-14 yr old. These elk were captured and equipped with global positioning system (GPS) collars (in accordance with IACUC 92-F-0004; Wisdom et al. 1993). All collars were programmed to collect locations at a 30-min or higher frequency during the hunting period. To simplify comparison of movement among animals, we subset all individuals to a 30min fix frequency. We then excluded data from collars with an acquisition rate lower than 0.85 to limit habitat-induced bias (resulting in n = 89 animal-years of location data from 62 unique individuals; Hebblewhite et al. 2007). Data on animal condition (e.g., body fat, lactation status) of these telemetered females were collected each winter (11 December-15 January) when animals moved to the winter area and could be captured and handled by biologists (Rowland et al. 1997, Cook et al. 2010. We used telemetry location data to evaluate female responses to hunting in terms of duration (change in movement rates) and magnitude (change in resource selection). We further used the telemetry data in combination with animal condition data to estimate the non-consumptive costs of elk behavioral responses to elk rifle hunters.

Duration
To quantify the phenology of changes in elk behavior associated with exposure to human hunters, we modeled the speed of elk movements (m/30 min) using quantile regression (Cribari-Neto and Zeileis 2010). We chose movement rates as our first focal metric because in ungulates it is more sensitive than resource selection as an indicator of behavioral change in response to hunting pressure (i.e., animals may be displaced without changing their resource selection behavior, but animals that change their resource selection must necessarily alter their spatial distribution; Brown et al. 2018). This approach allowed us to quantify changes in the magnitude and variation of elk movements across hunt periods and treatments. For each hunt period, we fit five models, estimating the 0.10, 0.25, 0.50, 0.75, and 0.90 quantiles, and quantified differences between the 5 d pre-hunt treatment and each subsequent 5-d treatment (i.e., scout, hunt, and post-hunt). Because recent research indicates strong diel temporal structure in elk behavior (Kohl et al. 2018, Spitz et al. 2018, all models were fit separately to nocturnal (one hour after sunset to one hour before sunrise) and diurnal (all remaining) locations. We chose these divisions so that diurnal data would include a 30-min buffer around legal hunting hours when hunters are likely actively moving on the landscape and may influence elk behavior. We interpret increases in speed as having two potential fitness consequences for elk: decreased foraging opportunity (opportunity cost) and increased metabolic demands (energetic cost; Lima and Bednekoff 1999).

Magnitude
To quantify hunters' effect on female elk behavior, we fit a series of population-level RSFs. As above, we fit separate models with diurnal and nocturnal data for each of the three hunt periods ❖ www.esajournals.org (resulting in six models total). Because we were primarily interested in quantifying the difference in behavior during hunt treatments compared to the pre-hunt baseline, in fitting these models we excluded data from other treatments (i.e., scout and post) to limit model complexity. We modeled each hunt period separately. Each model contained the same set of covariates (described in next paragraph) as well as an interaction between each covariate and a hunt indicator variable (0/1). Thus, these interactions quantified changes in elk resource selection behavior during the hunt treatment relative to the pre-hunt baseline. Because elk behavior across hunts may also change with plant phenology, we chose to define the pre-hunt baseline separately for each hunt period (using the corresponding pre-hunt treatment rather than a single consistent baseline across all models). This minimized the temporal separation of baseline and treatment (to ≤14 d), limiting potential for shifting phenology to confound the interpretation of our results. We structured the resulting models as where w is the relative probability of selection, x 1 -x i are covariates, b 1 -b i are the coefficients associated with these covariates, hunt is an indicator variable delineating the hunt treatment, and the h subscript indicates coefficients quantifying the magnitude of the behavioral change associated with the hunt treatment. We defined the whole of Starkey's main study area as available so that our models quantified which regions had higher relatively probability of selection within this broader area (second-order selection; Johnson 1980), and we used a 30-9 30-m pixel as our sampling unit. We divided our model terms into two groups, mechanistic variables, on which we focus for inference, and nuisance variables, which may be important to accounting for elk behavior in our system but resist direct interpretation. Mechanistic variables included canopy cover (0-100%); greenness (amplitude of the normalized difference vegetation index, i.e., the amount of photosynthetic activity relative to baseline; 1-100); distance to roads open to motorized use during the hunt (km); distance to all-terrain vehicle (ATV) trails (km); slope (°); and topographic positioning index (TPI, representing concavity, values <0 denote valleys, values >0 denote ridges). We interpret canopy cover as an index of hiding cover (Hillis et al. 1991), greenness as an index of available forage (Ranglack et al. 2017), and distance to roads and trails as an index of hunter access (Proffitt et al. 2013). Slope and TPI may also influence hunter access and hiding cover. The effect of roads and ATV trails was limited to 1 km based on prior estimates of elk's maximum flight distance relative to ATVs in our system (Preisler et al. 2013).
Nuisance variables included distance to roads only open to motor vehicles during the pre-hunt treatment (km; pre-hunt roads), distance to perimeter fence (km), and elevation (m). We considered pre-hunt roads a non-mechanistic parameter because their closure during the scout and hunt treatments prevented us from distinguishing the effects of their closure from the concomitant change in hunter activity. Similarly, we included distance to perimeter fence to account for possible edge effects resulting from our study system being bounded. We deemphasize the interpretation of this term because we were unable to distinguish between two divergent interpretations, which represent different underlying mechanisms: the selection by elk of areas near the fence that are subject to the no hunting buffer (proactive, suggesting partial refuge) or animals fleeing hunters being forced to alter course when pushed toward fence lines (reactive, no benefit of refuge). Elevation is difficult to interpret because it varies so little within our study area (<380 m), but likely captures information about local vegetative communities.
We tested covariates for correlation (|r| > 0.6, finding none) and centered and standardized continuous covariates by subtracting the mean and dividing by two standard deviations (Gelman 2008). Model parameters were estimated using logistic regression (Johnson et al. 2006). To evaluate model performance, we used k-fold cross-validation using 10-fold and~20% of our data (by individual), which we had withheld from model training (Boyce et al. 2002). We used Wald statistics (i.e., z = Est./SE; Hosmer et al. 2013) to evaluate the importance of individual covariates. A Wald statistics' magnitude represents a covariate's relative explanatory importance, and the statistic's sign indicates the relationship's direction.

Variation in behavioral cost
To test for metabolic costs of anti-predator resource selection behavior, we fit individuallevel RSFs, extracted coefficient estimates from these models, and then used linear regression to model early-winter body condition (i.e., body fat) as a function of individual resource selection behavior, life history, and climate. Previous work on fat cycles in elk indicates that reproductive efforts can have carry-over effects on body condition in subsequent years (Cook et al. 2013). We therefore restricted our analysis to individuals with GPS location data bounded by animal handling events in two consecutive winters (i.e., where body condition and reproductive status from the preceding winter were known). Because rates of pregnancy in our population are high and constant, we further restricted our analysis to individuals with similar reproductive histories by excluding animal-years of data that began with an animal not being pregnant.
For the remaining data, we modeled diurnal resource selection during the elk rifle hunt separately by individual, and year: where w jk is the relative probability of selection by animal j during year k, x 1 -x i are resource unit characteristics, b 1 -b i are the coefficients associated with these characteristics, and each unit represents a 30-9 30-m pixel. In contrast to our population-level models, each individual model excluded data from the pre-hunt treatment and thus did not contain interaction terms quantifying changes related to the hunt. Instead, these models only quantified an individual's diurnal resource selection behavior during the elk rifle hunt. We focused on this temporal period to extract an index of elk risk avoidance because we expect the largest risk response during the diurnal elk rifle hunt, and because this hunt is closest in time to the subsequent evaluation of elk body condition. Availability for individual models was defined using a match-case control design with one available location for each used location generated using the population distribution of steplengths and turn angles (Hosmer et al. 2013). These models were otherwise identical in structure to our populationlevel models (including identical covariates). We checked all covariates for correlation and centered and standardized each covariate by subtracting the population mean and dividing by two times the population standard deviation. We used the coefficients from these individuallevel models of resource selection as predictors of an individual's maximum rump fat thickness (MAXFAT; Cook et al. 2010) entering winter using linear mixed models with random effects to account for repeated measures. Maximum rump fat thickness is a useful index for total body fat because the two are linearly correlated across a wide range of body fat levels (Cook et al. 2010). Because we were specifically interested in how elk altered their behavior relative to security areas, we focused on two focal covariates: canopy cover and distance to roads open to motorized use. We interpret these covariates as respectively representing hiding cover and hunter access. The extremity of some beta estimates for distance to roads also led us to exclude from our model any value where | b| > 10 (all b values for canopy cover and greenness were between AE8). We also controlled for life history and climate by including terms for individual reproductive history (lactation status, age) and weather (August precipitation, cm; a key measure of seasonal forage quality; Johnson et al. 2013). We chose lactation status because for ungulates lactation is more energetically expensive than gestation (Clutton-Brock et al. 1989) and lactation in fall therefore indicates that an individual has borne the cost of lactation throughout the preceding summer. Because the effect we observed on body condition was large relative to the duration of the elk rifle hunt treatment, we performed a post hoc test for consistency in individual resource selection behavior across all three hunts. We did this by fitting individual-level resource selection models to the elk archery and deer rifle hunts using the same methods as above and then comparing individuallevel coefficients for distance to roads across hunts using Pearson's correlation coefficient. Our inference was based on the relative support for parameters within each model, which we evaluated using associated P-values. All analyses were conducted in R (R Core Team 2014).

Duration
Speed of female elk (m/30 min) varied by hunt period, treatment period, and diel intervals ❖ www.esajournals.org (Fig. 2). Diurnal scout treatments were characterized by an increase in the magnitude of extreme values (e.g., 0.90 quantile) which continued and peaked during the hunt period. Speeds during the post-hunt period generally approximated pre-hunt levels for all hunt and diel periods. Nocturnal speeds were much lower than diurnal speeds and remained consistent across treatments during the elk archery and deer seasons. During the elk rifle hunt, nocturnal speed did, however, still show a peaked response during the hunt treatment similar to diurnal intervals. Middle quantiles remained relatively flat throughout (Fig. 2). Overall, the hunt treatments for elk archery, deer rifle, and elk rifle resulted in a mean increase in speed of 21.3, 11.7, and 31.8 m/30 min, respectively (mean pre-hunt speed = 84.8 m/30 min), but the fastest diurnal movements in the 0.90 quantile increased substantially, by 64.4, 81.6, and 133.3 m/30 min, respectively.

Magnitude
Models of female elk resource selection generally performed with high precision (Appendix S1: Table S3), with lower precision for nocturnal elk archery and deer rifle models. Female elk resource selection behavior and the explanatory importance of covariates varied by time of day (day/ night), hunt period (elk archery/deer rifle/elk rifle), and treatment (pre-hunt/hunt; Table 1, Fig. 3). Greenness, canopy cover, and distance to roads and ATV trails were important to elk behavior, and elk displayed distinct relationships to each. Elk tended to select greenness across all periods, but it was more important during Fig. 2. Predicted speed quantiles and 95% confidence intervals of non-target female elk (Cervus canadensis) surrounding three hunt periods (elk archery, deer rifle-Odocoileus hemionus, and elk rifle) in northeast Oregon, USA. Predictions were divided by diel period (diurnal and nocturnal), and for each hunt period, elk behavior was compared across four 5-d treatments during which hunters were as follows: absent (pre); present, but not hunting (scout); actively hunting (hunt); and absent or inactive (post). The largest changes to elk speed were seen in the most extreme movements (i.e., highest quantiles; note that y-axes vary across plots). Responses to hunting typically varied across the diel cycle (diurnal > nocturnal), across treatments (hunt > scout > pre~post), and across hunt periods (elk rifle > elk archery > deer rifle).
pre-hunt vs. hunt treatments (Table 1). During the day, however, elk avoided greenness during both rifle hunt treatments as indicated by the change in slope from positive to negative (Fig. 3). During the day, elk also generally selected for higher canopy cover but avoided denser canopy at night during both rifle hunt periods. Patterns of selection for canopy during the pre-hunt closely paralleled those for the hunt treatment. Elk tended to avoid roads during the day, but selected areas closer to roads during the pre-hunt treatments of both rifle hunt periods. At night, elk generally moved closer to roads, except during elk hunts, when animals appeared indifferent to roads (elk archery) or continued to avoid them (elk rifle). Elk consistently avoided ATV trails, and this avoidance intensified during the diurnal rifle hunts, relative to pre-hunt periods (Table 1; Fig. 3). There was no clear pattern of elk behavior relative to slope and TPI.
Among nuisance variables, elk generally avoided the fence, except during the hunt treatment of the elk rifle period when they switched to selecting areas closer to the fence (as indicated by diurnal and nocturnal Wald statistics changing from positive to negative; Table 1). During the day, elk also avoided roads only open to vehicles during the prehunt, but this avoidance was more important during pre-hunt than during hunt treatments. There was no clear pattern in nocturnal elk behavior relative to these roads nor relative to elevation.
These relationships combined to alter the predicted pattern of elk use between and within treatment periods (Fig. 4). Predicted nocturnal patterns of relative probability of use appeared increasingly diffuse in later hunt periods, and nocturnal predictions appeared generally similar between pre-hunt and hunt treatments. By contrast, diurnal spatial predictions were most concentrated during the last hunt period (elk rifle) and showed dramatic shifts between the prehunt and hunt treatments. During the diurnal pre-hunt treatments, elk had a higher relative probability of using areas adjacent to roads, but Notes: Effects are divided by diurnal/nocturnal (A/B), hunt period (horizontal divisions: elk archery, deer rifle, elk rifle), and treatment (horizontal subdivisions: pre-hunt/hunt; see text for details). Elk behavior and the explanatory importance of covariates varied by diel period, hunt period, and treatment. Nuisance variables (i.e., those relationships we do not expect to extend beyond our study area) are shown in italics. Units are omitted because all variables were standardized. TPI, topographic positioning index; ATV, all-terrain vehicle. Fig. 3. Predicted changes in female elk (Cervus canadensis) resource selection behavior between pre-hunt (yellow) and hunt (blue) treatments. Plotted relationships are from population-level resource selection functions divided by diel (vertical) and hunt (horizontal: elk archery, deer rifle, elk rifle) periods. Predictions are shown during the hunt treatments, elk were predicted to shift their use to the NE and south, moving away from roads and trails and toward the fence.

Variation in behavioral costs
We found moderate support for the hypothesis that female body fat at the onset of winter was partially determined by individual differences in resource selection strategies during the final hunt treatment (i.e., elk rifle; Appendix S1: Table S4). Animals that showed stronger relative avoidance of roads during the elk rifle hunt treatment entered winter with lower fat reserves (P = 0.062; Fig. 5). As expected, we found strong support for lactating individuals entering winter with lower energetic reserves (P < 10 À5 ) and moderate support that older animals also entered winter with lower fat levels (P = 0.074). We did not, however, find support for an effect of August precipitation on body fat, nor did we find a significant relationship to individual differences in elk selection for canopy cover during the elk rifle hunt treatment (Appendix S1: Table S4). Our post hoc analysis of consistency in elk resource selection behavior across hunts indicated that the strength of an individual's avoidance of roads was positively related between the elk rifle and preceding two hunts (Pearson correlation coefficients 0.48 with P value = 0.0015 and 0.43 with P value = 0.0047, compared to elk archery and deer rifle hunts, respectively), indicating that risk-averse elk during one hunt were more likely to also be risk averse during other hunts.

DISCUSSION
Our results are novel in demonstrating that activity by human hunters can impose large nonconsumptive costs on non-target animals. As predicted, non-target female elk showed acute but short-lived behavioral changes in the presence of human hunters. Rather than the strength of elk responses declining chronologically (null hypothesis) or declining in archery compared to rifle hunts (hunting mode hypothesis), elk showed a stronger response to elk hunters than to deer hunters, and rifle hunters generally eliciting a stronger response than archery hunters. This pattern of responses supports the hypothesis that the magnitude of non-consumptive effects increases with predator lethality (Appendix S1: Tables S1 and S5) and further suggests that female elk distinguish between predators with similar coarse-scale cues and behaviors (i.e., deer and elk rifle hunters). This relative unimportance of hunting method (archery vs. rifle) compared to hunted species indicates that elk may be relying on finer-scale cues of predation risk (e.g., hunter movement characteristics) and that our understanding of how elk perceive risk could therefore be further improved through a closer examination of hunter behavior.

Duration
The comparatively short length of the hunting seasons we analyzed exposed elk to pulsed periods of perceived risk. The risk allocation hypothesis predicts that the magnitude of nonconsumptive effects (per unit time) should decline as exposure to risk is prolonged (Lima and Bednekoff 1999). Thus, we might expect the response in our study to be stronger than that of hunted elk in other systems (where elk are typically hunted 60-80 d, 6-8 times the season lengths we observed; Conner et al. 2001, Davidson et al. 2012, Ranglack et al. 2017. Nonetheless, less acute effects have the potential to be equally impactful on subsequent population performance when compounded over a more protracted period. The intensity of non-consumptive effects likely also varies as a function of predator density. Johnson et al. (2005) were able to show that the magnitude of ungulate movement responses increased with hunter density, and Davidson et al. (2012) observed that rates of pregnancy in elk at Starkey declined as the density of hunters during elk archery season increased. Hunter density in our study (0.32 tags/km 2 ) was lower than that with 95% confidence intervals and were generated across the observed range of each focal variable while all other covariates were held constant at mean values. The dashed lines at y = 0 demarcate non-significant (i.e., neither positive nor negative) effects. Non-significant pre-hunt effects were omitted (e.g., diurnal elk rifle canopy); hunt effects were included regardless of significance (e.g., nocturnal elk archery distance to roads).  observed in these studies (mean low-density treatment 0.53 tags/km 2 ; mean high-density treatment: 1.09 tags/km 2 ; Davidson et al. 2012), but was consistent with the range of hunting pressure found in other elk systems (e.g., 0.12 tags/km 2 , Proffitt et al. 2013;0.41 tags/km 2 , Vieira et al. 2003). Understanding how risk effects scale with predator density and longer temporal overlap between predators and prey should be a primary objective of future research.

Magnitude
Greenness, roads, and canopy cover provided a strong mechanistic basis for understanding the behavioral strategies elk used to respond to human hunting. These relationships largely conformed to expectations established by previous work that evaluated elk resource selection (Proffitt et al. 2013, Ranglack et al. 2017. Like changes in movement rates, changes in elk resource selection in response to hunting were largely diurnal, with some effects carrying over into nocturnal behavior. During hunts, elk generally displayed diurnal avoidance of greenness and roads, the reverse of their behavior during pre-hunt treatments. At night, elk continued to show stronger relative selection for greener areas, but this selection weakened during the hunt periods, paralleling the diurnal decline in selection for greenness. This suggests that elk responded to perceived risk by distancing themselves from roads and curtailing their use of nutritionally rich areas, similar to results observed by Ranglack et al. (2017). In contrast, diurnal selection for higher canopy cover was more consistent through time, indicating that elk selection for canopy cover may represent a baseline diurnal behavior rather than a reactive strategy of predator avoidance. Elk consistently avoided ATV trails during all treatments, and consistent with past research on anthropogenic disturbance (Preisler et al. 2013), this avoidance intensified during the rifle hunts when ATV use increased. The resource selection patterns of seeking areas of denser canopy, farther from motorized routes, and lower nutrition may serve to mediate energy expenditures associated with flight by decreasing an elk's probability of encountering human hunters. The opportunities forgone by spending time in areas with reduced forage quality (use of higher canopy and areas of lower greenness) may, however, limit individuals' ability to maintain or accrue fat, thus negatively impacting subsequent reproductive performance.
In part due to its small size and moderate canopy cover, our study area contains little to no elk security areas according to conventional published definitions (minimum of 1 km 2 of continuous forest >0.8 km from open roads, Hillis et al. 1991; ≥ 13% canopy cover ≥2.76 km from open roads, Ranglack et al. 2017). Nevertheless, our results affirm the importance of roads and vertical structure in defining areas perceived by elk as relatively secure. The spatiotemporal pattern of hunter avoidance we observed indicates that canopy cover is important to elk independent of hunter activity, but that avoidance of roads is integral to elk resource selection responses to human hunters. Thurfjell et al. (2017) reported that hunted female elk learned to successfully avoid human hunters, thus demonstrating the species capacity to perceive and successfully avoid risk. Elk can habituate to human activity (Berger 2007), but this is expected to be rare where hunters remain active (Frid and Dill 2002). The responses we observed further suggest that female elk overestimate the risk posed by human hunters, and thus incur non-consumptive costs associated with predator avoidance even when their true risk of predation is low or nonexistent.

Variation in behavioral cost
Despite the short duration of our hunt treatments, we found evidence that stronger avoidance of roads during the elk rifle hunt was associated with reduced energetic reserves with which individuals enter winter, presumably because stronger risk avoidance during this period is characteristic of broader patterns of risk avoidance by individual elk. This result is supported by our post hoc analysis and suggests that elk behavioral responses to hunters carry costs in energy expended and/or nutritional opportunities forgone. Ciuti et al. (2012) documented a spectrum of risk aversion in elk and stressed the importance of individual variation in baseline responses to understanding the population-level consequences of human disturbance. Unfortunately, this study was unable to control for the expected effect of individual body condition (i.e., the risk allocation hypothesis) and consequently the amount of this observed variation that can be attributed to personality remains unknown. We join these authors in emphasizing this as a compelling direction for future research. Cook et al. (2013) have shown that throughout the year elk employ a variety of compensatory mechanisms (including risk allocation), which lead to the convergence of energetic reserves among individuals with similar reproductive histories. Nonetheless, body condition entering winter is particularly important because it represents the end of the period in which individuals are able to accrue additional energetic reserves. The energetic reserves of temperate ungulates typically decline across winter even when individuals are fed ad libitum (Cook et al. 2004) and strong carryover effects are well documented wherein nutritional condition at the end of winter is often more strongly influenced by autumn nutritional condition than winter weather or winter range conditions (Cook et al. 2013). Thus, while earlier hunts leave individuals time to recoup energetic losses, late season hunts provide little opportunity for animals to build back fat stores before the onset of winter and may consequently have a greater impact on population performance. As more individuals face nutritional limitation, ungulate populations are expected to respond with a wellestablished pattern including sequential declines in juvenile survival, fecundity, and adult survival (Gaillard et al. 1998, Eberhardt 2002. Preisser and Bolnick (2008) argue that because risk effects can simultaneously operate through multiple pathways (e.g., by changing both foraging time and efficiency) studies that focus only on a single mechanistic path are likely to underestimate the cumulative magnitude of these effects. The changes in fat we observed are not, however, only a function of a single mechanism (e.g., opportunity costs attendant on changes in resource selection behavior), but represent a broader integration of these costs (also capturing, e.g., the energetic costs of related increases in movement). Thus, our results support Preisser and Bolnick's (2008) suggestion that risk effects may have a greater impact on performance than is commonly appreciated. Cook et al. (2004) suggest that reproduction may be affected when elk body fat falls below 9-13% ingesta-free body fat (IFBF). Given the range of behaviors we observed, our models predicted that the individuals in our study fall within this range (i.e., 7.6-9.7% IFBF, Cook et al. 2010) and that the individuals displaying the strongest responses to human hunters would enter winter with 2.1 points lower % IFBF compared to the least responsive individuals. This effect is nearly equivalent to our modeled prediction for the loss in body fat associated with the energetic costs of lactation (À2.2 points %IFBF). In studying hunted female elk, Davidson et al. (2012) showed that even early season hunts (elk archery) could have a negative effect on pregnancy rates in lactating elk, which are presumably limited by individual nutritional condition. Our results demonstrate that these negative effects on female elk nutrition and reproduction are possible even when only males are hunted and hunter density is relatively low.

Management implications
Unless the costs of seeking refuge can be mitigated, our results concur with other resource selection analyses (Ranglack et al. 2017) in indicating that elk will continue to redistribute from public to private lands or other areas with minimal human disturbance. Our research further suggests that managers can influence this trend in two primary ways. First, land managers can expand and enhance elk refuges from human hunters through road closures and by managing the landscape for a broader distribution of vertical structure that can provide hiding cover. A second complementary approach is to focus on vegetation manipulations to increase available nutrition (e.g., thinning of dense forest stands). Recent work in marine systems has shown that the non-consumptive effects prey incur in seeking refuge from predators decline with the foraging resources those refuges contain (Donelan et al. 2017). Together, these tools can help incentivize elk to stay on public lands by improving the quality of security areas, thus potentially increasing opportunity for hunters to pursue elk on public lands and reducing conflict between elk and private landowners. Although security areas for hunted species have typically been evaluated with respect to male vulnerability, risk effects on females during hunting seasons also have the potential to indirectly influence population productivity (Davidson et al. 2012). Following the risk allocation hypothesis, the demographic impacts of non-consumptive effects likely depend on a population's nutritional baseline. In populations known to be nutritionally stressed (Cook et al. 2013) or where juvenile recruitment is declining, managers should review the potential impact of non-consumptive effects on hunted females. In these cases, we encourage managers to consider road closures, either seasonally or permanently as well as the reduction of tags or season lengths to minimize the duration of human disturbance on the landscape. Sustaining robust elk populations on public lands requires that we further develop frameworks capable of integrating our understanding of elk behavior and spatial variation in nutrition and perceived risk (Gallagher et al. 2017).