Describing spatiotemporal memory patterns using animal movement modelling
Peter R. Thompson, Andrew E. Derocher, Mark A. Edwards, Mark A. Lewis
DDescribing spatiotemporal memory patternsusing animal movement modelling
Peter R. Thompson , Andrew E. Derocher , Mark A. Edwards ,and Mark A. Lewis Department of Biological Sciences, University of Alberta,Edmonton, AB, Canada Mammalogy Department, Royal Alberta Museum, Edmonton,AB, Canada Department of Renewable Resources, University of Alberta,Edmonton, AB, Canada Department of Mathematical and Statistical Sciences, Universityof Alberta, Edmonton, AB, Canada * Correspondence: [email protected] a r X i v : . [ q - b i o . Q M ] J a n Abstract
1. Spatial memory plays a role in the way animals perceive their environ-ments, resulting in memory-informed movement patterns that are observ-able to ecologists. Developing mathematical techniques to understandhow animals use memory in their environments allows for an increasedunderstanding of animal cognition.2. Here we describe a model that accounts for the memory of seasonal orephemeral qualities of an animal’s environment. The model captures mul-tiple behaviors at once by allowing for resource selection in the presenttime as well as long-distance navigations to previously visited locationswithin an animal’s home range.3. We performed a set of analyses on simulated data to test our model,determining that it can provide informative results from as little as oneyear of discrete-time location data. We also show that the accuracy ofmodel selection and parameter estimation increases with more locationdata.4. This model has potential to identify cognitive mechanisms for memoryin a variety of ecological systems where periodic or seasonal revisitationpatterns within a home range may take place. step selection function, cognitive map, spatial memory, grizzly bear,
Ursus arc-tos , animal movement, hidden Markov model
Animal movement modelling has rapidly emerged as a subfield of ecology (Nathanet al., 2008) due to advances in animal tracking (Skaug and Fournier, 2006) andcomputational technology (Kristensen et al., 2016). The products of these ad-vances have been widely applied to conservation and management (Fortin et al.,2005; Graham et al., 2012; Gerber et al., 2019). These models allow ecologiststo understand the size and shape of an animal’s home range (Worton, 1989) aswell as what habitat attributes animals prefer on a finer scale (Gaillard et al.,2010). To address the latter, ecologists have developed tools such as resourceselection functions (RSFs; Boyce and McDonald, 1999) and step selection func-tions (SSFs; Fortin et al., 2005). These allow for inference on an individual’shabitat preference in what is known as third-order selection (Johnson, 1980;Thurfjell et al., 2014). The fine temporal and spatial resolution of these modelsallows ecologists to draw inference about a variety of behavioral processes, suchas how an animal’s movement rates are affected by its environment (Avgar et al.,2016; Prokopenko et al., 2017) and how movement patterns change at different2emporal scales (Oliveira-Santos et al., 2016; Richter et al., 2020). And yet,even with the advances that have been made in animal movement modelling,some notable behavioral mechanisms are often not considered.Spatial memory, defined by Fagan et al. (2013) as memory of the spatialconfiguration of one’s environment, is one of the most important influenceson animal movement patterns. Many well-known behavioral processes, suchas home range emergence (Van Moorter et al., 2009), food caching (Claytonand Dickinson, 1998), and even migration (Bracis and Mueller, 2017; Merkleet al., 2019), require the ability to remember the spatial location of landmarksor regions. Animal species use spatial memory in different ways (Fagan et al.,2013), and the benefits an animal may receive from memory often depend on itsenvironment (Mueller and Fagan, 2008; Mueller et al., 2011). Theory on animalcognition has proposed that animals encode this spatial information in theirbrain as a cognitive map (Tolman, 1948). Ecologists have proposed multipletheories for the structure of these maps, with debate arising over whether aspatially explicit Euclidean map or a network-based topological map is moreaccurate (Bennett, 1996; Sturz et al., 2006; Normand and Boesch, 2009; Asensioet al., 2011). The true structure of these cognitive maps in animals is still unclearand may vary in different animal species.Many animals experience seasonal variation within their home ranges (Moreyet al., 2007; Wiktander et al., 2001), suggesting that a memory of these timingswould be beneficial to optimize foraging. A key tenet of optimal foraging theoryis that animals move to maximize their metabolic intake, and when the animaldoes decide to move, the timing of departure and the animal’s subsequent des-tination are both important (Charnov, 1976). For example, primates time theirjourneys to previously visited resource patches with optimal feeding conditionsas a means to maximize energetic intake (Janmaat et al., 2006). Sharks displayintra-population variation and plasticity in their partially migratory movements,highlighting the viability of these long-distance navigations as an efficient forag-ing tactic (Papastamatiou et al., 2013). These recursive movement patterns arenearly impossible without spatial (or spatiotemporal) memory, where the animalmust recall its previous interactions with the environment and associate thoseevents with both a location and a time (Fagan et al., 2013). Movement modelsthat incorporate spatial memory can provide insight on how human-animal con-flict (Buderman et al., 2018), habitat fragmentation (Marchand et al., 2017),and global warming (Mauritzen et al., 2001) affect memory-informed animals.Attempts to model these revisitations have proposed cognitive maps withspatial and temporal components, but have neglected to make inference aboutthe specific nature of these influences. While many such approaches exist(Dalziel et al., 2008; Avgar et al., 2013; Vergara et al., 2016; Harel and Nathan,2018), a common and simple technique involves integrating cognitive maps intoSSFs (Oliveira-Santos et al., 2016; Marchand et al., 2017). A notable exampleis the model developed by Schl¨agel and Lewis (2014), where cognitive mapsare based on time since last visit for each point in space. It is proposed thatanimals will only be encouraged to revisit locations when they have not visitedthem recently, as seen in some ecological systems (Davies and Houston, 1981).3his model was used to draw inference from gray wolf (
Canis lupus ) movementpatterns (Schl¨agel et al., 2017), but it does not provide detail on when animalschoose to revisit portions of their home range. The model only considers thelast visit to any point in space, disregarding any previous visits to that point.Time since last visit alone is insufficient to model the complex time-dependentspatial memory that inspires movement patterns described above, because wait-ing longer to revisit such locations may not always be beneficial for the animals(e.g., trees that lose their ripe fruit after too long).Here we describe a model that mathematically estimates the timing andprecision of these seasonally recursive movements, and thus can be used tomodel periodic revisitation patterns within an animal’s home range (Fig. 1).We employ innovative model fitting techniques (Kristensen et al., 2016; Fis-cher and Lewis, 2020) brought about by advances in computational methods todetect patterns in animal location data. Our modelling framework character-izes the movement of simulated or real animals according to four hypotheses:(N) the null hypothesis, assuming random walk behavior; (R) the resource-onlyhypothesis, assuming animals move entirely according to local resource gra-dients without memory; (M) the memory-only hypothesis, assuming animalsexhibit seasonal revisitation patterns within their home range with a prescribedperiodicity; and (RM) the resource-memory hypothesis, assuming animals aresimultaneously influenced by local resources and spatial memory.To test our model, we first simulated movement tracks according to themodel’s prescribed rules on simulated environments, subsequently analyzinghow sample size affects both model selection and parameter estimation. Wefound that even with data sizes equivalent to roughly one year of animal track-ing data, the model accurately identified movement patterns consistent with thefour different hypotheses and produced accurate parameter estimates. These re-sults improved when tracks with more locations were simulated. We then fitthe model to telemetry data from a population of Arctic grizzly bears (
Ursusarctos ) and performed the same simulation analysis with real landscape dataand movement parameters estimated for the bears. These bears live in a harshenvironment where food resources are seasonal (Edwards and Derocher, 2015)and sparsely distributed (Edwards et al., 2009). We found a heavy influence ofspatiotemporal memory in the bears’ movement patterns, although we deter-mined that more data may be required to analyze these populations than forsimulated movements.
Here we introduce a new modelling framework based on step selection functionsthat accounts for periodic or seasonal revisitations by animals that forage onephemeral resources (Fig. 1). We developed a nested structure of four modelsin discrete time and continuous space (see Table 1 for a summary of the param-eters and models) to address our four alternative hypotheses (N, R, M, RM).Our model fitting process, made possible through advanced automatic differ-4ntiation techniques, allows for further inference about the specific nature ofthese cognitive mechanisms. The novelty and complexity of the computationalprocesses used to analyze animal location data with our model motivated mul-tiple simulation-based studies to identify the statistical power and parameterestimability of our models.
We fit a hidden Markov model (HMM) to animal movement data to incorporatefor switching between stationary (or quasi-stationary) and non-stationary states.HMMs are a first-order Markov process, implying that the animal’s currentstate is entirely dependent on its most recent state. This approach is commonin movement ecology due to the multitude of behavioral strategies observed inforaging animals (Morales et al., 2004; Jonsen et al., 2013).An HMM consists of a Markov matrix A of state-switching probabilities aswell as conditional probability distributions of the animal’s spatial location foreach state (Jonsen et al., 2013). For a model with n different movement states, A maps from R n → R n , with each column summing to 1. Our model has twostates, so we can infer the structure of A from its diagonal. We denote theseentries λ and γ , representing the probability that the animal will stay in thestationary or non-stationary state, respectively, given it was just there.While our model is meant to be applied to continuous-space animal data, wemake an approximation by discretizing our landscape over a two-dimensionalsquare grid. Empirical landscape data is rarely continuous in space, and theresolution of this data can suggest a clear choice for the resolution of the domaingrid. We assume that the boundaries of this grid are reflective, so if the animalreaches its edge, it will simply “bounce back” and remain on the grid. We definepoints in continuous space as x (or x t to represent the animal’s location at time t ) and their corresponding grid cells as z or z t . Thus, x ∈ z is the animal’sinitial location.We define our conditional probability density functions for the stationaryand non-stationary state f s (which remains the same in all four models) and f ns , respectively. Each conditional probability distribution represents a first-order Markov process modelling the animal’s location x t and its bearing φ t over time, which depend only on x t − and φ t − from the previous time step.Due to observation error in animal tracking data, we assumed that the animal’sobserved location may change slightly even if it is not moving (Jonsen et al.,2013), so we allowed for small “movements” in our stationary state. The proba-bility distribution for bearings in the stationary state, g s ( φ t | φ t − ), is a uniformdistribution since we assume no directional autocorrelation here, so g s ( φ t | φ t − ) = 12 π , (1) f s ( x t , φ t | x t − , φ t − , ρ s ) = 2 π √ ρ s g s ( φ t | φ t − ) exp − (cid:107) x t − x t − (cid:107) πρ s . (2)5igure 1: Schematic describing our modelling framework. Given an animal’smovement track, quantified as a set of spatial coordinates, as well as landscapedata describing an animal’s environment, we fit four nested, competing modelsusing maximum likelihood estimation. The insight we gain from this processallowed us to make conclusions about the mechanistic drivers of animal behavior.6hen the animal is in the stationary state, the distance between x t and x t − has a half-Gaussian distribution with mean ρ s . The half-Gaussian distributionhas thinner tails than the more traditionally used exponential distribution, de-creasing the probability of longer movements from this state. We fix ρ s to reducemodel complexity, noting that it is fairly straightforward to do so based on theknown degree of observation error or the resolution of environmental data.In the non-stationary state, the animal may incorporate memory in its move-ments using a cognitive map (Tolman, 1948; Fig. 2). Our implementation ofa cognitive map expands on the concept of time since last visit (Davies andHouston, 1981; Schl¨agel and Lewis, 2014; Schl¨agel et al., 2017) by allowing forthe memory of more than just the last location to any point in space. Instead,we formulate the animal’s cognitive map as the set of times since previous visits(TSPVs) for any area in space. This form of episodic-like memory, where ani-mals remember a series of events with associated spatial locations and times, hasbeen displayed in birds (Clayton and Dickinson, 1998) and great apes (Martin-Ordas et al., 2010). These animals can associate a point in space with an event(or set of events) and time (or set of times) when they visited there (Martin-Ordas et al., 2010). We define this map Z t at each time t as a function overthe domain grid. At each grid cell z , Z t ( z ) is a linked list of integers, with eachelement of the list representing an animal’s visit to a point inside that cell. Z isa grid full of empty linked lists, except for z ; Z ( z ) is a list with one element,0. We can obtain Z t if we know Z t − as well as the animal’s location at time t .When t is incremented by 1, so is every entry on every linked list across the grid,and a new entry (0) is added to the linked list corresponding to the animal’snew location: Z t ( z ) = (cid:40) Z t − ( z ) + 1 x t / ∈ z [ Z t − ( z ) + 1 , x t ∈ z . (3)where [ Z t − ( z ) + 1 ,
0] implies adding 1 to every entry of the linked list Z t − ( z )and appending it with a new value 0.The function f ns , which models the animal’s location and bearing in the non-stationary state, resembles a step selection function (Fortin et al., 2005), withtwo main components: k , the resource-independent movement kernel; and W ,the environmental (or cognitive) weighting function. The function k describesthe animal’s locomotive capability while W , which may depend on the animal’scognitive map Z t − , describes how attractive the point is to the animal. Thisyields the following expression for f ns : f ns ( x t , φ t | x t − , φ t − , Z t − , Θ , Θ )= k ( x t | x t − , φ t − , Θ ) W ( x t | Z t − , Θ ) (cid:82) Ω k ( x (cid:48) | x t − , φ t − , Θ ) W ( x (cid:48) | Z t − , Θ ) dx (cid:48) . (4)The integral in the denominator serves as a normalization constant to ensurethat f ns integrates to 1. The parameter vector Θ represents parameters re-lated to the W and Θ represents the locomotive parameters associated with7igure 2: Diagram describing how an animal’s cognitive map Z changes overfour discrete time steps, given an animal’s movement track, which is illustratedin the shaded panels. Each cell of Z contains a linked list that starts out emptybut is iteratively appended as the animal traverses its environment. k , namely ρ ns which describes the animal’s mean step length and κ which de-scribes the degree of directional autocorrelation in the animal’s movements. Foreach of our four models (null, resource-only, memory-only, resource-memory),the animal’s resource-independent movement kernel k has the same formulation.The distance between x t and x t − , known as a step length, is exponentially dis-tributed with mean parameter ρ ns , and the bearing φ t is von Mises distributedaround φ t − with concentration parameter κ ≥ κ indicate straighter movement. We assume here that the animal’s steplengths and turning angles are independent. This modelling structure, knownmore generally as a correlated random walk, has been applied to a variety ofecological systems (Fortin et al., 2005; Codling et al., 2008; Auger-M´eth´e et al.,2015; Duchesne et al., 2015). We formulate k such that g ns ( φ t | φ t − ) = exp ( κ cos( φ t − φ t − ))2 πI ( κ ) , (5) k ( x t | x t − , φ t − , Θ ) = exp (cid:16) − (cid:107) x t − x t − (cid:107) ρ ns (cid:17) ρ ns g ns ( φ t | φ t − ) . (6)Notice that φ t is not explicitly included in the left side of Equation 6; it caninstead be calculated easily if x t and x t − are known. If we denote x t =8 x t, , x t, ) T ∈ R , then the following equation for φ t holds: φ t = tan − (cid:18) x t, − x t − , x t, − x t − , (cid:19) + πI ( x t, < x t − , ) , (7)where I ( x t, < x t − , ) is an indicator function that is 1 if x t, < x t − , and 0otherwise.The only mathematical difference between the four models is the formulationof W . To differentiate between these different formulations, we refer to themas W N , W R , W M , and W RM for the null, resource-only, memory-only, andresource-memory models, respectively. The set of parameters we estimate ineach model also varies, so we define Θ ,N , Θ ,R , Θ ,M , and Θ ,RM in a similarrespect. The null model describes an animal’s locomotive capability and directional au-tocorrelation based on its observed movement track. As a result, there is noextra weighting, so W N ( x t | Z t − , Θ ,N ) = 1 for all x t in space, and Θ ,N is theempty set. As a result, when considering the null model, f ns is equal to k . The resource-only model has the following key component:(R1) the animal’s movement is driven by third-order habitat selection in thepresent time within its immediate vicinity.This modelling component also exists in RSFs (Boyce and McDonald, 1999) andSSFs (Fortin et al., 2005), which have received ample attention in movementecology.As a result, W R resembles the weighting function from an SSF. If we areinterested in P different resource covariates (expressed mathematically at eachspatial location x as r ( x ) , ..., r P ( x )), we must estimate selection parameters β , ..., β P for each covariate. These parameters make up Θ ,R . The expressionfor our weighting function in the resource-only model is a linear combination ofthe covariates: W R ( x t | Z t − , Θ ,R ) = exp (cid:34) P (cid:88) p =1 β p r p ( x ) (cid:35) . (8)Notice that Z t − , the cognitive map, is not incorporated in the null or resource-only models since they do not include memory. The memory-only model contains the following key components:9M1) the animal uses a cognitive map to remember the timing of previous visitsto regions of its environment, and(M2) it will return to locations it previously visited at a prescribed and sched-uled time.This type of cognitive map has been supported in the literature (Normand andBoesch, 2009; Martin-Ordas et al., 2010; Schl¨agel and Lewis, 2014) as has thevalidity of path recursions and revisitations as a foraging strategy for animals(Berger-Tal and Bar-David, 2015; Schl¨agel et al., 2017).We calculate W M based on distance to previously visited points on the ani-mal’s track. Given some time lag τ , we can use the cognitive map Z t to find thepoint in space (or at least, the grid cell) where the animal was τ time indicesago. Tere is always exactly one grid cell z t − τ where τ is an element of the linkedlist Z t ( z t − τ ).For each time lag τ , we compute the distance between the animal’s currentlocation and z t − τ , (cid:107) x t − z t − τ (cid:107) , and transform it using an exponential decayfunction with decay parameter 10 α . α measures the order of magnitude onwhich the animal perceives its landscape spatially. Larger values of α indicatethat the animal perceives its landscape to be more homogeneous. If α is larger,then exp ( − α (cid:107) x − y (cid:107) ) is more likely to be small no matter how far apart thetwo points x and y are, suggesting that the animal cannot identify a differencebetween those distances. As α decreases, the mathematical difference betweena step 1000 m away and a step 2000 m away is amplified, suggesting that theanimal understands these differences in space on a wider scale.The animal’s periodic revisitation schedule, which is mediated by two pa-rameters, dictates the weights for each of these exponentially transformed dis-tances. The timing with which an animal navigates back to an existing locationcan be thought of as a random process, following a Gaussian distribution withmean parameter µ and standard deviation parameter σ . We can imagine thatthis timing reflects the state of the environment, with µ indicating the timescale at which resources may come and go and σ indicating the variability ofthese cycles. For any given time lag τ , the exponentially transformed distancebetween x t and z t − τ is weighted by the Gaussian probability distribution func-tion ϕ ( τ | µ, σ ). This produces a weighted mean of exponentially transformeddistances, following the hypothesis that animals will navigate towards pointsthey visited roughly µ time increments ago; the most “attractive” points forthe animal are closest to z µ . We introduce one final parameter, β d , a “selectioncoefficient” for memorized locations. This parameter can be thought of as therelative probability of revisiting a memorized location instead of moving ran-domly or selecting for present-time resources. We restricted β d ≥ . β d − β d > W M is as follows:10 M ( x t | Z t − , Θ ,M ) = exp (cid:32) log (cid:18) β d − β d (cid:19) × (cid:34) (cid:80) tτ =1 ϕ ( τ | µ, σ ) exp ( − α (cid:107) x t − z t − τ (cid:107) ) (cid:80) tτ =1 ϕ ( τ | µ, σ ) (cid:35)(cid:33) . (9) The resource-memory model incorporates both resource selection and memoryinto the animal’s movements, so (R1) and (M1) still remain as components inthis model. However, there is one additional component that is not present inthe resource-only or memory-only models:(RM1) the animal will return to locations it previously visited at a prescribedand scheduled time if habitat conditions there were favorable; otherwiseit will avoid these areas.Models combining resources and memory in some way have proven to be effectivein explaining movement patterns for many different animals (Dalziel et al., 2008;Merkle et al., 2014; Schl¨agel et al., 2017).The linear combination of resource covariates (cid:80) Pp =1 β p r p ( x ) is relative, so weintroduced an additional parameter β representing the relative probability ofvisiting a faraway location depending on its resource quality. If β = 1 then theanimal perceives all previously visited locations as “attractive” for revisitation.We transform this parameter with an inverse logistic function so it representsa pseudo-intercept (recall that traditional SSFs and are conditional models anddo not require an intercept; Fortin et al., 2005).The weighting function now includes present-time resource selection in thefirst sum and memorized information in the second term: W RM ( x t | Z t − , Θ ,RM ) = exp (cid:32) p (cid:88) p =1 β p r p ( x t ) + log (cid:32) β d − β d (cid:33) × (10) (cid:34) (cid:80) tτ =1 ϕ ( τ | µ, σ ) exp ( − α (cid:107) x t − z t − τ (cid:107) )(log (cid:16) β − β (cid:17) + (cid:80) Pp =1 β p r p ( z t − τ )) (cid:80) tτ =1 ϕ ( τ | µ, σ ) (cid:35)(cid:33) . The null model is a special case of both the resource-only and memory-onlymodels, which are both a special case of the resource-memory model. Setting β i = 0 for i = 1 , , ..., P and log (cid:16) β − β (cid:17) = 1 in the resource-memory modelyields the memory-only model, while setting β d = 0 yields the resource-onlymodel. Nesting models is advantageous for many mathematical reasons, includ-ing the ability to conduct likelihood ratio tests between models (Burnham andAnderson, 2004). 11 nits Description N R M RM ρ ns distancetime Mean movement speed in non-stationary state X X X X κ N/A Degree of directional autocorrelation X X X X β N/A Probability of revisitation X β i r i units Resource selection coefficient(s) X X β d N/A Strength of selection for memorized areas X X µ time Mean time lag between revisitations X X σ time Standard deviation in time between revisitations X X α Degree of perceptual resolution X X λ N/A Probability of staying in stationary state X X X X γ N/A Probability of staying in non-stationary state X X X XTable 1: Description of model parameters, including units (N/A implies thatthe parameter is unitless) and models (N = null; R = resource-only; M =memory-only; RM = resource-memory) in which the parameters are estimated.
We fit the four models to discrete-time, continuous-space animal movementdata and used information theory to identify which corresponding hypothesiswas most likely to be true. We identified the optimal set of parameters for agiven track using maximum likelihood estimation, and used likelihood profilingto obtain accurate confidence intervals for our parameters.
The likelihood of a set of model parameters for one step is a weighted sum ofthe conditional likelihood functions ( f s and f ns ), weighted by the probabilityof being in each state. These state probabilities depend on probabilities for theprevious step, so for the first point we fit (there is no previous step), we fixed δ s , the probability of being in the stationary state right before the data begins,as the proportion of steps shorter than ρ s .The likelihood function for the entire track is a product of the likelihoodsfor each step included in model fitting. We omitted all animal locations beforesome time t ∗ , since our model (or at least, the memory-only and resource-memory models) relies on past information to explain where the animal maygo. We left the portion of the track that happened before t = t ∗ to “train”the model on what the animal remembers. Thus, our iterative formula for thelikelihood function begins at t = t ∗ . We define Φ t ∈ R as the vector of stateprobabilities for time t ≥ t ∗ , and we calculate our likelihood using the iterativeequations below: 12 t ∗ = ( δ s , − δ s ) T , (11) P t = (cid:18) f s ( x t | x t − , ρ s ) 00 f ns ( x t , φ t | x t − , φ t − , Z t − , Θ , Θ ) (cid:19) , (12)Φ t = Φ Tt − P t − (cid:107) P t − Φ t − (cid:107) A . (13)Then, following Whoriskey et al. (2017), the overall likelihood for the model is (cid:81) t m axt = t ∗ Φ Tt P t , where = (1 , T .We approximate the denominator of Equation 4 with a sum so we do not haveto integrate every time we evaluate the likelihood function. As is commonly donewith SSFs (Thurfjell et al., 2014), we calculated W at a set of “control points”for each observed point x t . If x t , the endpoint of a step from x t − , is a randomvariable conditional on Z t − , Θ , and Θ , the integral in the denominator ofEquation 4 is E ( W ( x t )). Thus, we can approximate it by estimating the meanvalue of W at a set of simulated draws from x t , which has probability densityfunction k . This gives us the following approximation for f ns :˜ f ns ( x t , φ t | x t − , φ t − , Z t − , Θ , Θ )= k ( x t | x t − , φ t − , Θ ) W ( x t | Z t − , Θ ) K (cid:80) Kj =1 W ( x t,j | Z t − , Θ ) , (14)where x t,j represents the j th control point (a simulated step starting at x t − )and K is the number of control points per observed step. We fit the model to data using maximum likelihood estimation, with the Tem-plate Model Builder (TMB) R package (Kristensen et al., 2016) improving nu-merical accuracy for this complex problem. TMB has been used to fit complexanimal movement models, including HMMs (Albertsen et al., 2015; Auger-M´eth´eet al., 2017; Whoriskey et al., 2017). TMB uses automatic differentiation tocalculate the gradient of a multidimensional likelihood function using pseudo-analytical methods, as opposed to traditional finite-difference methods that areslow and frequently result in numerical errors (Skaug and Fournier, 2006). Wewrote a likelihood function for each model in C++, which TMB compiles andreturns as a callable function in R (Kristensen et al., 2016). This allowed usto use an R optimizer of our choice while also benefiting from C++’s superiorprogramming speed.We used the R nlminb function to obtain maximum likelihood estimatesfor the negative log of our likelihood function. To prevent our model fromproducing errors or unrealistic results, we imposed various bounds on some ofthe parameters. We bounded the estimation for µ at t ∗ because if µ > t ∗ , wewould not be able to identify a signal due to a lack of training data. We also put13 lower bound on σ ; when this parameter was small, the partial derivative of ourlikelihood function with respect to µ became noisy, leading to computationalerrors in optimization. We found that a lower bound of approximately 20 timeindices eliminated this problem. We additionally required estimates for α < − log (¯ ρ ), where ¯ ρ is the animal’s empirical mean step length (for context, weexpect ¯ ρ to be close to but slightly smaller than ρ ns ). Values of α above thisbound imply that the animal cannot perceive a difference between a few steplengths, which is unreasonable biologically. For parameters with fairly restrictivebounds ( λ , γ , β d , and β , which are bounded between 0 and 1), we performedlogit transformations (˜ λ = log − λ , for example) so the optimizer would moreeffectively traverse the parameter space.We tested two “initial values” for µ for each dataset we fit the model to,picking the fit that gave us the best likelihood function value. When profilingthe likelihood surface with respect to this parameter, we often found many lo-cal optima, so we fit the model with initial values of t ∗ and t ∗ . This incursadditional computational time (we are effectively running the optimization al-gorithm twice) but is necessary due to the importance of picking a good initialvalue for each parameter (Pan and Wu, 1998). Using a different number ofinitial values for µ may be advantageous for some datasets.For a model as complicated as this one, obtaining confidence intervals (CIs)using traditional Wald-type methods does not always produce accurate results.We frequently found this to be true for our model in practice so we used thelikelihood profiling from Fischer and Lewis (2020). Given a multidimensionalobjective function with a known optimum, this algorithm finds confidence in-tervals for one parameter at a time by performing a binary search algorithmfor a target function value (typically, the optimum minus some small confidencethreshold). The algorithm starts searching at the optimal parameter value, andtries an initial step, fixing the parameter in question at this value and optimizingthe rest of the function parameters. This process is repeated subsequently untilthe lengths of each step in parameter space are small enough for the algorithmto converge (Fischer and Lewis, 2020).We used the Bayesian Information Criterion (BIC) to rank the four modelsby their likelihood and identify the hypothesis that was most likely to be true.BIC has a stronger penalization for model complexity than the more commonlyused Akaike Information Criterion (AIC), and is a more useful criterion formodel selection when one is interested in the truth of a hypothesis rather thanthe predictability of a model (Burnham and Anderson, 2004). Before applying our model to an ecological system, we simulated data and usedit to test the model. These simulations are individual-based representations ofour model that produce movement patterns associated with our four hypotheses.At each time index, we used our Markov matrix A to decide whether the animalwould change its behavioral state. If the animal was in the stationary state wesimulated a random step from f s (half-Gaussian step length, uniform turning14ngle). For the non-stationary state, we simulated from f ns using a MonteCarlo method (Parzen, 1961). We first calculated W for the entire grid, thenwe simulated a large number of random steps from k (Equations 5 and 6). Thissimulation process resembles the generation of control points in Equation 14,but we simulated N r = 10000 steps at each point in time. Making N r verylarge did not greatly affect computational time, so we did so in the interestof accurately approximating Equation 4. We then randomly choose one of thesteps based on the values of W at each step, with the probability of any step x t,i being chosen described below: W ( x t,i | Z t − , Θ ) (cid:80) N r j =1 W ( x t,j | Z t − , Θ ) . (15)For models that incorporate memory, we simulated memoryless training data( W M = W N for the memory-only model, and W RM = W R for the resource-memory model) for t < t ∗ . As expected, these initial points are omitted frommodel fitting. We simulated tracks on artificial landscapes with preset model parameters, thenfit the model to these tracks to explore parameter estimability and model se-lection accuracy. We varied the length of these tracks, T = t max − t ∗ , as wellas K , the number of controls points per step, to evaluate the amount of datarequired for accurate inference. Specifically, we tested four “treatment groups”: T = 600 , K = 10; T = 600 , K = 50; T = 1200 , K = 10; and T = 1200 , K = 50.We used the R NLMR package (Sciaini et al., 2018) to simulate spatially au-tocorrelated Gaussian random fields representing our resource covariates. Foreach treatment group, we simulated 50 random movement tracks for each hy-pothesis. Each group of 50 tracks had the same set of parameters. In oursimulations, we simulated environments for P = 3 resource covariates per trackusing the nlm gaussianfield function in R. We then fit all four models to eachtrack individually, then used BIC to identify how often the “correct” modelwas selected for each movement track. We compared these results with AIC toconfirm that BIC is the most suitabile information criterion for our modellingframework. We also estimated the bias and mean squared error (MSE; the meansquared difference between the parameter estimate and the true value) for eachparameter with each model. We applied the model to grizzly bears in the Canadian Arctic, and then re-peated the simulation study with data and model parameters from this system.Bears were captured from 2003 to 2006 and released with global positioning sys-tem (GPS) collars. Collars returned a location every four hours while the bearwas not hibernating, and remained on the bears for up to four years (Edwardset al., 2009). The University of Alberta Animal Care and Use Committee for15iosciences approved all animal capture and handling procedures, which werein accordance with the Canadian Council on Animal Care. The bears were col-lared in the Mackenzie River Delta region in the Northwest Territories (Edwardset al., 2009). Resources in the region are sparse and heterogeneous both in spaceand time (Shevtsova et al., 1995; Edwards and Derocher, 2015). To survive andforage optimally, bears take advantage of ephemeral, unpredictable, or season-ally available resources through a variety of foraging strategies (Edwards et al.,2009, 2011; Edwards and Derocher, 2015).We analyzed grizzly bear habitat selection using multiple sources of envi-ronmental data describing the Mackenzie Delta region. Vegetation class datafor the region assigned a one of 47 classes (indicating the dominant plant typeor terrain) to each 30x30 m cell. A digital elevation model for the region (with30x30 m cell resolution) provided information on elevation and slope. We alsoused an RSF layer estimating resource use for Arctic ground squirrels (
Urocitel-lus parryii ), a common grizzly bear prey item (Barker and Derocher, 2010;Edwards and Derocher, 2015). We considered P = 6 resource covariates: berrydensity, represented as a likelihood of having berries for each vegetation class;distance to turbid water, an indicator of broad whitefish ( Coregonus nasus ;a grizzly bear prey item; Barker and Derocher, 2009) density; Arctic groundsquirrel density, taken directly from the RSF; sweetvetch (
Hedysarum alpinum ;a key grizzly bear food item; Edwards and Derocher, 2015) density, estimatedby the vegetation class data; distance to the nearest of two towns in the region;and distance to six remote industrial camps (likely with little human activity).Of the 21 bears with enough data for model fitting (at least two years ofGPS collar data), we selected the eight with the most GPS fixes (these bearshad at least three years of collar data). We set ρ s = 30 meters, correspondingto the length of one grid cell for the environmental raster data, and we set t ∗ = 365 days. We used K = 50 control points when fitting the models. Foreach of these bears, we fit the models to the entire track as well as each yearindividually, comparing model selection between years. We then replicated thatanalysis using simulated bear tracks; for each bear, we simulated 100 movementtracks using the optimal parameters for each bear and the Mackenzie Deltaenvironmental data. We simulated tracks of length T = 600 (approximatelyone year of grizzly bear GPS data, accounting for missed fixes and hibernation)and T = 1200 to evaluate how model selection accuracy changed with samplesize. We used BIC to identify the hypothesis that most accurately explainedeach movement track, and also conducted likelihood ratio tests for each pair ofnested models to determine the significance of specific behavioral signals. Our modelling structure allows ecologists to explain movement patterns iden-tified from location data according to a set of four hypotheses, of which twoincorporate complex time-dependent spatial memory. For animals that appearto use memory, our parametric approach evaluates the timing and precision of16avigations to previously visited locations in an animal’s home range. By fit-ting the model to simulated data we showed that the accuracy of the model isimproved by sample size, and ecologists can also increase parameter estimabil-ity by simulating additional control points. Still, the amount of data requiredto draw accurate inference from the model is not large, as we show both withsimulated environments and real-life landscape data (where the model is slightlyless accurate).
The model’s ability to accurately characterize each type of movement behav-ior increased with the amount of location data ( T ) but not with the numberof control points ( K ; Table 2). The model identified null and resource-onlymovements accurately at all treatment levels, but the model’s ability to identifymemory-only and resource-memory movement increased for longer simulatedtracks. As a whole, increasing K does not improve model selection accuracyfor either choice of T . The most common misidentification at all sample sizeswas mistaking resource-memory movement for resource-only or memory-onlymovement. K = 10 K = 50N R M RM N R M RM T = N 48 0 0 2 47 0 0 3R 0 45 0 5 0 46 0 4M 7 0 40 3 4 0 44 2RM 0 5 7 38 0 8 7 35 T = N 49 0 0 1 45 0 0 5R 0 46 0 4 0 47 0 3M 2 0 47 1 0 0 50 0RM 0 3 2 45 0 7 2 41Table 2: Breakdown of model selection counts using BIC for the simulatedtracks. The row represents the “true” model that the tracks were simulatedfrom (N = null; R = resource-only; M = memory-only; RM = resource-memory),while the column represents the model that was identified as the most parsimo-nious explanation of the data using BIC. Treatment groups are identified by theouter left and upper portions of the table and are separated by shading.Using AIC instead of BIC resulted in a higher rate of “false positives” formemory (i.e., the resource-memory or the memory-only model was identified asthe most parsimonious explanation for memoryless simulated tracks), and mademodel selection less accurate overall (Appendix: Table S1). Likelihood ratiotests on the same dataset for each pair of nested models revealed a similar trend;the likelihood ratio test often identified memory when it was not incorporatedinto the simulated tracks (Appendix: Table S2).17igure 3: Violin plot of parameter estimates for β in the resource-memorymodel for our four treatment groups (listed on the x-axis), with 50 simulationsper plot. The true value of 7.5 is denoted by a horizontal red line.The model produced more accurate parameter estimates with larger valuesof T and K (Table 3). When focusing on β in the resource-memory model, wecan see that bias does not change as much with different treatment groups asMSE (Fig. 3). For the simpler movement parameters ( ρ ns , κ , λ , γ ), parameterestimates were consistent even with smaller values of T and K (Table 3). Memory played a significant part in the movements of five of the eight bearsin the population (Table 4). When the data were broken up into one-year in-crements, model selection results varied annually, and sometimes differed evenfrom the full dataset. For three of the bears (GF1008, GF1016, GM1046), themodel identified as most explanatory of the bears’ movement behaviors by BICwas different for the full dataset, the first subset, and the second subset. Theresource-memory model was the most parsimonious explanation of the move-ment patterns of four bears, while the resource-only (2), memory-only (1), andnull (1) models were also identified as most parsimonious in some cases. Four ofthe five memory-informed bears exhibited seasonal memory timescales close toone year ( µ >
320 days), while GF1016 had a µ value of 3 days. The six bearswith resource selection included in their “best model” displayed similar resourceselection patterns: significant selection for areas indicative of berries and Arctic18rue value T = 600 T = 1200 T = 600 T = 1200 K = 10 K = 10 K = 50 K = 50Bias MSE Bias MSE Bias MSE Bias MSE ρ ns κ β β β -7.5 -1.2 7.8 -1.1 5.5 0.2 3.8 0.5 1.5 β β d µ
500 -8 795 -13 543 -22 7819 -25 7277 σ
25 2 341 -3 32 1 197 -1 217 α -1.78 -0.45 2.64 -0.28 1.30 -0.14 1.81 -0.67 2.17 λ γ T = 600. With T = 1200, this improved to 89. Whenwe used likelihood ratio tests to compare the resource-memory model with theresource-only model (a special case of the resource-memory model) for GF1008,we found that at T = 600, 76 of our 100 simulated tracks registered a p-value19elow 0.05, indicating that the resource-memory model was significantly moreexplanatory than the resource-only model 76% of the time. With T = 1200,this increased to 95. We observed similar trends for the other three resource-memory bears (GF1004, GF1041, GM1046) but not as strongly. It should benoted that GF1008 had the smallest estimate for β d (2.3) of these bears. Whenwe performed BIC model selection on simulated tracks based on GF1086, a“resource-only” bear, a false memory signal was identified more frequently withlarger T (from 4 to 12 out of 100). This trend was not replicated for GF1107,the other “resource-only” bear (decrease from 10 to 7). Our model builds on existing literature to identify unique behavioral and cogni-tive mechanisms from animal movement data. Using advanced computationaltechniques, this novel and complex modelling framework can provide statisticalinference for a variety of ecological systems. Our simulation studies providedinsight on the viability of the model for different amounts of data.We formulated a model that expresses parameters with clear biological im-plications to aid in the interpretation of our results, but we had to do so carefullyto ensure that these parameters could be estimated accurately. Finding a set ofbiologically meaningful parameters with low mean squared error (Table 3) re-quired a degree of trial and error, especially for β and β d . We chose to expressthem in a way that makes sense both biologically (where they represent relativeprobabilities) and mathematically (where they can easily be estimated with lesserror). While we can redefine these parameters without actually changing ourlikelihood function, we made sure to define parameters that are easy to estimateand biologically meaningful.Our results provided support for a positive effect of the amount of locationdata and control points on parameter estimation, with the number of controlpoints having a negligible effect on model selection accuracy. However, at alltreatment groups, parameter estimates were occasionally inaccurate (Fig. 3, Ta-ble 3), and the model occasionally mistook the movement mechanisms drivingthe simulated tracks (Table 2). These outliers may be due to the stochasticityof simulated movement tracks; for example, a resource-only simulated animalmay happen to visit similar portions of its landscape at a coincidentally regularinterval, which the model might mistake for memory-informed movement. Con-versely, an animal following the “memory-only” rules may coincidentally visitlocations that happen to be particularly high (or low) in specific resource values,resulting in the movement track being best explained by the resource-memoryor even the resource-only model.Increasing the number of observed animal locations ( T ) improved our results,but we are more encouraged by the positive effects of simulating additional con-trol points ( K ). While increasing T may require such costly tasks as usinglonger-lasting tracking devices, re-capturing animals and equipping them withnew tracking devices periodically, or increasing the temporal resolution of track-20ng devices, increasing K is easy to do post-hoc. While increasing K may notyield benefits as large as increasing T , the cost of increasing K is much smaller.Our simulated tracks consistently underestimated ρ ns and κ in the resource-only and resource-memory models (Table 3), which is an artifact of the waywe simulated the data. In these models, the animal “chooses” a step from N r proposed steps, which are simulated from k , which depends on ρ ns and κ . Oursimulated landscapes are spatially autocorrelated, so if the simulated animalfound itself in a resource-rich patch, it would be very likely to stay put. Thesemovements are also less directionally autocorrelated than would be suggestedby κ for similar reasons. Using an integrated step selection function (Avgaret al., 2016) could remedy these issues but for our purposes, it adds additionalcomplexity to the model and is not our primary concern.Our estimates of bias and MSE for α did not consistently decrease withincreases in the amount of location data or the number of control points, poten-tially because of an odd bimodal distribution of parameter estimates (Appendix:Fig. S1). The larger portion of this bimodal distribution is clearly centeredaround the true value of approximately -1.78 for all four treatment groups, butcuriously the “second” mode, which appears to be centered around -4.5, seemsto account for more of the estimates T and K increase. These smaller esti-mates for α would imply that the hypothetical organisms moving according toour simulation rules occasionally behave with a much wider understanding oftheir environment, which they perceive to be spatially heterogeneous. The exactcause of these patterns requires further investigation.When we applied the analysis to field data, we notice that the model’s ef-fectiveness, especially when it comes to identifying a memory signal, increasedgreatly with sample size. Our simulations revealed that the model may miss amemory signal with inadequate data, which could explain the disparity betweensubsets of the data in Table 4. While it is possible that grizzly bears, especiallyfemales that take on different reproductive roles in different years, would changetheir movement strategies between years, it is also likely that the model mayhave not had enough data to identify a memory signal in an individual year.With twice data, the simulations accurately identified memory more often, sug-gesting that the memory signals identified for the entirety of each bear’s trackare legitimate.We occasionally observed “false positives” (the model identified memory asa driver of movement from memoryless simulated data) that increased in longeranimal tracks ( K = 1200 vs. K = 600) when simulating tracks with the grizzlybear data. This trend may be an artifact of how the Mackenzie Delta landscapedata influenced simulated tracks, since false positives were much less frequentin the simulation study with artificial landscapes. When comparing this resultwith the real-life subsetting for the bears, we saw examples of subsetted dataregistering memory when the full data set did not, but we also saw examples ofthe opposite.Our modelling framework operates under the assumption that resources varyin time, forcing animals to exhibit seasonal movement patterns within theirhome ranges. We do not explicitly model temporal resource variation, instead21ssuming that our resource covariates r i ( x ) are constant in time. Temporallyvarying resource data are difficult to obtain (for example, they were unavailablefor the Mackenzie Delta region), and transforming our data with some periodicfunction of time (e.g., ˜ r i ( x, t ) = r i ( x ) g ( t )) would either require arbitrary choiceson our part or the addition of parameters to an already complex model. Instead,we model resource variation in time through the animal’s behavior; in a way, weassume that µ is informative about how long resources take to re-appear, andas a result how long animals take to return to them.While we used grizzly bears as a case study, the model was designed tobe general and can be applied to a variety of different systems. Many ani-mals, including turkey vultures ( Cathartes aura ; Holland et al., 2017), blackvultures (
Coragyps atratus ; Holland et al., 2017), caribou (
Rangifer tarandus ;Lafontaine et al., 2017), and eastern indigo snakes (
Drymarchon couperi ; Bauderet al., 2016), perform seasonal movements within their home ranges. For datawith higher temporal resolution, it would be possible to model complex time-dependent recursive movements on a diel scale, since many animals exhibitrepetitive day-to-day movements within their home range (Christiansen et al.,2016; Herbig and Szedlmayer, 2016). Collecting data at finer temporal resolu-tions would be beneficial for inference on memory-informed movement, assum-ing observation errors are accounted for. Even patrolling predators, which weremodeled by Schl¨agel and Lewis (2014), could be modeled using our framework,although we may expect estimates for µ to be smaller than in the grizzly bears.Schl¨agel et al. (2017) displayed the importance of time since last visit for graywolves, but insight on when exactly wolves deem parts of their home range“re-visitable” could be interesting. Of course, migration is also seasonal andpredictable, and although it is typically difficult to obtain environmental datafor an animal’s entire migratory route, spatial memory has been identified asa key driver of migration in many instances (Mueller and Fagan, 2008; Muelleret al., 2011; Fagan et al., 2013; Bracis and Mueller, 2017; Merkle et al., 2019).Fitting this model to migratory populations could provide insights on how toquantify or potentially even predict these mechanisms. Our model uses patterns in animal movement data to obtain information oncomplex time-dependent spatial memory patterns. Made possible by advancedcomputational techniques, we expand on existing literature from animal move-ment modelling as well as animal cognition to generate a model that can beapplied to a variety of ecological systems. The model can estimate the tim-ing of periodic revisitation patterns observed in an animal, which is novel, andalso allows for the interaction of present-time resource selection and memory-informed navigation. We verify our model fitting process using simulated databefore testing its utility on GPS collar data from grizzly bears, finding thatthis very complex model can be effective without need for immense data col-lection. We hope to apply this model more broadly to animals with different22oraging strategies as a means to compare the nature of time-dependent memorymechanisms in different ecological systems.
PRT was funded by the Ashley and Janet Cameron Graduate Scholarship,with support from UAlberta North, as well as the Alberta Graduate Excel-lence Scholarship. PRT was also supported by a University of Alberta Master’sRecruitment Award as well as a University of Alberta Doctoral RecruitmentAward. MAL gratefully acknowledges support from a Canada Research Chairand NSERC Discovery grant. The authors declare no conflict of interest.
PRT designed the analysis with recommendations from MAL, AED, and MAE.MAE and AED collected and supplied grizzly bear data. PRT wrote the draftof the manuscript. All authors have reviewed and provided modifications to themanuscript.
We intend to upload grizzly bear movement data to Dataverse, which is availablethrough the University of Alberta. 23 eferences
Albertsen, C. M., Whoriskey, K., Yurkowski, D., Nielsen, A., and Flemming,J. M. (2015). Fast fitting of non-gaussian state-space models to animal move-ment data via template model builder.
Ecology , 96:2598–2604.Asensio, N., Brockelman, W. Y., Malaivijitnond, S., and Reichard, U. H. (2011).Gibbon travel paths are goal oriented.
Animal Cognition , 14(3):395–405.Auger-M´eth´e, M., Albertsen, C. M., Jonsen, I. D., Derocher, A. E., Lidgard,D. C., Studholme, K. R., Bowen, W. D., Crossin, G. T., and Mills Flemming,J. (2017). Spatiotemporal modelling of marine movement data using templatemodel builder (tmb).
Marine Ecology Progress Series , 565:237–249.Auger-M´eth´e, M., Derocher, A. E., Plank, M. J., Codling, E. A., Lewis, M. A.,and B¨orger, L. (2015). Differentiating the L´evy walk from a composite corre-lated random walk.
Methods in Ecology and Evolution , 6(10):1179–1189.Avgar, T., Deardon, R., and Fryxell, J. M. (2013). An empirically parameter-ized individual based model of animal movement, perception, and memory.
Ecological Modelling , 251:158–172.Avgar, T., Potts, J. R., Lewis, M. A., and Boyce, M. S. (2016). Integrated stepselection analysis: bridging the gap between resource selection and animalmovement.
Methods in Ecology and Evolution , 7(5):619–630.Barker, O. E. and Derocher, A. E. (2009). Brown bear (
Ursus arctos ) predationof broad whitefish (
Coregonus nasus ) in the mackenzie delta region, northwestterritories.
Arctic , 62(3):312–316.Barker, O. E. and Derocher, A. E. (2010). Habitat selection by arctic groundsquirrels (
Spermophilus parryii ). Journal of Mammalogy , 91(5):1251–1260.Bauder, J. M., Breininger, D. R., Bolt, M. R., Legare, M. L., Jenkins, C. L.,Rothermel, B. B., and McGarigal, K. (2016). Seasonal variation in east-ern indigo snake (
Drymarchon couperi ) movement patterns and space use inpeninsular florida at multiple temporal scales.
Herpetologica , 72(3):214–226.Bennett, A. T. (1996). Do animals have cognitive maps?
Journal of Experi-mental Biology , 199:219–224.Berger-Tal, O. and Bar-David, S. (2015). Recursive movement patterns: reviewand synthesis across species.
Ecosphere , 6(9).Boyce, M. S. and McDonald, L. (1999). Relating populations to habitats usingresource selection functions.
TREE , 14(7):268–272.Bracis, C. and Mueller, T. (2017). Memory, not just perception, plays an im-portant role in terrestrial mammalian migration.
Proc Biol Sci , 284(1855).24uderman, F. E., Hooten, M. B., Alldredge, M. W., Hanks, E. M., and Ivan,J. S. (2018). Time-varying predatory behavior is primary predictor of fine-scale movement of wildland-urban cougars.
Mov Ecol , 6:22.Burnham, K. P. and Anderson, D. R. (2004). Multimodel inference.
SociologicalMethods & Research , 33(2):261–304.Charnov, E. (1976). Optimal foraging, the marginal value theorem.
TheoreticalPopulation Biology , 9:129–136.Christiansen, F., Esteban, N., Mortimer, J. A., Dujon, A. M., and Hays, G. C.(2016). Diel and seasonal patterns in activity and home range size of greenturtles on their foraging grounds revealed by extended fastloc-gps tracking.
Marine Biology , 164(1).Clayton, N. S. and Dickinson, A. (1998). Episodic-like memory during cacherecovery by scrub jays.
Nature , 395(6699):272–274.Codling, E. A., Plank, M. J., and Benhamou, S. (2008). Random walk modelsin biology.
J R Soc Interface , 5(25):813–34.Dalziel, B. D., Morales, J. M., and Fryxell, J. M. (2008). Fitting probabilitydistributions to animal movement trajectories: using artificial neural networksto link distance, resources, and memory.
Am Nat , 172(2):248–58.Davies, N. and Houston, A. (1981). Owners and satellites: The economicsof territory defence in the pied wagtail,
Motacilla alba . Journal of AnimalEcology , 50(1):157–180.Duchesne, T., Fortin, D., and Rivest, L. P. (2015). Equivalence between stepselection functions and biased correlated random walks for statistical inferenceon animal movement.
PLoS One , 10(4):e0122947.Edwards, M., Derocher, A. E., Hobson, K., Branigan, M., and Nagy, J. (2011).Fast carnivores and slow herbivores: differential foraging strategies amonggrizzly bears in the canadian arctic.
Oecologia , 165:877–889.Edwards, M. A. and Derocher, A. E. (2015). Mating-related behaviour of grizzlybears inhabiting marginal habitat at the periphery of their north americanrange.
Behavioural Processes , 111:75–83.Edwards, M. A., Nagy, J. A., and Derocher, A. E. (2009). Low site fidelity andhome range drift in a wide-ranging, large arctic omnivore.
Animal Behaviour ,77(1):23–28.Fagan, W. F., Lewis, M. A., Auger-Methe, M., Avgar, T., Benhamou, S., Breed,G., LaDage, L., Schlagel, U. E., Tang, W. W., Papastamatiou, Y. P., Forester,J., and Mueller, T. (2013). Spatial memory and animal movement.
EcologyLetters , 16(10):1316–29. 25ischer, S. M. and Lewis, M. A. (2020). A robust and efficient algorithm to findprofile likelihood confidence intervals.
ArXiv , pages 1–18.Fortin, D., Beyer, H., Boyce, M. S., Smith, D., Duchesne, T., and Mao, J.(2005). Wolves influence elk movements: Behavior shapes a trophic cascadein yellowstone national park.
Ecology , 86(5):1320–1330.Gaillard, J. M., Hebblewhite, M., Loison, A., Fuller, M., Powell, R., Basille, M.,and Van Moorter, B. (2010). Habitat-performance relationships: finding theright metric at a given spatial scale.
Philos Trans R Soc Lond B Biol Sci ,365(1550):2255–65.Gerber, B. D., Hooten, M. B., Peck, C. P., Rice, M. B., Gammonley, J. H.,Apa, A. D., Davis, A. J., and Lemaˆıtre, J. (2019). Extreme site fidelityas an optimal strategy in an unpredictable and homogeneous environment.
Functional Ecology , 33(9):1695–1707.Graham, R. T., Witt, M. J., Castellanos, D. W., Remolina, F., Maxwell, S.,Godley, B. J., and Hawkes, L. A. (2012). Satellite tracking of manta rayshighlights challenges to their conservation.
PLoS One , 7(5):e36834.Harel, R. and Nathan, R. (2018). The characteristic time-scale of perceivedinformation for decision-making: Departure from thermal columns in soaringbirds.
Functional Ecology , 32(8):2065–2072.Herbig, J. L. and Szedlmayer, S. T. (2016). Movement patterns of gray trigger-fish,
Balistes capriscus , around artificial reefs in the northern gulf of mexico.
Fisheries Management and Ecology , 23(5):418–427.Holland, A. E., Byrne, M. E., Bryan, A. L., DeVault, T. L., Rhodes, O. E.,and Beasley, J. C. (2017). Fine-scale assessment of home ranges and activitypatterns for resident black vultures (
Coragyps atratus ) and turkey vultures(
Cathartes aura ). PLoS One , 12(7):e0179819.Janmaat, K. R., Byrne, R. W., and Zuberbuhler, K. (2006). Primates takeweather into account when searching for fruits.
Curr Biol , 16(12):1232–7.Johnson, D. (1980). The comparison of usage and availability measurements forevaluating resource preference.
Ecology , 61(1):65–71.Jonsen, I. D., Basson, M., Bestley, S., Bravington, M. V., Patterson, T. A.,Pedersen, M. W., Thomson, R., Thygesen, U. H., and Wotherspoon, S. J.(2013). State-space models for bio-loggers: A methodological road map.
DeepSea Research Part II: Topical Studies in Oceanography , 88-89:34–46.Kristensen, K., Nielsen, A., Berg, C. W., Skaug, H., and Bell, B. M. (2016).TMB: Automatic differentiation and Laplace approximation.
Journal of Sta-tistical Software , 70(5). 26afontaine, A., Drapeau, P., Fortin, D., and St-Laurent, M. H. (2017). Manyplaces called home: the adaptive value of seasonal adjustments in range fi-delity.
J Anim Ecol , 86(3):624–633.Marchand, P., Garel, M., Bourgoin, G., Duparc, A., Dubray, D., Maillard,D., and Loison, A. (2017). Combining familiarity and landscape featureshelps break down the barriers between movements and home ranges in a non-territorial large herbivore.
J Anim Ecol , 86(2):371–383.Martin-Ordas, G., Haun, D., Colmenares, F., and Call, J. (2010). Keepingtrack of time: evidence for episodic-like memory in great apes.
Anim Cogn ,13(2):331–40.Mauritzen, M., Derocher, A. E., and Wiig, Ø. (2001). Space-use strategies offemale polar bears in a dynamic sea ice habitat.
Canadian Journal of Zoology ,79(9):1704–1713.Merkle, J. A., Fortin, D., and Morales, J. M. (2014). A memory-based foragingtactic reveals an adaptive mechanism for restricted space use.
Ecol Lett ,17(8):924–31.Merkle, J. A., Sawyer, H., Monteith, K. L., Dwinnell, S. P. H., Fralick, G. L., andKauffman, M. J. (2019). Spatial memory shapes migration and its benefits:evidence from a large herbivore.
Ecol Lett , 22(11):1797–1805.Morales, J. M., Haydon, D. T., Frair, J., Holsinger, K. E., and Fryxell, J. M.(2004). Extracting more out of relocation data: Building movement modelsas mixtures of random walks.
Ecology , 85(9):2436–2445.Morey, P. S., Gese, E. M., and Gehrt, S. D. (2007). Spatial and temporal vari-ation in the diet of coyotes in the Chicago metropolitan area.
The AmericanMidland Naturalist , 158(1):147–161.Mueller, T. and Fagan, W. F. (2008). Search and navigation in dynamic en-vironments - from individual behaviors to population distributions.
Oikos ,117:654–664.Mueller, T., Fagan, W. F., and Grimm, V. (2011). Integrating individual searchand navigation behaviors in mechanistic movement models.
Theoretical Ecol-ogy , 4(3):341–355.Nathan, R., Getz, W., Revilla, E., Holyoak, M., Kadmon, R., Saltz, D., andSmouse, P. E. (2008). A movement ecology paradigm for unifying organismalmovement research.
PNAS , 105(49):19052–19059.Normand, E. and Boesch, C. (2009). Sophisticated Euclidean maps in forestchimpanzees.
Animal Behaviour , 77(5):1195–1201.Oliveira-Santos, L. G., Forester, J. D., Piovezan, U., Tomas, W. M., and Fer-nandez, F. A. (2016). Incorporating animal spatial memory in step selectionfunctions.
J Anim Ecol , 85(2):516–24.27an, L. and Wu, L. (1998). A hybrid global optimization method for inverseestimation of hydraulic parameters: Annealing-simplex method.
Water Re-sources Research , 34(9):2261–2269.Papastamatiou, Y. P., Meyer, C. G., Carvalho, F., Dale, J. J., Hutchinson,M. R., and Holland, K. N. (2013). Telemetry and random-walk models revealcomplex patterns of partial migration in a large marine predator.
Ecology ,94(11):2595–2606.Parzen, E. (1961). On estimation of a probability density function and mode.
The Annals of Mathematical Statistics , 33:1065–1076.Prokopenko, C. M., Boyce, M. S., Avgar, T., and Tulloch, A. (2017). Charac-terizing wildlife behavioural responses to roads using integrated step selectionanalysis.
Journal of Applied Ecology , 54(2):470–479.Richter, L., Balkenhol, N., Raab, C., Reinecke, H., Meißner, M., Herzog, S.,Isselstein, J., and Signer, J. (2020). So close and yet so different: the im-portance of considering temporal dynamics to understand habitat selection.
Basic and Applied Ecology , 43:99–109.Schl¨agel, U. E. and Lewis, M. A. (2014). Detecting effects of spatial memoryand dynamic information on animal movement decisions.
Methods in Ecologyand Evolution , 5(11):1236–1246.Schl¨agel, U. E., Merrill, E. H., and Lewis, M. A. (2017). Territory surveillanceand prey management: Wolves keep track of space and time.
Ecol Evol ,7(20):8388–8405.Sciaini, M., Fritsch, M., Scherer, C., Simpkins, C. E., and Golding, N. (2018).Nlmr and landscapetools: An integrated environment for simulating and mod-ifying neutral landscape models in R.
Methods in Ecology and Evolution ,9(11):2240–2248.Shevtsova, A., Ojala, A., Neuvonen, S., Vieno, M., and Haukioja, E. (1995).Growth and reproduction of dwarf shrubs in a subarctic plant community:Annual variation and above-ground interactions with neighbours.
Journal ofEcology , 83(2):263–275.Skaug, H. J. and Fournier, D. A. (2006). Automatic approximation of themarginal likelihood in non-gaussian hierarchical models.
ComputationalStatistics & Data Analysis , 51(2):699–709.Sturz, B. R., Bodily, K. D., and Katz, J. S. (2006). Evidence against integrationof spatial maps in humans.
Anim Cogn , 9(3):207–17.Thurfjell, H., Ciuti, S., and Boyce, M. S. (2014). Applications of step-selectionfunctions in ecology and conservation.
Movement Ecology , 2(4):1–12.28olman, E. C. (1948). Cognitive maps in rats and men.
The PsychologicalReview , 55(4):189–208.Van Moorter, B., Visscher, D., Benhamou, S., B¨orger, L., Boyce, M. S., andGaillard, J.-M. (2009). Memory keeps you at home: a mechanistic model forhome range emergence.
Oikos , 118(5):641–652.Vergara, P. M., Soto, G. E., Moreira-Arce, D., Rodewald, A. D., Meneses,L. O., and Perez-Hernandez, C. G. (2016). Foraging behaviour in Magellanicwoodpeckers is consistent with a multi-scale assessment of tree quality.
PLoSOne , 11(7):e0159096.Whoriskey, K., Auger-Methe, M., Albertsen, C. M., Whoriskey, F. G., Binder,T. R., Krueger, C. C., and Mills Flemming, J. (2017). A hidden Markovmovement model for rapidly identifying behavioral states from animal tracks.
Ecol Evol , 7(7):2112–2121.Wiktander, U., Olsson, O., and Nilsson, S. J. (2001). Seasonal variation in home-range size, and habitat area requirement of the lesser spotted woodpecker(
Dendrocopos minor ) in southern sweden.
Biological Conservation , 100:387–385.Worton, B. (1989). Kernel methods for estimating the utilization distributionin home-range studies.
Ecology , 70(1):164–168.29 K = 10 K = 50N R M RM N R M RM T = N 20 7 2 21 27 4 3 16R 0 24 0 26 0 26 0 24M 3 0 24 23 0 0 25 25RM 0 2 5 43 0 0 2 48 T = N 28 1 1 20 29 2 1 18R 0 35 0 15 0 32 0 18M 0 0 15 35 0 0 20 30RM 0 0 1 49 0 0 0 50Table S1: Breakdown of model selection counts using AIC for different types ofsimulated tracks. The row represents the “true” model that the tracks were sim-ulated from (N = null; R = resource-only; M = memory-only; RM = resource-memory), while the column represents the model that was identified as themost parsimonious explanation of the data using AIC. Treatment groups areidentified by the outer left and upper portions of the table and are separatedby shading. K = 10 K = 50N-R N-M N-RM R-RM M-RM N-R N-M N-RM R-RM M-RM T = N 5 2 26 25 9 3 3 17 18 27R 50 24 50 22 50 50 22 50 20 50M 4 45 46 45 22 3 48 50 50 23RM 50 47 50 47 45 50 49 50 50 48 T = N 0 0 21 24 27 1 1 19 19 24R 50 31 50 13 50 50 29 50 16 50M 2 48 50 50 33 0 50 50 50 28RM 50 50 50 50 49 50 49 50 50 50Table S2: Breakdown of likelihood ratio test results for different types of simu-lated tracks. The row represents the “true” model that the tracks were simulatedfrom (N = null; R = resource-only; M = memory-only; RM = resource-memory),while the column represents the two models that were compared using a likeli-hood ratio test. Counts represent the number of simulated tracks (50 per cell)that registered a p-value below 0.05 for the indicated likelihood ratio test (forexample, the “N-R” column indicates likelihood ratio tests determining whetherthe resource-only model is significantly more explanatory than the null model).Treatment groups are identified by the outer left and upper portions of the tableand are separated by shading. 30igure S1: Violin plot of parameter estimates for α parameter in the resource-memory model for our four treatment groups (detailed on the x-axis), with 50simulations per treatment. The true value of log
10 160 ≈ − ..