https://cpw.cvlcollections.org/files/original/26518640dc85dec1efd9a23bfd21034a.pdf 7a9a762b1d7227a8f085fdf01f729ff3 PDF Text Text The research in this publication was partially or fully funded by Colorado Parks and Wildlife. Dan Prenzlow, Director, Colorado Parks and Wildlife • Parks and Wildlife Commission: Marvin McDaniel, Chair • Carrie Besnette Hauser, Vice-Chair Marie Haskett, Secretary • Taishya Adams • Betsy Blecha • Charles Garcia • Dallas May • Duke Phillips, IV • Luke B. Schafer • James Jay Tutchton • Eden Vardy �Buderman et al. Movement Ecology (2018) 6:22 https://doi.org/10.1186/s40462-018-0140-6 RESEARCH Open Access Time-varying predatory behavior is primary predictor of fine-scale movement of wildland-urban cougars Frances E. Buderman1* , Mevin B Hooten2, Mathew W Alldredge3, Ephraim M Hanks4 and Jacob S Ivan3 Abstract Background: While many species have suffered from the detrimental impacts of increasing human population growth, some species, such as cougars (Puma concolor), have been observed using human-modified landscapes. However, human-modified habitat can be a source of both increased risk and increased food availability, particularly for large carnivores. Assessing preferential use of the landscape is important for managing wildlife and can be particularly useful in transitional habitats, such as at the wildland-urban interface. Preferential use is often evaluated using resource selection functions (RSFs), which are focused on quantifying habitat preference using either a temporally static framework or researcher-defined temporal delineations. Many applications of RSFs do not incorporate time-varying landscape availability or temporally-varying behavior, which may mask conflict and avoidance behavior. Methods: Contemporary approaches to incorporate landscape availability into the assessment of habitat selection include spatio-temporal point process models, step selection functions, and continuous-time Markov chain (CTMC) models; in contrast with the other methods, the CTMC model allows for explicit inference on animal movement in continuous-time. We used a hierarchical version of the CTMC framework to model speed and directionality of fine-scale movement by a population of cougars inhabiting the Front Range of Colorado, U.S.A., an area exhibiting rapid population growth and increased recreational use, as a function of individual variation and time-varying responses to landscape covariates. Results: We found evidence for individual- and daily temporal-variability in cougar response to landscape characteristics. Distance to nearest kill site emerged as the most important driver of movement at a population-level. We also detected seasonal differences in average response to elevation, heat loading, and distance to roads. Motility was also a function of amount of development, with cougars moving faster in developed areas than in undeveloped areas. Conclusions: The time-varying framework allowed us to detect temporal variability that would be masked in a generalized linear model, and improved the within-sample predictive ability of the model. The high degree of individual variation suggests that, if agencies want to minimize human-wildlife conflict management options should be varied and flexible. However, due to the effect of recursive behavior on cougar movement, likely related to the location and timing of potential kill-sites, kill-site identification tools may be useful for identifying areas of potential conflict. Keywords: Animal movement, Hierarchical model, Individual variation, Population-level, Predation, Telemetry, Wildlandurban interface * Correspondence: franny.buderman@colostate.edu 1 Colorado State University, Departments of Fish, Wildlife, and Conservation Biology, 1484 Campus Delivery, Fort Collins, CO 80523, USA Full list of author information is available at the end of the article © The Author(s). 2018 Open Access This article is distributed under the terms of the Creative Commons Attribution 4.0 International License (http://creativecommons.org/licenses/by/4.0/), which permits unrestricted use, distribution, and reproduction in any medium, provided you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made. The Creative Commons Public Domain Dedication waiver (http://creativecommons.org/publicdomain/zero/1.0/) applies to the data made available in this article, unless otherwise stated. �Buderman et al. Movement Ecology (2018) 6:22 Background Individual-level movement decisions are one of the underlying processes that give rise to population-level patterns such as species distributions or their density and abundance on the landscape [1]. Movement decisions are a function of a number of variables, including the current location of the individual and the alternative available landscape [1]. Therefore, a central theme of animal ecology is the assessment of an individual’s selection for habitat, given what is available [2]. Habitat selection is typically characterized using resource selection functions (RSF), which are often fit using logistic regression to compare the locations used by an individual or population to a random sample taken across some area defined as “available” [3]. Use that is disproportionate to habitat availability implies that the individual selects for, or avoids, the given habitat [3]. However, inference on selection depends on what components are considered available to the animal [2]. For example, an animal may use a resource disproportionately less than is available in its home range, however it may have chosen its home range because the resource was abundant [2]. In addition, availability is constrained by an individual’s range of movement. To account for dynamic availability, spatio-temporal point process models simultaneously estimate the resource selection function and time-varying availability kernels, which is the area an individual is capable of moving to over a given period of time [4–6]. The more commonly used method, a step selection function, approximates the availability kernel by using conditional logistic regression and a sample of “available” steps that an individual could have taken (e.g., [7–9]). Recent methods have used conditional logistic regression to separately approximate the movement and time-varying availability kernels, in the vein of spatio-temporal point process models (e.g., [10]). However, because all of these methods are formulated in discrete time, inference is made only when data were observed and not on the unobserved path. In addition, aside from the spatio-temporal point process of [6], none of these methods account for measurement error in the observed locations. In contrast to many resource selection studies, one of the primary goals of continuous-time movement models is to estimate the true path of an individual when it was unobserved [6, 11–14]. Continuous-time movement models can also incorporate measurement error and irregular observations in time. However, movement models are typically time consuming and computationally intensive to fit, making it difficult to obtain inference on multiple individuals [15]. If inference on multiple individuals is attainable, it may be possible to identify a population-level response that is consistent across individuals, which would provide a rigorous link between individual choices and population-level patterns [1]. In Page 2 of 16 addition, understanding individual variability may help identify individuals that associate more strongly with certain features of the landscape [16]. A recently developed method, continuous-time Markov chain (CTMC) modeling, incorporates an explicit movement model to obtain information on travel speeds and directionality. Travel speeds may provide indirect inference on resource selection [17] and avoid absolute statements about selection [2]. The CTMC method [18, 19] is fit in two stages, where the first stage uses a continuous-time movement model to obtain inference on where the individual was when it was unobserved and account for measurement error, while the second stage allows for evaluation of landscape drivers of animal movement. The second stage of the analysis uses a Poisson specification with an offset to model transition rates; therefore, statistical software based on a Poisson likelihood can implement the CTMC movement model [19]. The flexibility of the CTMC framework can account for time-varying responses to landscape drivers by allowing coefficients to vary temporally [19], and it can also be implemented in a Bayesian hierarchical framework, allowing for inference on individual- and population-level drivers. Previous applications of the CTMC framework focused on inference for single individuals and did not make inference across multiple individuals [19, 20]. Quantifying individual variability in habitat selection, while simultaneously estimating population-level patterns, can be important for management and conservation issues where resources are heterogeneous or cause points of conflict [21]. Some large carnivores, such as cougars (Puma concolor), have undergone recent range expansions into human-modified landscapes [22], but they rarely use the heavily modified landscapes in urban and suburban areas, instead relying on the rural and exurban areas at the wildland-urban interface [21, 23]. Along with increased risk from human interactions [23], human-modified landscapes may contain greater numbers of both primary (ungulates, e.g., [24]) and secondary (domestic animals, e.g., [25]) prey for large carnivores compared to adjacent wild-land areas. As early as 1998, the frequency of human-cougar interactions along portions of the Front Range, a mountain range extending north-south from Casper, Wyoming to Pueblo, Colorado, have increased due to encroaching residential development, increasing cougar populations, and increasing prey densities near human populations [26]. The Front Range Urban Corridor runs along the eastern edge of the Front Range, while the Front Range itself contains a matrix of towns and areas that are managed for recreational use by county, state, and federal agencies. Human-cougar interactions have remained high in recent years (Mat Alldredge, Colorado Parks and Wildlife, personal communication), and cougars have been observed using developed areas in the Front Range as a hunting �Buderman et al. Movement Ecology (2018) 6:22 ground [27, 28]. In addition, due to their desirable qualities, regions adjacent to protected areas have higher human population growth compared to growth in rural, non-protected areas [29], increasing the potential for human-wildlife conflict [30]. Given the increasing potential for human-wildlife conflict as development permeates rural and wildland areas along the Front Range and elsewhere in the West, we sought to extend previous work by explicitly modeling fine-scale cougar movement to identify key drivers of their behavior, and in doing so, better understand their use of the wildland-urban landscape in both space and time. Many cougar studies do not explicitly model movement, and instead focus on resource selection; the animal locations used for inference were sometimes obtained only during daylight (e.g., [31–33]), obtained during night and day but were treated equivalently (e.g., [34]), or obtained at unspecified times (e.g., [35, 36]). Inference on time-varying behavior has been limited to separate analyses on discretized temporal periods (e.g., [17, 22]). Some studies have also focused exclusively on kill site and hunting locations (e.g., [27]) or non-kill site locations (e.g., [17, 22]). We used the CTMC framework to model individualand population-level cougar responses to landscape features in continuous time, which allowed for direct inference on how behavior varied at a temporally fine scale, given what was available. In addition, by using a hierarchical modeling framework, we accounted for individual-level variation, which may be a function of the spatial distribution of prey items or the behavioral flexibility of a generalist predator Page 3 of 16 [21, 37], while still obtaining population-level inference across a suite of individuals. Methods Data collection and study area As part of an ongoing study by Colorado Parks and Wildlife (CPW), cougars were trapped and fit with global positioning system (GPS) collars and released along the Front Range of Colorado (CPW ACUC 01–2008; Fig. 1). We focused on 19 adult individuals (M = 5, F = 14) that were monitored April 1–15 2011, 21 adult individuals (M = 7, F = 14) that were monitored during June 16–30 2011, and 21 adult individuals (M = 3, F = 18) that were monitored October 1–15 2011. The time periods of interest were chosen for two reasons: first, because observations were available for a large number of individuals, which is critical for making population-level inference across individuals, and second, because we were interested in examining seasonal differences in cougar movement due to temporal variability in the landscape-level covariates. For example, we expected a strong response to prey-based covariates year round but with seasonal shifts to reflect seasonal changes in prey availability. In June, mule deer fawns are born, and form a primary prey source for cougars [38] and are at a disproportionately high risk for predation [39]. However, cougars have been observed relying on smaller prey items in April, potentially due to competition with other species, and by October fawn predation has decreased and their diet switches to deer and elk (Mat Alldredge, CPW, personal observation). Fig. 1 Map of Colorado counties, with the cougar movement study area plotted in gray (Fig. 1a). Elevation (m; Fig. 1b) and land classified as developed (dark gray is < 10 acres/unit; Fig. 1c) is shown for the study area and surrounding area �Buderman et al. Movement Ecology (2018) 6:22 Cougars were trapped using cage traps, hounds, and foothold snares, and the minimum difference between trapping and the analysis window was 12 days, with the average time from trapping to analysis being 114 days for April, 130 days for June, and 235 days for October. One individual (AF69) was darted and relocated mid-analysis (April), however we generated the CTMC data separately for each period before making inference on movement drivers across the observation period. All individuals were monitored with Vectronics collars (Vectronics GmbH, Berlin, Germany) programmed to obtain fixes every 3 hours. Our study area comprised a 2,700 km2 region in the Colorado Front Range to the north-west of Denver (Fig. 1). The study area consisted of a matrix of private (43%) and public (57%) land [27]. Private land included areas of rural, exurban, suburban development, and small towns. Public land was managed by federal, state, and municipal governments for recreational activities or as open space. Road density was, on average 1.6 km/km2, but ranged from no roads to 16.6 km/km2. Elevation ranged between 1,522 and 4,328 m, generally increasing east to west, with development decreasing along a similar gradient. Ponderosa pine, Douglas fir (Pseudotsuga menziesii), lodgepole pine (Pinus contorta), and spruce-fir (Picea engelmannii-Abies lasiocarpa) were the predominant non-agriculture vegetation types, in order of dominance from east to west. Cougar density was approximately 2.4 individuals per 100 km2 (Mat Alldredge, CPW, personal communication). Native herbivores (elk and deer) form the predominant component of cougar diet in the region, with a smaller contribution from domestic species (pets and livestock) and synanthrophic wildlife [40]. However, the proportion of domestic and synanthropic species in the diet varies with cougar location along the wildland-urban gradient [28]. Stage 1: Continuous-time Markov chain model We used a Bayesian hierarchical CTMC model to evaluate drivers of cougar movement; this model is an extension of the model proposed by [19] and allows for inference on movement rates and directional bias, as opposed to use, in continuous time and discrete space. The initial step in the CTMC framework is to estimate a continuous movement path from the observed data points. We used the functional movement model developed by [13] to account for measurement error and predict locations every 10 minutes for the selected 2 weeks of each month. The functional movement model can Page 4 of 16 be implemented using the fmove.bayes function in the ctmcmove package [41]. The smoothness of the imputed paths can be controlled using the arguments associated with the precision matrix and the shape and scale parameters of the inverse gamma prior for the partial sill parameter of the spline basis coefficients. The details of the functional movement model are beyond this manuscript, but for the purpose of reproducibility we note that we used prior knowledge to fix the standard deviation pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi of the measurement error to logð10=4Þ , and modeled the variance of the basis functions using a prior consisting of a lag two conditional autoregressive precision matrix with the partial sill modeled as an inverse gamma with shape and scale parameters equal to one. Additional details of the functional movement model can be found in [13, 14, 41, 42]. Although the CTMC method is computationally efficient in terms of speed, there are trade-offs between the duration of the time-period of interest, the temporal resolution of the path interpolation, and the spatial resolution of the rasters, which together can create extremely large and difficult to store data frames. In addition, when the latent variable formulation is used, the discrete cell sequence must be contiguous, meaning that the spatial and temporal resolutions must match (e.g., an individual cannot move further than one grid cell between time points). In our analysis, the time-period of interest (2 weeks), temporal resolution (10 min intervals), and spatial resolution (100-m by 100-m or 1 ha), were selected to focus on fine-scale movement decisions. Additionally, crossing a 100 m2 area encompasses a cougar’s ability to move over a ten-minute interval [17]. To account for uncertainty in the true movement path, a random subset of imputed paths from the posterior predictive distribution of the movement model were spatially discretized to a latent variable formulation with a cell size of 100-m by 100-m, which was the lowest resolution among the available covariates. The CTMC model consists of a product of two components: the time an individual spends in a grid cell and the direction that an individual moves when it leaves a grid cell. The time an individual spends in a grid cell (motility) is exponentially distributed, such that a large rate parameter corresponds to fast movement out of the cell. When an individual leaves a grid cell, the probability that they move to a particular neighboring grid cell (directionality) is the ratio between the movement rate into that cell and the sum of the movement rates into all neighboring grid cells. Therefore, higher �Buderman et al. Movement Ecology (2018) 6:22 proportional rates indicate directional bias in movement. Thus, movement rate parameters, which are a function of covariates (i.e., landscape variables that correspond to the position of the cell on the landscape), control both motility and directionality. Hanks et al [19] showed that the likelihood of the CTMC model (i.e., the product of the motility and directional components) for movement can be expressed as a Poisson GLM using a latent variable formulation. In the latent variable formulation, each transition corresponds to four data points (the four neighboring grid cells); the response variable is equal to one if the neighboring grid cell is the cell that the individual transitioned into and zero otherwise. Modeling the latent variables (zeros and ones) as Poisson random variables with an offset for the amount of time an individual spends in a grid cell results in a likelihood that is equivalent to the CTMC likelihood. This allows inference to be obtained using standard GLM software, and the R package ctmcmove facilitates creation of the CTMC latent variable formulation [41]. Full CTMC details are available in Additional file 1. Using multiple imputed paths accounts for the uncertainty in the true path of the individual and is a process version of multiple imputation [15, 18, 19, 43], a method frequently used for missing data [44]. Process imputation is more computationally efficient than using the entire posterior distribution, but still approximates the uncertainty associated with the unobserved path. We generated 30 imputations for each individual [19, 43], using 20 imputations to fit the models for individual-and population-level inference on transition rates and 10 imputations to calculate the posterior predictive score that was used to select regularization terms. Regularization shrinks the effect of unimportant covariates toward zero to prevent over-fitting and, in a Bayesian context, this is achieved by using an informative prior for the coefficients [45]. Stage 2: Poisson models for movement inference We used a hierarchical generalized linear model (H-GLM) for individual- and population-level inference on average cougar behavior, as measured by movement rates and directional bias, as a function of landscape features. Because cougars and humans are active at different times throughout the day, we proposed an additional model, a hierarchical generalized additive model (H-GAM), to account for individual- and population-level diel time-varying behavior. Covariates were centered and scaled to the individual, meaning that the coefficients are relative to the mean and standard deviation of the values that each individual encountered during a given two-week period. This is similar to the idea proposed by [2], where selection was Page 5 of 16 determined by comparing some measure of usage and availability of a landscape feature on an individual basis. The hierarchical component of the model allows individual-level responses to vary around a population-level mean response, where both the individual and population-level estimates are obtained simultaneously. In the CTMC framework, the response variables, zij, were a sequence of zeros and ones, where zij~Poisson(λij), for i = 1, …, T and j = 1, …, J, where T was the total number of cell transitions, and J was the number of individuals. Landscape covariates were incorporated using the log link function, such that logðλij Þ ¼ logðτ ij Þ þ x0ij β j . The residence times were represented by the constants τij, and the landscape variables by xij. The parameter βj was a vector of P individual-level coefficients that arose from the population-level distribution β j � N ðμβ ; Σβ Þ. The covariance matrix, Σβ ≡ σ 2β diagðϕÞ , where the vector ϕ scaled the value σ 2β to each coefficient. The vector of scaling parameters consisted of a one for p = 1 (ϕ1 = 1) and was modeled as logðϕ p Þ � N ð0; 0:04Þ for p = 2, …, P. The population-level distribution had a mean that was modeled with a multivariate normal distribution μβ � N ð0; σ 2μ IÞ , where I is the identity matrix. Both σ 2β and σ 2μ were used as regularization terms, where σ 2β was selected a priori and σ 2μ ¼ 0:1 , to shrink the coefficients toward zero; this prevented over-fitting and allows for correlated predictors [45]. The H-GAM was formulated as a varying coefficient model [46], where the response to covariates varied over space or time. By expanding the landscape covariates with a basis function [47], we created a new vector, vij, that was the Kronecker product of the P length vector of covariates, x0ij , and the Q length vector of the values of the basis at the time of transition i, w(i). For diel movement, we used cubic cyclic spline basis functions (w(i)), because they constrain the start and end points of the varying coefficients to be equal, which is an important property for time spans that are cyclic in nature. The GAM for hourly movement was similar to the GLM, except logðλij Þ ¼ logðτ ij Þ þ v 0ij α j , where αj was a vector of length PQ. Each parameter in αj was the collective effect of the basis function and the corresponding covariate at the time of transition i. Using the vector w(i), αj can be back-transformed to obtain the time-varying effect of the covariate. In the hierarchical framework α j � N ðμα ; Σα Þ, where Σα ≡ diagðσ 2α ϕÞ . The vector ϕ again reduced the number of parameters we need to select a priori by scaling the σ 2α term to each parameter, and μα � N ð0; σ 2μ IÞ. Both σ 2α and σ 2μ served as regularization terms, where σ 2α was selected a priori and σ 2μ ¼ 0:1. Finally, to assess whether males and females exhibited different amounts of temporal variation in their response �Buderman et al. Movement Ecology (2018) 6:22 to potential movement drivers, we fit the GLM and GAM models to males and females separately for each time period. This resulted in four models: 1.) a GLM fit to all individuals, 2.) a GAM fit to all individuals, 3.) a GLM fit to females and a GAM fit to males, 4.) and a GAM fit to females and a GLM fit to males. We calculated the posterior predictive score for each model (i.e., the sum of the posterior predictive score for the models fit to males and females separately) and compared the scores across models within each month. Models were fit using a Markov Chain Monte Carlo (MCMC) algorithm written in R [48]. We performed adaptive tuning over an initial 50,000 MCMC iterations. We used the selected tuning parameters as constants in the subsequent 50,000 iterations that were used to calculate the posterior predictive score for the a priori regularization parameter grid-search. The final models were fit using 100,000 MCMC iterations with a burn-in period of 10,000 iterations. Landscape covariates Each covariate can be included as either a motility or directional driver of movement in the CTMC model. Motility covariates are based on the value of the grid cell that the individual is in currently and control the absolute rate of movement; positive coefficients indicate faster movement with increasing values of the covariate (and slower movement with decreasing values), and negative coefficients correspond to faster movements with decreasing values of the covariate (and slower movement with increasing values). Directional covariates account for the correlation between movement and the gradient of a covariate and contribute to the probability that an individual moves toward a grid cell. The directional drivers were calculated such that a positive coefficient indicates that individuals move predominantly in the direction that the covariate decreases (decreasing distance, such that they orient toward a feature), whereas a negative coefficient indicates that individuals move in the direction that the covariate increases (increasing distance, orient away from a feature). All rasters were aggregated to a 100-m by 100-m resolution, which is within the distance that a cougar might typically move over a ten-minute interval [17]; individuals cannot skip grid cells (enter a non-neighboring cell), therefore the spatial resolution of the rasters should reflect our prior knowledge about movement speeds. We hypothesized that a number of landscape covariates may contribute to transition rates and directional bias of cougars: mule deer (Odocoileus hemionus) utilization (as a proxy for availability), distance to nearest potential kill site, distance to nearest structure, distance to nearest road, elevation, heat insolation load index, and topographic wetness. We also used an Page 6 of 16 autoregressive parameter to account for directional persistence, or an individual’s tendency to move in the direction in which it was already moving [19]. Prey availability is a driving factor in cougar habitat selection. For example, cougars in western Washington used areas where suspected prey availability was high, such as low-elevation, early successional forests, and areas near water [21], and [27] observed cougars foraging in areas with high mule deer utilization. We approximated prey availability using two covariates: annual mule deer utilization and nearest potential kill site. The model averaged prediction for mule deer utilization [49] approximates prey availability given a suite of landscape covariates. We hypothesized that cougars would move slower in areas with high values for mule deer utilization and orient toward areas of high mule deer use during crepuscular and nocturnal movements [21, 27, 50]. Blake et al [51] found that many of the landscape variables that contribute to the location of predation events were the same as those contributing to non-predation habitat use, which led them to determine that cougars spend the majority of their time moving across the landscape in hunting mode. Including the location of a potential kill site may act as a proxy for unmeasured landscape variables and non-mule deer prey presence. In addition, potential kill sites represent known spatially recursive behavior based on memory and perception of the landscape [52–54]. Memory and recursively used locations have been incorporated into resource-selection analyses using individual-level intensity distributions [55, 56] and model-based Dirichlet processes [57]. Potential kill sites were determined using a clustering algorithm on the GPS points, where a location was classified as a potential kill site if two or more GPS locations, occurring between the average time of sunset and sunrise for each two-week period, were found within 200 m of the site within a six-day period (modified from [50, 58]). We calculated the distance (m) to nearest potential kill site identified within the two-week period to account for dependence in the movement process due to the known temporary activity centers induced by the potential kill sites. Up to nine potential kill sites were identified for each individual during the observation window (two-weeks). Although the data were used to generate the clusters, the CTMC model is not assessing resource use, but is determining whether speed and directionality vary as a function of the location of the clusters (e.g., we would not detect a response if they were not correlated with variation in movement). We expected individuals to move faster as distance to potential kill site increased, because decreasing distance may correspond to an individual returning to a cached kill, and caches are more often located in areas of high vegetation cover [59]. �Buderman et al. Movement Ecology (2018) 6:22 We calculated distance to nearest structure (m) as the Euclidean distance to the nearest man-made roofed structure [49]. Distance to road was calculated using major roads data (i.e., a major highway primarily for through traffic usually on a continuous route and streets whose primary purpose is to serve the internal traffic movement within an area) obtained from Colorado Department of Transportation. Due to increased human activity around structures and roads, we expected cougars to move faster when closer to roofed structure and distance to nearest road [17, 33, 60]. However, females may respond less to structures and roads than males, given that there may be additional factors, such as food limitation and offspring, which drive them to tolerate human-modified landscapes [37, 61]. We also expected there to be high temporal variability in the response to structures, because individuals have been observed avoiding areas of anthropogenic activity less at night, while avoiding contiguous forest habitat less during the day [22]. We used a digital elevation model (Fig. 1) to characterize elevation. Blecha et al [27] found that cougars avoided foraging in higher elevations, but [37] observed cougars selecting for higher elevations in developed areas. We expected cougars to show high temporal variability in their directional response to elevation, with cougars moving toward lower elevations when they are hunting (main prey is concentrated in lower elevations) and toward increasing elevations at other times. We used a raster based on the continuous heat insolation load index ([62], modified from [63]), to measure the accumulation of solar radiation at that location over the course of a year (MJ/cm2/yr). Heat insolation is higher on south-facing slopes that are more xeric and open than north-facing slopes [64]. Cougars have been observed using less rugged terrain for travel [17], selecting for south-facing slopes containing shrubs [22], and avoiding foraging on north-facing slopes [28]. Therefore, we expected that cougars may orient toward areas of high heat insolation, but move quickly through them. The topographic wetness plus metric (TWI+) predicts soil moisture based on slope, as originally described by [65], and aspect, as modified by [66]. Because cougars have been observed selecting for and hunting in riparian areas [21, 33, 60, 61], we expected cougars to move slowly in areas of high topographic wetness and demonstrate temporal variability in their directional response (toward areas of increasing topographic wetness when hunting). We also analyzed a subset of individuals and the interaction between housing density and their response to deer utilization and distance to nearest kill site. Despite cougars demonstrating avoidance of high housing densities while foraging (locations preceding a successful kill Page 7 of 16 and following previous prey handling), kill sites were positively related to housing density [27]. In addition, the temporal variability in the response to anthropogenic structures that was observed by [22] was stronger for cougars in rural, rather than wilderness, areas. Therefore, these are the two variables that we expected to vary most with housing density due to the potential trade-offs between increased prey abundance but increased mortality risk. To determine the effect of housing density on the response of cougars to deer utilization and potential kill sites, we discretized the landscape into developed (< 10 acres/unit) and undeveloped areas (Fig. 1). Only 13, 15, and 17 individuals for April, June, and October, respectively, were used in the secondary analyses because the remaining individuals did not spend time in developed areas in the selected two-week periods. We were unable to evaluate an interaction in the H-GAM (time-varying) framework due to high variability in the percentage of locations for each individual that were classified as occurring in undeveloped areas. Results There was no detectable effect of many of the landscape covariates on average motility or directionality at a population-level (Fig. 2). However, distance to potential kill site emerged as the primary driver of both motility and directionality in the GLM framework (Figs. 2 and 3). As individuals increased their distance from a potential kill site, their transition rate increased (Fig. 3a). In addition, individuals oriented movement toward their potential kill site (Fig. 3b). We also detected significant directional persistence (95% CI for April: 1-1.15, June: 1.03-1.16, and October: 1.07-1.20), or residual autocorrelation, indicating that individuals tended to continue moving in the direction they had previously been moving, after accounting for landscape features. Of the remaining potential drivers of movement, the largest seasonal differences were observed in the effect of heat loading, elevation, distance to nearest roofed structure, and distance to nearest road; however, the 95% credible intervals consistently overlapped zero (Fig. 2). Based on the posterior mean, individuals were observed moving slower at higher elevations in April, but faster in June and October (Fig. 4a). In contrast, individuals moved faster than average in areas where heat loading was high in April, and to a lesser degree, June, but moved slower with higher heat loading in October (Fig. 4b). We detected a consistent seasonal effect of distance to structure, but the posterior mean was negative, meaning that individuals moved faster as distance to structure decreased (Fig. 4c). A similar pattern was detected with distance to roads, however the effect became positive in June, with individuals moving faster as distance to roads increased (Fig. 4d). In most cases, individual-level �Buderman et al. Movement Ecology (2018) 6:22 Page 8 of 16 Fig. 2 The mean and 95% credible intervals for the population-level mean effects of landscape covariates on movement rates (M) and directionality (D) of cougar movement in the Colorado Front Range for two-week periods in April, June, and October 2011 uncertainty tended to be high, with a few individuals showing statistically significant responses despite a non-significant population-level response (Figs. 3, 4). Our results suggest that distance to nearest potential kill site was also the predominant motility and directionality driver in the diel time-varying framework (H-GAM; Figs. 5 and 6). However, the strength of the motility response to distance to nearest potential kill site varied over time and with seasons. The strongest motility response occurred around dawn, decreased steadily during daylight hours, and then increased around dusk (Fig. 5). The magnitude of this variation was strongest in June and October, and weakest in April (Fig. 5). The strength of the directional bias toward potential kill sites also varied through time, but was consistent across seasons (Fig. 6). The evidence for an effect of potential kill sites on directionality suggests that individuals orient less toward their kill site during daylight hours, and may even orient away from their potential kill sites during late afternoon (Fig. 6). While the 95% credible intervals overlapped zero for much of the day, we detected modest temporal responses in both motility and directionality to elevation and distance to nearest structure (Fig. 7). The motility response to elevation varied seasonally, as in the GLM framework. The average negative response to elevation observed in April (Fig. 4a) was reflected in a negative response to elevation around dawn (individuals move slower as elevation increases), with a slightly positive response later in the day (Fig. 7a). We observed little time variation in June, but the pattern observed in October was the opposite of April, with individuals moving faster with increasing elevation around dawn, with a decreasing effect through the rest of the day (Fig. 7a). Individuals moved toward higher elevations mid-day and toward lower elevations at other times, a pattern that was consistent across seasons (Fig. 7b). The strongest negative effect of distance to structure on motility (individuals move faster as distance decreases) occurred around dawn and dusk for all seasons (Fig. 7c). The effect on directionality was less consistent, with orientation toward roofed structures just after dawn, followed by orientation away from structures, in April and June; this pattern shifted toward pre-dawn in October (Fig. 7d). In addition, we did not see evidence for an interaction between development and deer utilization, which remained a statistically insignificant driver of cougar movement rates and directionality in both the H-GLM and H-GAM models. The positive effect of distance to potential kill site on speed (faster as distance to kill site increases) and directional bias (more orientation toward the kill site) was consistent between developed and undeveloped areas (Fig. 8a). However, we detected a difference in average movement rate between the two areas, with individuals in each month moving faster in developed areas (Fig. 8b). Finally, the H-GAM for both sexes was the best model in terms of predictive performance across all months, whereas the GLM performed the worst (Fig. 9). The models that were a mixture of a GAM and GLM, varying by sex, were generally equivalent (Fig. 9). The largest difference between the two sex-varying models was observed in October, when the better of the two models �Buderman et al. Movement Ecology (2018) 6:22 (a) Page 9 of 16 (b) Fig. 3 Posterior means and 95% credible intervals for the individual- and population-level static effects of distance to nearest potential kill site on motility (Fig. 3a) and directionality (Fig. 3b) of cougar movement in the Colorado Front Range for two-week periods in April, June, and October 2011 included time variation for males and no time variation for females (Fig. 9). Discussion Fine-scale cougar movement The observed response to distance to nearest potential kill site over a short time period is potentially due to cougars returning to the carcasses (e.g., spatial memory; [52–54]) and unmeasured fine-scale covariates related to landscape features that increase the likelihood of a successful hunting attempt. Blecha et al [27] found that kill sites, compared to preceding locations, occurred more frequently in areas with higher housing densities and lower topographic positions, such as drainage areas, despite drainage areas having lower prey availability. We did not measure hunting success, but we did find that cougars moved toward lower elevations at dusk, when cougars are likely to hunt or return to a carcass. We found that individual response to nearest potential kill site was variable within the two-week period and across individuals; this is likely a function of timing of successful kills and the size of the prey item, with stronger positive responses being correlated with larger prey (as individuals return to the site over a longer period of time). In addition, the variation in motility across months may also be a function of the available prey �Buderman et al. Movement Ecology (2018) 6:22 (a) Page 10 of 16 (b) (c) (d) Fig. 4 Posterior means and 95% credible intervals for the individual- and population-level static effects of elevation (Fig. 4a), heat loading (Fig. 4b), distance to nearest roofed structure (Fig. 4c), and distance to nearest road (Fig. 4d) on cougar motility in the Colorado Front Range for two-week periods in April, June, and October 2011 April 0.50 June October 0.25 0.00 −0.25 0 2 4 6 8 10 12 14 16 18 20 22 24 Hour Fig. 5 Posterior means and 95% credible intervals for the population-level diel time-varying effect of distance to nearest potential kill site on cougar motility in the Colorado Front Range for two-week periods in April, June, and October 2011. The gray box represents 0630 h to 1930 h. �Buderman et al. Movement Ecology (2018) 6:22 Page 11 of 16 April 0.50 June October 0.25 0.00 −0.25 0 2 4 6 8 10 12 14 16 18 20 22 24 Hour Fig. 6 Posterior means and 95% credible intervals for the population-level diel time-varying effect of distance to nearest potential kill site on directionality of cougar movement in the Colorado Front Range for two-week periods in April, June, and October 2011. The gray box represents 0630 h to 1930 h April 0.25 June October Elevation: Motility Elevation: Direction 0.00 (a) −0.25 0.25 (b) Distance to Structure: Direction Distance to Structure: Motility 0.00 (c) (d) −0.25 0 2 4 6 8 10 12 14 16 18 20 22 24 0 2 4 6 8 10 12 14 16 18 20 22 24 Hour Fig. 7 Posterior means and 95% credible intervals for the population-level diel time-varying effect of elevation (Fig. 7a, b) and distance to nearest roofed structure (Fig. 7c, d) on motility and directionality of cougar movement in the Colorado Front Range for two-week periods in April, June, and October 2011. The gray box represents 0630 h to 1930 h �Buderman et al. Movement Ecology (2018) 6:22 Page 12 of 16 (a) (b) Fig. 8 Posterior means and 95% credible intervals for the population-level effect of distance to nearest potential kill site on motility and directionality of cougar movement (Fig. 8a) and average movement rate (Fig. 8b) as a function of development (developed being < 10 acres/unit) in the Colorado Front Range for two-week periods in April, June, and October 2011 Our results indicated that cougars moved slightly faster in areas with a higher heat insolation load index in April compared to June and October. These areas correspond to xeric, south-facing slopes, which, in the montane zone of the Front Range, mostly consist of open stands of ponderosa pine, compared to the more dense items during a given season. We observed the weakest response to kill sites during April, before mule deer fawns are born [38]. We also observed an increasingly strong positive response to distance to kill site from dusk to dawn, implying that, from dusk to dawn, individuals moved increasingly faster the farther away they were from a potential kill site. Model 0.00100 GLM GAM FGAM+MGLM MGAM+FGLM Density 0.00075 0.00050 0.00025 0.00000 1202000 1206000 1210000 April 1214000 1055000 1057500 1060000 1062500 June 1065000 1067500 1195000 1200000 October Fig. 9 Posterior predictive score distributions for the model set, where smaller values indicate better predicting models. The GAM fit to all individuals performed the best in terms of predictive ability for all months 1205000 �Buderman et al. Movement Ecology (2018) 6:22 north-facing slopes [64]. The more open forest floor may facilitate cougars using south-facing slopes as travel corridors, leading to greater transition rates. Similarly, [17] found that cougars used less rugged terrain than the surrounding area while traveling, while [22] found that cougars selected for south-facing slopes and areas with shrub habitat. The monthly difference in effect size for the response to heat loading may be related to seasonal changes in vegetation (shrub cover may be denser in June, reducing speed) or a product of unobserved weather patterns (e.g., more snow on north-facing slopes could lead to a greater tendency to use south-facing slopes as corridors). Similarly, we also observed slower movements at high elevations in April compared to June and October; high elevation areas of the Front Range may still contain snow in April, which could result in slower movement rates at high elevations. The time-varying directional response to elevation indicates that individuals are moving to higher elevations during the day, and then toward lower elevations at night. Blecha et al [27] found that predation events, which typically occur at night, occurred at lower elevations, which may explain the temporal pattern we observed. The difference among seasons for the dawn motility response (slower [faster] at high [low] elevations in April, faster [slower] at high [low] elevations in October) might also be a response to unobserved fine-scale vegetation changes or dietary shifts. The documented response of cougars to disturbed and developed landscapes varies in the literature, and is likely a function of the level of disturbance encountered and how disturbance was quantified. For example, [21] found no difference in cougar movement rates in wildland and residential areas throughout the day. However, [22] observed that cougars avoided developed landscapes, while also noting a temporal shift in usage of those areas, with cougars avoiding areas near development more during the day. On average, we observed a more negative relationship to distance from structures at dawn and dusk compared to mid-day and evening (i.e., individuals moved faster when closer to structures during dawn and dusk than mid-day and evening, when there was a slight positive relationship between speed and distance to structure), which could be explained by increased human activity caused by the start and end of the workday. However, the uncertainty was fairly large for the diel effect of distance to nearest roofed structure, despite subtle positive and negative shifts. These minor differences, and an overall lack of consistent response, could be explained by unmeasured spatial and temporal relationships, such as individual interactions, fine-scale human disturbance (e.g., recreational activities, noise, and construction), and individual risk-avoidance strategies [67]. For example, [37] found that cougars showed stronger Page 13 of 16 avoidance of more consistent sources of anthropogenic disruption, such as neighborhoods, than intermittent sources, such as low-traffic roads. We detected a faster average movement rate in areas with more development, which is consistent with work by [17] and [68], who found that individuals expended more calories and moved further in developed areas, and [69] also found that large mammal movement rates vary as a function of local conditions. Distance to roads and structures may not be adequate proxies for how cougars perceive anthropogenic disturbance, and alternative measures (e.g. density or categorical variables), or responses (actual road crossings, e.g., [70]), may be more relevant for cougar movement behavior. Other studies have detected significant individual variation [21, 37], and [61] and [37] found that selection differed between males and females. We did not see consistent sex-specific responses to covariates, which could be due to the timing of the observations; for example, females may respond differently to males when breeding, but similarly at other times. In addition, males were underrepresented in our sample. Some of the unexplained individual variation could be due to the amount of anthropogenic landscape features each individual was likely to encounter in their movements [22, 61], as opposed to the amount of development in the immediate vicinity during a given movement. Benson et al [61] also hypothesized that the amount of development in many studies of cougar habitat selection has been too low to cause cougar behavioral changes. We propose that our findings regarding a lack of evidence for significant landscape drivers of movement may have three potential biological causes. First, cougars are generalists, therefore, they are expected to demonstrate less habitat selection at the landscape scale than a habitat specialist would [71]. Though we were assessing movement, generalists may likewise demonstrate less variation in speed and directionality as a function of the landscape than would a specialist, or the risks and rewards present in the Front Range are not significant enough to cause a detectable response in behavior. In addition, individuals moving within an established home range, such as in this study, may be acclimated (demonstrating minimal change in behavior) to the disturbances that they encounter during daily movements. Second, significant individual-variation within a 24-h period would make determining a consistent population-level response difficult. Individual variation can occur across and within individuals, and may be a function of the unmeasured internal state of the animal (breeding status, body condition, and energetics), or external factors (interactions with other individuals, fine-scale landscape features, and unmeasured brief disturbances). Finally, cougar movement may correspond to, or interact with, �Buderman et al. Movement Ecology (2018) 6:22 lower frequency environmental variation (e.g., weather patterns and food availability). Comparing behavior among years could be used to assess seasonal consistency in observed patterns; however, due to the difficulty of performing multiple comparisons of time-varying effects, a study of among year differences would likely need to focus on a particular season of interest. Modelling framework Many of the hypothesized movement drivers did not have a consistent statistically significant relationship with movement, while other studies on cougars in the same geographic area have found strong effects for landscape variables on cougar resource use while hunting [27]. However, unlike RSFs and traditional SSFs (not integrated, e.g. [10]), the CTMC framework is measuring the effect of landscape variables on speed and directionality, not habitat selection. For example, cougars may select for areas with high mule deer use [27], but cougars may not alter their speed based on the amount of mule deer usage. Although there has been a recent development that uses the limiting distribution of a CTMC to estimate utilization density (home-range; [20]), the relationship between changes in speed and directionality and habitat selection may vary by species. In addition, although multiple imputation appropriately accounts for the uncertainty in the unobserved true path, it does introduce an additional source of variation that is not accounted for when using only the observed locations. In addition, the variation among imputed paths for each individual will increase with measurement error, because there is less information to constrain the realizations of the true, unobserved path. Although measurement error was minimal in our study, the added uncertainty introduced by the multiple imputation framework may make it difficult to detect statistically significant effects. The varying coefficient modeling framework, implemented in this study as a GAM, can reveal hidden process dynamics [72] and allows for complex nonlinear patterns that would be difficult to model in a traditional framework (e.g., [73]). While we expanded each parameter in to the temporal space, one could make each covariate a function of another parameter, such as a different temporal predictor (e.g., time since kill) or another parameter in the model (e.g., distance to structure). Allowing the coefficients to vary in time (or another covariate space) can also improve the predictive ability of the model, as it did in our study. GLMs can mask time-varying responses to covariates (e.g., [74]), because the response variable is aggregated over the time period of interest. Therefore, if the response of an individual switches between positive and negative (faster or slower movement rates), the estimated response will Page 14 of 16 be approximately zero. Studies have found that cougars use a broader range of habitats for nocturnal movements than for daybed locations [17] and demonstrate temporal variability in their response to anthropogenic landscape features [22]. Therefore, restricting analysis of locations to a particular temporal subset (e.g. day vs night) may not be indicative of all behavior [75]. The time-varying CTMC framework represents an important step forward in detecting latent temporal patterns in animal movement and is especially useful when behavior is known to vary in time. Conclusions Recursive events, as measured by potential kill site locations, were identified as the primary driver of motility and directionality for cougars in the Front Range of Colorado. Observed cougars also moved faster, on average, in developed areas compared to undeveloped areas. Many other landscape features, including proxies for anthropogenic development, did not have a strong population-level effect on cougar movement, potentially due to unexplained individual-level variation. We did not detect a link between cougars that had reports of human conflict and response to development. Cougars have demonstrated different second- and third-order selection to roads in previous studies [33], therefore, nuisance individuals may select for, or end up in, home ranges near human development, but do not respond differentially to areas closer to development within their home range [76]. The high degree of individual variation suggests that, if agencies want to minimize human-wildlife conflict, a “one size fits all” approach to cougar management and conflict abatement will likely be unsuccessful and management options should be varied and flexible. A potential proactive mitigation of cougar conflict is to identify potential kill sites and orient recreation (e.g., trails, camp sites) away from those areas. Kill site identification tools, such as those developed by [27], could prove useful in this regard. We observed temporal variation in the population-level response to some landscape features (potential kill site, elevation, distance to structure), which highlights the importance of considering time-varying effects of covariates on movement behavior. Time-varying effects may be particularly important to consider when animal behavior is known to vary in time and when temporally static covariates may contain uncharacterized temporal variation (e.g., roads that vary in traffic load according to time of day and season). This study is also the first hierarchical application of the CTMC modeling framework, and demonstrates its ability to provide computationally efficient inference on individual- and population-level drivers of animal movement behavior. �Buderman et al. Movement Ecology (2018) 6:22 Additional file Additional file 1: Continuous-Time Markov Chain Model Details. Provides additional details on the CTMC model specification. (PDF 133 kb) Abbreviations CPW: Colorado Parks and Wildlife; CTMC: Continuous-time Markov chain; H-GAM: Hierarchical generalized additive modelH-GLMHierarchical generalized linear modelGPSGlobal positioning system; RSF: Resource selection function; SSF: Step selection function Acknowledgements Any use of trade, firm, or product names is for descriptive purposes only and does not imply endorsement by the U.S. Government. We would like to thank the reviewers, whose comments helped to improve the manuscript. Funding Funding was provided by Colorado Parks and Wildlife (1304), the National Park Service (P12AC11099), Colorado Department of Transportation, NSF DMS 1614392, and NSF EEID 1414296. Availability of data and materials The datasets analysed during the current study are available from Mat Alldredge (mat.alldredge@state.co.us) on reasonable request. Authors’ contributions MA oversaw initial data collection. FEB analyzed the data, MBH and EMH provided statistical support, and JSI and MA provided ecological knowledge and assisted FEB with interpretation of results. FEB took the lead in writing the manuscript, and all authors contributed to, and approved, the final manuscript. Ethics approval Data were collected according to protocols approved by the Colorado Parks and Wildlife Animal Care and Use Committee (CPW ACUC #01–2008). Page 15 of 16 4. 5. 6. 7. 8. 9. 10. 11. 12. 13. 14. 15. 16. Consent for publication Not applicable. 17. Competing interests The authors declare that they have no competing interests. 18. 19. Publisher’s Note 20. Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations. 21. Author details 1 Colorado State University, Departments of Fish, Wildlife, and Conservation Biology, 1484 Campus Delivery, Fort Collins, CO 80523, USA. 2U.S. Geological Survey, Colorado Cooperative Fish and Wildlife Research Unit, Departments of Fish, Wildlife, and Conservation Biology and Statistics, Colorado State University, 1484 Campus Delivery, Fort Collins, CO 80523, USA. 3Colorado Parks and Wildlife, 317 W Prospect Road, Fort Collins, CO 80526, USA. 4 Pennsylvania State University, W-250 Millennium Science Complex, University Park, State College, PA 16802, USA. 22. 23. 24. Received: 19 June 2018 Accepted: 26 September 2018 25. References 1. Wiens JA, Stenseth NC, Van Horne B, Ims RA. Ecological mechanisms and landscape ecology. Oikos 1993;66:369–380. 2. Johnson DH. The comparison of usage and availability measurements for evaluating resource preference. Ecology. 1980;61:65–71. 3. Manly B, McDonald L, Thomas D, McDonald T, Erickson W. Resource selection by animals: statistical analysis and design for field studies. New York (NY): Kluwer academic publishers; 2007. 26. 27. 28. Christ A, Ver Hoef J, Zimmerman DL. An animal movement model incorporating home range and habitat selection. Environ Ecol Stat. 2008;15:27–38. Johnson DS, Thomas DL, Ver Hoef JM, Christ A. A general framework for the analysis of animal resource selection from telemetry data. Biometrics. 2008; 64:968–76. Brost BM, Hooten MB, Hanks EM, Small RJ. Animal movement constraints improve resource selection inference in the presence of telemetry error. Ecology. 2015;96:2590–7. Boyce MS, Mao JS, Merrill EH, Fortin D, Turner MG, Fryxell J, Turchin P. Scale and heterogeneity in habitat selection by elk in Yellowstone National Park. Ecoscience. 2003;10:421–31. Fortin D, Beyer HL, Boyce MS, Smith DW, Duchesne T, Mao JS. Wolves influence elk movements: behavior shapes a trophic cascade in Yellowstone National Park. Ecology. 2005;86:1320–30. Forester JD, Im HK, Rathouz PJ. Accounting for animal movement in estimation of resource selection functions: sampling and data analysis. Ecology. 2009;90:3554–65. Avgar T, Potts JR, Lewis MA, Boyce MS. Integrated step selection analysis: bridging the gap between resource selection and animal movement. Methods Ecol Evol. 2016;7:619–30. Johnson DS, London JM, Lea MA, Durban JW. Continuous-time correlated random walk model for animal telemetry data. Ecology. 2008; 89:1208–15. Patterson TA, Thomas L, Wilcox C, Ovaskainen O, Matthiopoulos J. State-space models of individual animal movement. Trends Ecol Evol. 2008;23:87–94. Buderman FE, Hooten MB, Ivan JS, Shenk TM. A functional model for characterizing long-distance movement behaviour. Methods Ecol Evol. 2016; 7:264–73. Hooten MB, Johnson DS. Basis function models for animal movement. J Am Stat Assoc. 2017;112:578–89. Hooten MB, Buderman FE, Brost BM, Hanks EM, Ivan JS. Hierarchical animal movement models for population-level inference. Environmetrics. 2016;27: 322–33. Aune KE. Increasing mountain lion populations and human-mountain lion interactions in Montana. In: Braun CS, editor. Mountain lion-human interaction symposium and workshop. Fort Collins (CO): Colorado Division of Wildlife; 1991. p. 86–94. Dickson BJ, Jenness JS, Beier P. Influence of vegetation, topography, and roads on cougar movement in southern California. J Wildl Manag. 2005;69: 264–76. Hooten MB, Johnson DS, Hanks EM, Lowry JH. Agent-based inference for animal movement and selection. J Agric Biol Environ Stat. 2010;15:523–38. Hanks EM, Hooten MB, Alldredge MW. Continuous-time discrete-space models for animal movement. Ann Appl Stat. 2015;9:145–65. Wilson K, Hanks E, Johnson D. Estimating animal utilization densities using continuous-time Markov chain models. Methods Ecol Evol. 2018;9:1232–40. Kertson BN, Spencer RD, Marzluff JM, Hepinstall-Cymerman J, Grue CE. Cougar space use and movements in the wildland-urban landscape of western Washington. Ecol Appl. 2001;21:2866–81. Knopff AA, Knopff KH, Boyce MS, Clair CCS. Flexible habitat selection by cougars in response to anthropogenic development. Biol Conserv. 2014;178: 136–45. Burdett CL, Crooks KR, Theobald DM, Wilson KR, Boydston EE, Lyren LM, Fisher RN, Vickers TW, Morrison SA, Boyce WM. Interfacing models of wildlife habitat and human development to predict the future distribution of puma habitat. Ecosphere. 2010;1:1–21. Torres RT, Carvalho J, Panzacchi M, Linnell J, Fonseca C. Comparative use of forest habitats by roe deer and moose in a human-modified landscape in southeastern Norway during winter. Ecol Res. 2011;6:781–9. Torres SG, Mansfield TM, Foley JE, Lupo T, Brinkhaus A. Mountain lion and human activity in California. Wildl Soc Bull. 1996;24:451–60. Manfredo MJ, Zinn HC, Sikorowski L, Jones J. Public acceptance of mountain lion management: a case study of Denver, Colorado, and nearby foothills areas. Wildl Soc Bull. 1998;26:964–70. Blecha KA, Boone RB, Alldredge MW. Hunger mediates apex predator’s risk avoidance response in wildland–urban interface. J Anim Ecol. 2018;87:609–22. Moss WE, Alldredge MA, Logan KA, Pauli JN. Human expansion precipitates niche expansion for an opportunistic apex predator (Puma concolor). Sci Rep. 2016;6:39639. �Buderman et al. Movement Ecology (2018) 6:22 29. Wittemyer G, Elsen P, Bean WT, Burton ACO, Brashares JS. Accelerated human population growth at protected area edges. Science. 2008;321:123–6. 30. White PC, Ward AI. Interdisciplinary approaches for the management of existing and emerging human-wildlife conflicts. Wildl Res. 2011;37:623–9. 31. Anderson AE, Bowden DC, Kattner DC. The puma on Uncompahgre plateau, Colorado: technical publication 40. Fort Collins (CO): Colorado Division of Wildlife; 1992. 32. Beier P. Dispersal of juvenile cougars in fragmented habitat. J Wildl Manag. 1995;59:228–37. 33. Dickson BG, Beier P. Home-range and habitat selection by adult cougars in southern California. J Wildl Manag. 2002;66:1235–45. 34. Hemker TP, Lindzey FG, Ackerman BB. Population characteristics and movement patterns of cougars in southern Utah. J Wildl Manag. 1984;48: 1275–84. 35. Ruth TK, Logan KA, Sweanor LL, Hornocker MG, Temple LJ. Evaluating cougar translocation in New Mexico. J Wildl Manag. 1998;62:1264–75. 36. Sweanor LL, Logan KA, Hornocker MG. Cougar dispersal patterns, metapopulation dynamics, and conservation. Conserv Biol. 2000;14:798–808. 37. Wilmers CC, Wang Y, Nickel B, Houghtaling P, Shakeri Y, Allen ML, KermishWells J, Yovovich V, Williams T. Scale dependent behavioral responses to human development by a large predator, the puma. PLoS One. 2013;8: e60590. 38. Pojar TM, Bowden DC. Neonatal mule deer fawn survival in west-Central Colorado. J Wildl Manag. 2004;68:550–60. 39. Hornocker MG. An analysis of mountain lion predation upon mule deer and elk in the Idaho primitive area. Wildl Monogr. 1970;21:1–39. 40. Moss WE, Alldredge MW, Pauli JN. Quantifying risk and resource use for a large carnivore in an expanding urban–wildland interface. J Appl Ecol. 2016; 53:371–8. 41. Hanks E. ctmcmove: Modeling Animal Movement with Continuous-Time Discrete-Space Markov Chains. 2017. 42. Buderman FE, Hooten MB, Ivan JS, Shenk TM. Large-scale movement behavior in a reintroduced predator population. Ecography. 2018;41:126–39. 43. Scharf HR, Hooten MB, Johnson DS. Imputation approaches for animal movement modeling. J Agric Biol Environ Stat. 2017;22:335–52. 44. Rubin DB. Multiple imputation for nonresponse in surveys. New York (NY): John Wiley & Sons; 1987. 45. Hooten MB, Hobbs N. A guide to Bayesian model selection for ecologists. Ecol Monogr. 2015;85:3–28. 46. Hastie T, Tibshirani R. Varying-coefficient models. J R Stat Soc Ser B Methodol. 1993;55:757–96. 47. Hefley TJ, Broms KM, Brost BM, Buderman FE, Kay SL, Scharf HR, Tipton JR, Williams PJ, Hooten MB. The basis function approach for modeling autocorrelation in ecological data. Ecology. 2017;98:632–46. 48. R Core Team. R: a language and environment for statistical computing. R Foundation for statistical computing. Austria: Vienna; 2016. 49. Blecha KA. Risk-reward tradeoffs in the foraging strategy of cougar (Puma concolor): prey distribution, anthropogenic development, and patch selection. MSc [thesis]. Fort Collins (CO): Colorado State University; 2015. 50. Anderson CR Jr, Lindzey FG. Estimating cougar predation rates from GPS location clusters. J Wildl Manag. 2003;67:307–16. 51. Blake LW, Gese EM. Resource selection by cougars: influence of behavioral state and season. J Wildl Manag. 2016;80:1205–17. 52. Fagan WF, Lewis MA, Auger-Méthé M, Avgar T, Benhamou S, Breed G, LaDage L, Schlägel UE, Tang WW, Papastamatiou YP, Forester J. Spatial memory and animal movement. Ecol Lett. 2013;16:1316–29. 53. Berger-Tal O, Bar-David S. Recursive movement patterns: review and synthesis across species. Ecosphere. 2015;6:1–12. 54. Piper WH. Making habitat selection more “familiar”: a review. Behav Ecol Sociobiol. 2011;65:1329–51. 55. Oliveira-Santos LGR, Forester JD, Piovezan U, Tomas WM, Fernandez FA. Incorporating animal spatial memory in step selection functions. J Anim Ecol. 2016;85:516–24. 56. Avgar T, Baker JA, Brown GS, Hagens JS, Kittle AM, Mallon EE, McGreer MT, Mosser A, Newmaster SG, Patterson BR, Reid DE. Space-use behaviour of woodland caribou based on a cognitive movement model. J Anim Ecol. 2015;84:1059–70. 57. Brost BM, Hooten MB, Small RJ. Leveraging constraints and biotelemetry data to pinpoint repetitively used spatial features. Ecology. 2017;98:12–20. Page 16 of 16 58. Knopff KH, Knopff AA, Warren MB, Boyce MS. Evaluating global positioning system telemetry techniques for estimating cougar predation parameters. J Wildl Manag. 2009;73:586–97. 59. Husseman JS, Murray DL, Power G, Mack C, Wenger C, Quigley H. Assessing differential prey selection patterns between two sympatric large carnivores. Oikos. 2003;101:591–601. 60. Nicholson KL, Krausman PR, Smith T, Ballard WB, McKinney T. Mountain lion habitat selection in Arizona. Southwest Nat. 2014;59:372–80. 61. Benson JF, Sikich JA, Riley SP. Individual and population level resource selection patterns of mountain lions preying on mule deer along an urbanwildland gradient. PLoS One. 2016;11:e0158006. 62. Theobald DM, Harrison-Atlas D, Monahan WB, Albano CM. Ecologicallyrelevant maps of landforms and physiographic diversity for climate adaptation planning. PLoS One. 2015;10:e0143619. 63. McCune B, Keon D. Equations for potential annual direct incident radiation and heat load. J Veg Sci. 2002;13:603–6. 64. Veblen TT, Donnegan JA. Historical range of variability for forest vegetation of the national forests of the Colorado front range. Golden (CO): USDA Forest Service, Rocky Mountain Region; 2005. 65. Beven K, Kirkby MJ. A physically based, variable contributing area model of basin hydrology. Hydrol Sci J. 1979;24:43–69. 66. Theobald D. LCaP v1. 0: landscape connectivity and pattern tools for ArcGIS. Fort Collins (CO): Colorado State University; 2007. 67. Scrafford MA, Avgar T, Heeres R, Boyce MS. Roads elicit negative movement and habitat-selection responses by wolverines (Gulo gulo luscus). Behav Ecol. 2018;29:534–42. 68. Wang Y, Smith JA, Wilmers CC. Residential development alters behavior, movement, and energetics in an apex predator, the puma. PLoS One. 2017; 12:e0184687. 69. Avgar T, Mosser A, Brown GS, Fryxell JM. Environmental and individual drivers of animal movement patterns across a wide geographical gradient. J Anim Ecol. 2013;82:96–106. 70. Prokopenko CM, Boyce MS, Avgar T. Characterizing wildlife behavioural responses to roads using integrated step selection analysis. J Appl Ecol. 2017;54:470–9. 71. Katnik DD, Wielgus RB. Landscape proportions versus Monte Carlo simulated home ranges for estimating habitat availability. J Wildl Manag. 2005;69:20–32. 72. Fan J, Zhang W. Statistical methods with varying coefficient models. Stat Interf. 2008;1:179–95. 73. Polansky L, Robbins MM. Generalized additive mixed models for disentangling long-term trends, local anomalies, and seasonality in fruit tree phenology. Ecol Evol. 2013;3:3141–51. 74. Cheng MY, Zhang W, Chen LH. Statistical estimation in generalized multiparameter likelihood models. J Am Stat Assoc. 2009;104:1179–91. 75. Comiskey EJ, Bass OL Jr, Gross LJ, McBride RT, Salinas R. Panthers and forests in South Florida: an ecological perspective. Conserv Ecol. 2002;6:18. 76. Linnell JDC, Odden J, Smith ME, Aanes R, Swenson JE. Large carnivores that kill livestock: do “problem individuals” really exist? Wildl Soc Bull. 1999;27: 698–705. � https://cpw.cvlcollections.org/files/original/9d61ca9d47a5c4d37ccc868196ef3beb.pdf 9efd70b7e6948e78af400ac4cdf053c9 PDF Text Text Appendix A Continuous-Time Markov Chain Model Details We spatially discretize a posterior predictive continuous path from the movement model to the resolution of the rasters of interest and decompose it into two elements: c, a state sequence consisting of the sequential grid cells (of N possible grid-cells) visited by the individual, and τ , a vector of residence times that describe how long the individual spent in each grid cell. It is important to note that individuals cannot skip grid cells (enter a non-neighboring cell), therefore the spatial resolution of the rasters should mirror our prior knowledge about movement speeds. We describe the cell sequence in terms of the transition rates α where αij is a parameter controlling movement from cell i to cell j that can be a function of spatial covariates: 0 αij = exij β (1) If we designate t as the tth observation in the state-sequence (t ∈ T ), then the residence time τt is exponentially distributed with a rate equal to the sum of all αij (the total transition rate): [τt |β] = N X ! αij e−τk PN j=1 αij . (2) j=1 In the above notation, [τt |β] represents the probability distribution of the random variable τt given the parameters β; this notation will appear again. We assume that it is impossible to move directly to non-neighboring cells, and therefore αij = 0 for all j except for the cells adjacent to cell i. When an individual transitions to a neighboring cell, the probability of transitioning to cell ct+1 = l is αil [ct+1 = l|ct = i] = PN j=1 αij . (3) Assuming independence, the joint likelihood is the product of the transition probabilities 1 �and the residence times in the state sequence c is: N X αil [τt , ct+1 = l|ct = i, β] = PN j=1 αij = αil e−τt PN ! e−τt αij PN j=1 αij (4) j=1 j=1 αij (5) Using a latent variable representation, where zij =   1, if j = ct+1 (6)  0, if j 6= ct+1 and z [zij , τt |β] ∝ αijij e−τt αij , (7) then the product of [zct k , τt |β] over all N is proportional to the likelihood of the observed transition: [zij , τt |β] ∝ T X N X z e−τt αij αijij . (8) t=1 j=1 Additional details can be found in Hanks et al. 2015. The above process is parameterized with a single realization from the movement model. To avoid computational storage limitations, we use multiple imputation to account for the uncertainty in the path and make approximate posterior predictive inference on transition rates. Literature Cited Hanks, E.M., M.B. Hooten, and M.W. Alldredge. 2015. Continuous-time discrete-space models for animal movement. Annals of Applied Statistics 9:145165. 2 � Dublin Core The Dublin Core metadata element set is common to all Omeka records, including items, files, and collections. For more information see, http://dublincore.org/documents/dces/. Title A name given to the resource Journal Articles Description An account of the resource CPW peer-reviewed journal publications Text A resource consisting primarily of words for reading. Examples include books, letters, dissertations, poems, newspapers, articles, archives of mailing lists. Note that facsimiles or images of texts are still of the genre Text. Dublin Core The Dublin Core metadata element set is common to all Omeka records, including items, files, and collections. For more information see, http://dublincore.org/documents/dces/. Title A name given to the resource Time-varying predatory behavior is primary predictor of fine-scale movement of wildland-urban cougars Rights Information about rights held in and over the resource <a href="http://rightsstatements.org/vocab/InC-NC/1.0/" target="_blank" rel="noreferrer noopener">In Copyright - Non-Commercial Use Permitted</a> Date Created Date of creation of the resource. 2018-11-02 Subject The topic of the resource Animal movement Hierarchical model Individual variation Population-level Predation Telemetry Wildland-urban interface Description An account of the resource <h4 class="c-article__sub-heading">Background</h4> <p>While many species have suffered from the detrimental impacts of increasing human population growth, some species, such as cougars (<i>Puma concolor</i>), have been observed using human-modified landscapes. However, human-modified habitat can be a source of both increased risk and increased food availability, particularly for large carnivores. Assessing preferential use of the landscape is important for managing wildlife and can be particularly useful in transitional habitats, such as at the wildland-urban interface. Preferential use is often evaluated using resource selection functions (RSFs), which are focused on quantifying habitat preference using either a temporally static framework or researcher-defined temporal delineations. Many applications of RSFs do not incorporate time-varying landscape availability or temporally-varying behavior, which may mask conflict and avoidance behavior.</p> <h4 class="c-article__sub-heading">Methods</h4> <p>Contemporary approaches to incorporate landscape availability into the assessment of habitat selection include spatio-temporal point process models, step selection functions, and continuous-time Markov chain (CTMC) models; in contrast with the other methods, the CTMC model allows for explicit inference on animal movement in continuous-time. We used a hierarchical version of the CTMC framework to model speed and directionality of fine-scale movement by a population of cougars inhabiting the Front Range of Colorado, U.S.A., an area exhibiting rapid population growth and increased recreational use, as a function of individual variation and time-varying responses to landscape covariates.</p> <h4 class="c-article__sub-heading">Results</h4> <p>We found evidence for individual- and daily temporal-variability in cougar response to landscape characteristics. Distance to nearest kill site emerged as the most important driver of movement at a population-level. We also detected seasonal differences in average response to elevation, heat loading, and distance to roads. Motility was also a function of amount of development, with cougars moving faster in developed areas than in undeveloped areas.</p> <h4 class="c-article__sub-heading">Conclusions</h4> <p>The time-varying framework allowed us to detect temporal variability that would be masked in a generalized linear model, and improved the within-sample predictive ability of the model. The high degree of individual variation suggests that, if agencies want to minimize human-wildlife conflict management options should be varied and flexible. However, due to the effect of recursive behavior on cougar movement, likely related to the location and timing of potential kill-sites, kill-site identification tools may be useful for identifying areas of potential conflict.</p> Creator An entity primarily responsible for making the resource Buderman, Frances E. Hooten, Mevin B. Alldredge, Mathew W. Hanks, Ephraim M. Ivan, Jacob S. Language A language of the resource English Is Part Of A related resource in which the described resource is physically or logically included. Movement Ecology Format The file format, physical medium, or dimensions of the resource application/pdf Extent The size or duration of the resource. 16 pages Bibliographic Citation A bibliographic reference for the resource. Recommended practice is to include sufficient bibliographic detail to identify the resource as unambiguously as possible. <p>Buderman, F. E., M. B. Hooten, M. W. Alldredge, E. M. Hanks, and J. S. Ivan. 2018. Time-varying predatory behavior is primary predictor of fine-scale movement of wildland-urban cougars. Movement Ecology 6:22. <a href="https://doi.org/10.1186/s40462-018-0140-6" target="_blank" rel="noreferrer noopener">https://doi.org/10.1186/s40462-018-0140-6</a></p> <div></div> <p></p> Type The nature or genre of the resource Article