Beyond habitat: Memory versus environment in shaping animal space use

For nearly half a century, ecologists have sought to explain animal space use through characteristics of the environment (i.e., habitat). Recent evidence suggests animals also use memory of previous experiences to decide when and where to move. Yet, the relative influence of the two in explaining animal space use has not been resolved. Using six large ungulate species in the Rocky Mountains (USA), we evaluated the performance of a habitat selection model with 16 environmental variables and another with two variables representing previous use (i.e., memory). While memory outperformed the environment for two species and the environment outperformed memory for four species, the influence of memory and the environment was overall comparable. The environment best explained space use of specialists, while memory best explained species with strong site fidelity. Our work challenges traditional habitat selection theory, showcasing that animals build their spatial preferences through experience just as much as merely responding to their environment.

Methods
Animal location data

We compiled GPS-collar datasets from large herbivores throughout Wyoming and the Greater Yellowstone Ecosystem (USA). We included 211 individuals monitored for ≥3 consecutive years with ≥8 GPS locations per day from four mule deer populations, three elk (Cervus canadensis) populations, two populations each of bison (Bison bison), bighorn sheep (Ovis canadensis), and moose (Alces alces), and one pronghorn (Antilocapra americana) population (Figure 1).

Study area

The study area was characterized by long, cold winters and short warm summers. Mean daily maximum temperatures ranged from -4 to 2° C in the winter and 20 to 29° C in the summer (1904 – 2016; https://wrcc.dri.edu). Mean annual precipitation varied widely, with low elevations (the lowest point is approximately 1300 m) receiving only 10-20 cm of precipitation per year, while the higher elevations (the highest point is approximately 4100 m) can receive approximately 10 times as much (150-250 cm; PRISM Climate Group, Oregon State University, http://prism.oregonstate.edu). Sagebrush (Artemisia spp.) dominated low elevations, while herbaceous meadows and mixed forest (lodgepole pine [Pinus contorta], Douglas fir [Pseudotsuga menziesii], and aspen [Populus tremuloides]) dominated median elevations, and high elevations were characterized by herbaceous meadows and mixed conifer forests (spruce [Picea engelmannii], subalpine fir [Abies lasiocarpa], and whitebark pine [Pinus albicaulis]).

Habitat selection framework 

To evaluate the relative importance of environment or previous experience in animal space use, we assessed individual habitat selection within seasonal home ranges using an RSF framework (Boyce et al., 2002; Manly et al., 1993). We assessed whether habitat selection within a seasonal home range can be better explained by the contemporary distribution of environmental characteristics within an animal’s home range or the animal’s previous use of the landscape.

Starting with the third or fourth year of data (depending on whether the animal had a full three or four years of data), we separated the data into 4 seasons. We defined seasons by calendar date (winter = December-February, spring = March-May, summer = June-August, fall = September-November) based on migration periods and weather patterns. We only included data for a given season for a given individual if there were at least 60 days of data. Most collars collected data between 2- and 3-hour intervals. For those collars collecting data more frequently, we randomly subsampled so that those individuals included no more than 12 relocations per day.

For each season of observed GPS data for a given individual, we drew a sample of available points across the home range. The home range was defined as the 99.99% volume contour of a Brownian bridge movement model (Horne et al., 2007) fit to the GPS data using a 50m resolution grid. We used such a large contour to encompass extra available area that the animal likely knew about due to sight or smell. We sampled twice as many available points as used points.

Memory models

We quantified previous use (modifying the methods of Oliveria-Santos et al. 2016) using two metrics: 1) footprint of previous use (binary: whether an animal visited a location), and 2) occurrence distribution (continuous: intensity of use), hereafter referred to as memory models. Our approach primarily captures individual knowledge accumulated through personal experience, though spatial information can also be acquired socially through interactions with conspecifics or transmitted across generations from parents to offspring (Gil et al. 2018, Jesmer et al. 2018). We indexed these previous use variables by season and year in the following different ways: the same season in the previous year, the same season two years prior, the same season three years prior, the same season in the first two years combined, the same season in all three years combined, the previous year’s use only, experience only from 2 years into the past, use only from 2 and 3 years into the past, use from years 1 and 2 combined, use from years 2 and 3 combined, and all previous use combined (i.e., year 1, 2, and 3). For each season and year of data in the past for each individual, we calculated an occurrence distribution using a BBMM using a 50m resolution grid (Horne et al., 2007). The footprint of their previous use was based on the 99% contour of the occurrence distribution. We built separate models representing the 6 to 11 different previous use indices (depending on whether the animal had 2 or 3 years’ worth of previous experience). We allowed for selection for the previous occurrence distribution to be non-linear by including a second-order polynomial term.

While these indices cannot capture the full complexity of spatial memory, they provide tractable measures of how previous space use influences current habitat selection. The previous footprint identifies whether an animal has ever visited a given location, representing spatial familiarity. The previous occurrence distribution quantifies the intensity of use, representing the strength of spatial association. By including both metrics in our models, we can detect whether animals select for or against areas they have previously visited while simultaneously selecting for or against areas of high previous space use, allowing for the detection of more nuanced patterns that may reflect both positive and negative past experiences.

Environment models

We extracted a suite of environmental variables from the used and available locations that have been known to influence the habitat selection of large herbivores. Variables included 16 biotic and abiotic factors such as elevation, distance to roads, land cover (including % cover of shrubs, trees, bare ground), and indices of overall productivity or biomass of perennial and annual forbs and grasses (complete list in Table 1 of the Supplement). We built a single environment model for each id-year-season that included uncorrelated environmental variables and squared terms of those variables (to allow for non-linear responses) with relative empirical support. We began by building a full model (with all 16 variables). We then used a modified version of the iterative variable selection process described by Dunk et al. (2019), first computing variance inflation factors (VIF) for all variables and sequentially removing individual variables with VIF ≥ 4 until all remaining variables exhibited minimal multicollinearity (i.e., VIF < 4). Once a final suite of uncorrelated environmental variables was identified, we iteratively included a squared term for each variable and tested whether the model still converged. Our final model contained all uncorrelated habitat variables and squared terms for those variables that did not cause model convergence issues.

Fitting models and model performance

For each id-year-season, we fit one RSF with environmental variables only, and then 6 to 11 RSFs indexing how the different levels (season vs. annual, and number of years) of previous use influence habitat selection. While we also explored combined models incorporating both environmental variables and previous use variables (results in the Supplement, Fig. S8), our primary analytical framework focused on comparing these components separately to directly assess their relative influence on habitat selection. We fit all RSFs using generalized linear models with a binomial distribution, where the response variable was a used versus available GPS location. We assessed relative empirical support for each RSF using the following framework. We calculated the predicted values from the available locations for each RSF and then partitioned the values into 10 bins based on equally spaced quantiles of the predicted values. Then, for each RSF, we calculated the proportion of used GPS locations falling in the top three predicted bins (bins 7-10, top 30%). In other words, we calculated the proportion of observed locations falling within areas predicted to be the highest relative probability of use. We did not assess relative empirical support using Akaike’s Information Criterion because it is biased toward models with more predictor variables (Link & Barker, 2006; Symonds & Moussalli, 2011), and our habitat models had up to 28 terms, while all of our memory models had only three terms.

Summarizing results across species, seasons, etc.

Our initial summaries include determining: 1) the importance of the amount of time considered as previous experience, 2) if there was a clear distinction in memory model performance between models that indexed previous use in the same season in prior years versus the full prior year(s) of data, and 3) whether memory models performed better in certain seasons (e.g., migratory versus summer or winter). These methods and results can be found in the Supplement.

We also compared the proportion of top models that were memory versus environment and the performance of memory and environment models across species. For these comparisons, we included individuals with up to four years of da, ta but limited each individual to only four seasons (i.e., each individual could only have models for one summer, fall, winter, and spring). First, we compared the proportion of top models that were either environment or memory. To do this, we first identified the top model for each id-year-season (n = 858). We then used linear regression to evaluate if the proportion of top models that were based on previous use was influenced by species, accounting for repeated measures by including id-population-species as a random effect. Next, we compared the performance of the top environment models to the top memory models (n = 858). We evaluated whether model performance was influenced by model type and species using linear regression, including id-year-season-population as a random effect. Finally, we compared the model performance of the environment models with the top memory model for each id-year-season (n = 1716) using the same method described in the previous sentence.

Lastly, we assessed whether general patterns of space use could explain variation in model performance across species. Specifically, we sought to evaluate whether: 1) memory models perform better for species with high versus low site fidelity, and 2) environment models perform better for habitat specialists versus generalists. We implemented a non-parametric, randomization-based statistical approach to identify significant differences in model performance between species. For each model type (memory and environment), we conducted pairwise Fisher-Pitman permutation tests (Hollander and Wolfe 1999) between all possible species pairs (e.g., mule deer vs. elk, mule deer vs. pronghorn, elk vs. pronghorn, etc.) using the coin package (Hothorn et al., 2008) with 100,000 permutations.

For the memory model analysis, we first identified the best-performing memory model for each id-year-season. We then compared the performance of the best memory model between species pairs to determine if species with higher documented site fidelity showed significantly better memory model performance. We applied the Bonferroni correction to adjust p-values for multiple comparisons. We interpreted these memory model results in relation to documented species-specific differences in site fidelity. Previous research by Morrison et al. (2021) reported that mule deer exhibit the highest site fidelity, followed by moose, bighorn sheep, elk, and pronghorn. Although bison were not included in Morrison et al. (2021), we included them in our analysis to provide a comprehensive comparison.

Similarly, for the environment model analysis, we compared the performance of the environment model for each individual between species pairs using the same permutation-based statistical approach. We interpreted these results in the context of habitat specialization, where we considered species with strong preferences for specific habitat types (e.g., bighorn sheep with escape terrain, moose with riparian areas) to be specialists, compared with species that utilize a broader range of habitats, such as elk that can be associated with various elevations and vegetation types.

Extended methods

Data Source: Topographic variables were derived from a 30-meter resolution digital elevation model (DEM). All spatial processing was performed in R using the terra package (Hijmans, 2022).

Derived Variables

Aspect

Aspect represents the compass direction that a slope faces, measured in degrees from north (0°). Aspect was calculated using the terrain() function in the terra package, which implements standard GIS algorithms for deriving terrain attributes from elevation data.

Transformed Aspect (TRASP)

We transformed aspect into a biologically meaningful metric (TRASP), which provides a continuous measure of topographic heat load. TRASP assigns a value of 0 to cool, moist north-northeast facing slopes (typically 30° aspect) and a value of 1 to hot, dry south-southwest facing slopes (typically 210° aspect). TRASP was calculated using the following equation:

TRASP = [1 - cos(π/180 × (aspect - 30))] / 2

This transformation linearizes the circular aspect variable into a measure that better represents the moisture and temperature conditions experienced at different slope exposures.

Slope

Slope gradient was calculated in degrees using the terrain() function in the terra package. This calculation determines the maximum rate of change in elevation between each cell and its neighbors.

Terrain Position Index (TPI)

Terrain Position Index quantifies the difference between the elevation of a cell and the mean elevation of its surrounding cells, effectively identifying topographic positions such as ridges (positive values) and valleys (negative values). We calculated TPI using a 9×9 cell moving window (270m × 270m), which captures landform features at an ecologically relevant scale for our study organisms. TPI was calculated as:

TPI = z₀ - z̄

Where z₀ is the elevation of the focal cell and z̄ is the mean elevation of all cells in the moving window (excluding the focal cell).

Model evaluation technique:

Regarding our binning approach and the choice of top bins, we conducted additional analyses to ensure our results weren't sensitive to this methodological choice. We re-evaluated our models using various thresholds: the top 10% (1 bin), top 20% (2 bins), top 30% (3 bins, current threshold), top 40% (4 bins), and top 50% (5 bins) of available area. This sensitivity analysis revealed that the mean proportions of used locations in top bins were consistent across different thresholds.

When examining which model type (Environment or Memory) performed best across our 858 individual RSFs, the results were fairly similar across thresholds, although we noted that using the top 40% or lower slightly favored environment models. This consistency suggests our findings are robust to the specific threshold chosen.