EXAMINING THE INTRINSIC AND EXTRINSIC DIMENSIONS OF UNGULATE
MOVEMENT AND RESOURCE SELECTION
By
Kyle Martin Redilla

A THESIS
Submitted to
Michigan State University
in partial fulfillment of the requirements
for the degree of
Fisheries and Wildlife – Master of Science
2017

ABSTRACT
EXAMINING THE INTRINSIC AND EXTRINSIC DIMENSIONS OF UNGULATE
MOVEMENT AND RESOURCE SELECTION
By
Kyle Martin Redilla
Large mammalian herbivores are extensively-studied worldwide, often to gain new
insights into the relationship between these animals and their environment. Elucidating a
mechanistic basis of processes such as movement and resource selection can inform conservation
and management practice. Therein, relating these processes to intrinsic (both state- and stagerelated) and extrinsic dimensions (abiotic and biotic characteristics of the environment) is of
paramount importance. I center my thesis on the role of these two dimensions in ungulate
movement and resource selection. In chapter 1, I focused on the extrinsic dimension by
employing both linear and non-linear regression techniques to evaluate the plausibility of
“critical temperatures” in movement of moose (Alces alces) in Norway, using a rich dataset of
GPS-location data. I found weak evidence for these thresholds in the movement of moose, and I
discuss this finding in light of a changing climate. In Chapter 2 I studied the intrinsic dimension
via an examination of individual variation in resource selection of elk (Cervus elaphus) in the
Ozark Mountains of Missouri. I investigated the consequences of prevailing practice, whereby
individual information is pooled to fit an aggregate-level model, by fitting models at the level of
each individual elk and making comparisons. I found that important inferences can be missed if
resource selection is only considered at the aggregate level. In summary, my research
demonstrates the importance of using wild-living individuals and multiple modeling perspectives
to develop functional population- or species-level inferences regarding the roles of intrinsic and
extrinsic factors in animal-environment relationships of ungulates.

ACKNOWLEDGEMENTS

I am deeply indebted to all those who made this research possible, and who believed in
me along the way. I am grateful to my major advisor, Dr. Robert Montgomery, for seeing
potential in me as a directionless undergraduate and granting me the opportunity to conduct
research in the RECaP laboratory. I am thankful to my current thesis committee members, Drs.
Joshua Millspaugh and Jerry Urquhart, for invaluable guidance in forming useful research
questions and interpreting findings in light of established theory, current applied perspectives,
and broader global ecological processes. I would also like to thank Dr. Brian Maurer for advice
on maximizing my potential as a quantitative ecologist, and consultation on my research as well.
I also thank staff members in the department of fisheries and wildlife, who assisted me with the
logistics of my degree, particularly Jill Cruth and Marcia Barr, whose willingness to understand
and answer any of my questions saved me a great deal of time and energy.
I would like to acknowledge my colleagues in RECaP, whose hard work and dedication
to their craft provided an invigorating climate for work, especially Steven Gray, Arthur Muneza,
and Remington Moll, for being role models that I learned a great deal from as I acquired my
footing and adapted to this discipline. I am especially grateful to Rem for his encouragement,
practical advice, and contagious passion for discussing any and all topics related to ecology,
uncertainty, and learning.
I owe a tremendous amount of gratitude to my parents, Martin Redilla and Donna Gee,
for always encouraging me to pursue my passion and giving me everything I ever needed to be
successful.
Finally, I would like to thank my close friends, who made my own self-belief possible
and taught me some crucial life lessons and skills: Karl Supal, you taught me how to grind (an

iii

invaluable skill in the RECaP lab); Seth Rockafellow, you taught me the value of my own
intellect; Tyler Trepanier, you taught me that I can handle more than I think; Sara Tanis, you
taught me how to treat the ones I love; Benjamin Goheen, you taught me a lot about everything,
and you give me hope for humanity; Mitchell Goheen, you taught me how to love everyone, and
you give me hope for finding balance in this life; Molly Robinett, you taught me how to love
myself, and you give me hope for the future.

iv

TABLE OF CONTENTS

LIST OF TABLES ........................................................................................................................ vii
LIST OF FIGURES ..................................................................................................................... viii
INTRODUCTION .......................................................................................................................... 1
REFERENCES ............................................................................................................................ 5
CHAPTER 1 ................................................................................................................................... 9
INVESTIGATING SIGNATURES OF HEAT STRESS IN WILD-LIVING MOOSE ................ 9
Abstract ....................................................................................................................................... 9
1.1.
Introduction .................................................................................................................. 11
1.2.
Methods ........................................................................................................................ 14
1.2.1. Study area................................................................................................................ 14
1.2.2. GPS-tracking and temperature ................................................................................ 15
1.2.3. Movement Metric Calculation ................................................................................ 15
1.2.4. Mixed-effects modeling .......................................................................................... 17
1.2.5. Probabilistic Movement Metrics ............................................................................. 19
1.3.
Results .......................................................................................................................... 20
1.4
Discussion .................................................................................................................... 23
1.5.
Conclusions .................................................................................................................. 27
Acknowledgements ................................................................................................................... 29
APPENDIX ............................................................................................................................... 30
REFERENCES .......................................................................................................................... 46
CHAPTER 2 ................................................................................................................................. 53
EVALUATING THE CONSEQUENCES ON ECOLOGICAL INFERENCE OF
AGGREGATING WILDLIFE TELEMETRY DATA WHEN ESTIMATING RESOURCE
SELECTION FUNCTIONS ......................................................................................................... 53
Abstract ..................................................................................................................................... 53
2.1.
Introduction .................................................................................................................. 55
2.2.
Methods ........................................................................................................................ 58
2.2.1. Resource selection data ........................................................................................... 58
2.2.2. Environmental variables ......................................................................................... 59
2.2.3. RSF modeling ......................................................................................................... 59
2.2.4. Comparison between aggregate- and individual-level RSFs .................................. 62
2.3.
Results .......................................................................................................................... 64
2.3.1. Number and composition of parameters ................................................................. 65
2.3.2. Parameter estimates ................................................................................................ 65
2.3.3. Predicted relative probabilities ............................................................................... 66
2.4.
Discussion .................................................................................................................... 67
Acknowledgements ................................................................................................................... 73
APPENDIX ............................................................................................................................... 74

v

REFERENCES .......................................................................................................................... 80
CONCLUSION ............................................................................................................................. 86
REFERENCES .......................................................................................................................... 89

vi

LIST OF TABLES

Table 1.1. Review of articles published between 2010-2015 that have directly referenced the
temperature thresholds suggested by Renecker and Hudson (1986) to develop models
or discuss climate related impacts on moose. .............................................................. 31
Table 1.2. Counts of temperatures recorded by GPS collars of GPS-tracked moose in central
Norway from July through August of each year during 2006-2010. Counts during
crepuscular period (20:00 – 06:00 next day) above, counts of temperatures recorded
during daytime period (08:00 – 18:00) middle, differences of temperature (ΔT)
between crepuscular and daytime periods below. ........................................................ 33
Table 1.3. Comparison of ARMA error structure for linear mixed models of crepuscular (left
half) and daytime (right half) movements of moose, sorted by ascending AIC score. 34
Table 1.4. Comparison of best approximating LMMs. The top three linear mixed models of
movement rate (Rate) and mean turning angle (Turn) response variables during
crepuscular and daytime movement metrics, sorted by increasing AIC score. ........... 35
Table 1.5. Comparison of best approximating AMMs. The top three additive mixed models of
movement rate (Rate) and mean turning angle (Turn) response variables during
crepuscular and daytime movement metrics, sorted by increasing AIC score. ........... 36
Table 1.6. Best approximating LMMs, daytime hours. Parameter estimates for best
approximating linear mixed models of daytime movements, as determined by AIC.
Max. daytime collar temperature (Tday) was the temperature covariate used in these
models. ......................................................................................................................... 37
Table 1.7. Best approximating LMMs, crepuscular hours. Parameter estimates for best
approximating linear mixed models of crepuscular movements, as determined by AIC.
The temperature covariate used in each model is indicated in parentheses. Tcrep = max.
crepuscular collar temperature, Tday = max. daytime collar temperature, and ΔT = Tday
– Tcrep. ........................................................................................................................... 38
Table 2.1. WAIC values for top 10 aggregate-level resource selection function models. +
indicates an environmental variable included in the model, spaces left blank if not
included. Environmental variable abbreviations as follows: D2t = distance to twotrack road; Dpav = Distance to paved road (m); Slope = slope (degrees); Aspect =
aspect (degrees); Rd dens = road density (km/km2); D edge = Distance to nearest
wooded edge; % Can = percent forested canopy cover; Rx burn = years since
prescribed fire; Habitat = indicator variables for 8 cover types. .................................. 75

vii

LIST OF FIGURES

Figure 1.1. Summer locations of GPS-tracked moose in Central Norway during 2006-2010.
Norway is shaded light blue, other countries colored olive. ....................................... 39
Figure 1.2. Mean ln-transformed speeds (Âą SE) of moose movement in Central Norway, during
the months of July and August of each year from 2006-2010 (n = 152 moose). Speeds
were estimated from continuous-time correlated random walk movement models fit to
each summer track using the package crawl in R. ...................................................... 40
Figure 1.3. Boxplots of daytime (a) and crepuscular (b) movement metrics calculated from GPS
data on moose tracked during July-August of 2006-2010. Movement metrics are
binned by maximum collar temperature recorded during the respective period for
which the movement metric was calculated. .............................................................. 41
Figure 1.4. Predictions based on best approximating AMMs (Âą 95% CI shaded grey) for
movement rates of male moose at median elevation during daytime (a) and
crepuscular periods (b-d). Daytime movement rates predicted as a function of daytime
temperature, and crepuscular movement rates predicted as a function of crepuscular
temperature (b), daytime temperature (c), and the difference between daytime and
crepuscular temperature (d). ....................................................................................... 42
Figure 1.5. Predictions based on best approximating AMMs (Âą 95% CI shaded grey) for mean
turning angles of male moose at median elevation during daytime (a) and crepuscular
periods (b-d). Daytime turning angles predicted as a function of daytime temperature,
and crepuscular turn angles predicted as a function of crepuscular temperature (b),
daytime temperature (c), and the difference between daytime and crepuscular
temperature (d). ........................................................................................................... 43
Figure 1.6. Predictions based on best approximating AMMs of movement rates calculated from
probabilistic movement paths. Predictions for daytime (a) and crepuscular (b-d)
movement rates of male moose. Daytime movement rates predicted as a function of
daytime temperature, and crepuscular movement rates predicted as a function of
crepuscular temperature (b), daytime temperature (c), and the difference between
daytime and crepuscular temperature (d). ................................................................... 44
Figure 1.7. Predictions based on best approximating AMMs of mean turning angles calculated
from probabilistic movement paths. Predictions for daytime (a) and crepuscular (b-d)
movement rates of male moose. Daytime movement rates predicted as a function of
daytime temperature, and crepuscular movement rates predicted as a function of
crepuscular temperature (b), daytime temperature (c), and the difference between
daytime and crepuscular temperature (d). ................................................................... 45
Figure 2.1. Inclusion of parameters at the individual level in top models (a panel) and top 5% of
models (b panel) as ranked by WAIC for elk reintroduced into the Missouri Ozarks

viii

(2011-2014). Green boxes (panel a) indicate the selection parameters in line with the
row was included in the top model as ranked by WAIC. Shade of purple (panel b) is
proportional to the number of times the corresponding parameter was included in the
top 5% of sub-models, where the darkest shade corresponds to all models and lightest
shade indicates that parameter was not included in any of the top models.. .............. 76
Figure 2.2. Aggregate-level random effects distributions (boxes) and individual-level RSF
parameter estimates (blue dots) with 95% Bayesian credible intervals (CIs, in
translucent blue) for reintroduced elk in the Missouri Ozarks (2011-2014). The
middle horizontal bars within the boxes are the point estimates of the mean
hyperparameters of the random effects distributions. The light gray horizontal bars
within the boxes are point estimates for the standard deviation hyperparmeters, and
the ends of the boxes represent 95% Bayesian credible intervals on mean + standard
deviation point estimates............................................................................................. 77
Figure 2.3. Principal components analysis (PCA) ordination of the 88 individual elk reintroduced
into the Missouri Ozarks (2011-2014) based on individual-level parameter estimates
for selection of 9 resource covariates. Grouped by intrinsic factors (panel a) and by
GPS-collaring and release year (cohort; panel b). ...................................................... 78
Figure 2.4. Predicted relative probabilities of use of the elk restoration zone based on and
aggregate-level RSF (panel a) and individual-level RSFs for 20 randomly selected elk
(panel b) reintroduced to the Missouri Ozarks (2011-2014). ..................................... 79

ix

INTRODUCTION

The collective ability to measure and map various ecological phenomena has expanded
dramatically over the last two decades, as evidenced by a myriad of novel data types and
volumes, ushering in a new era of inquiry that has been likened to the bioinformatics movement
of the early 2000’s (Jones et al., 2006; Michener and Jones, 2012). Unprecedented rates of data
volume and data acquisition have now become commonplace due to advancements in technology
and networking, leading to a “data deluge” (Baraniuk, 2011). Novel data collection methods
including sensor networks, such as the National Ecological Observatory Network (NEON),
automated sensing (including camera trap technology), and both aircraft- and satellite-borne
remote sensing platforms are now regularly used to collect environmental observations that have
great utility for ecological research (Keller et al., 2008; Pfeifer et al., 2012; Porter et al., 2009,
2012). Not surprisingly, developments that facilitate the quantitative analysis of big ecological
data have emerged in tandem, from best practices in data management to statistical methods and
software tools (Borer et al., 2009; Purves et al., 2013; Steiniger et al., 2009).
These advancements have altered the trajectory of scientific discourse in virtually all
types of applied ecological research, but one area of particular growth has been animal
movement ecology. This previously marginal sub-discipline only received cursory attention until
recent years, when a paradigm shift has been facilitated by decreasing costs and increasing
flexibility of tracking methods, most notably those based on Global Positioning System (GPS)
technologies (Hussey et al., 2015; Kays et al., 2015; Nathan et al., 2008). A central theme of this
new paradigm is a shift in focus from Eulerian (i.e. place- or population-based) approaches to
Lagrangian (i.e. individual-based) approaches for studying movement and space use (Nathan et
al., 2008; Smouse et al., 2010). When combined with data derived from sensor networks and

1

remote sensing systems, high-resolution tracking data holds promise for revealing mechanistic
connections between individual animals and their environment (Cooke et al., 2004; Schick et al.,
2008). Recent advancements in Global Positioning System (GPS) telemetry technology provide a
venue to study the mechanisms underlying fundamental concepts in animal ecology, such as
home range estimation and resource use (Cagnacci et al., 2010).
A mechanistic understanding of species movement through human-dominated landscapes
is of critical importance to both conservation and management (Fleishman et al., 2011). Indeed,
animal-borne telemetry and Geographic Information System (GIS) software have been
indispensable tools to wildlife biologists and managers for at least 25 years (Millspaugh and
Marzluff, 2001). Barriers to a mechanistic mode of inquiry have traditionally been driven by lack
of data however, reflecting the advancements outlined above, animal movement data has recently
become voluminous and ubiquitous enough to earn the “big” data title (Kays et al., 2015; Urbano
et al., 2010). Databases specific to animal movement data have been created to archive data and
facilitate collaboration. For example, as of January 2017, the open-access repository known as
Movebank included over 387 million locations on 635 taxa (Wikelski and Kays, 2017). Coupled
with the growth across these data collection/management spheres, there has been ongoing
development of statistical and software tools necessary to accommodate these increasingly
resolute data (Calabrese et al., 2016; Gurarie et al., 2016). Many methods have been proposed for
modelling the behavioral mechanisms of animal movement using both Fisherian and Bayesian
techniques, but the Bayesian framework offers an intuitive way to accommodate uncertainty of
latent states in these types of models (Calabrese et al., 2016; Clark, 2005; Gurarie et al., 2016;
Hobbs and Hooten, 2015; Johnson et al., 2008).

2

Studies of ungulate ecology can particularly benefit from sophisticated analyses of
telemetry and remotely-sensed data. Movements and behaviors of these animals typically occur
over vast spatiotemporal scales, and their relatively large body size permits the attachment of
tracking devices that can collect precise data over long periods of time (Kays et al., 2015). Given
the ecological and cultural significance of many ungulate species across the globe, there is
growing concern around the ways in which these species will respond to continued
anthropogenic disturbance and changing climate patterns. Herein, my thesis investigates
questions in ungulate ecology arising from the complex interaction between individuals, their
movement and resource selection behaviors, and the intrinsic and extrinsic factors affecting these
processes. In chapter 1, I explored the relationship between Scandinavian moose movements and
a changing climate, specifically rising ambient temperatures. Using five years’ worth of
relocation data from over 150 GPS-tracked moose in Norway, I derived movement metrics over
daytime and crepuscular periods and employed linear and non-linear regression techniques to test
established theory pertaining to the “critical temperatures” beyond which moose behavior is
expected to change considerably (Renecker and Hudson, 1986). In discussing this analysis, I
address the support for these established thresholds and the importance of using free-roaming
animals to delineate threshold responses to extrinsic factors. In chapter 2, I explored the
consequences of analyzing resource selection by elk in an aggregate manner, as is conventionally
done to achieve statistical inference at the population level. With four years of relocation data on
elk introduced into southern Missouri, along with associated GIS layers mapping environmental
features in the landscape, I fit discrete-choice models at the aggregate and individual levels and
compared inferences between these two broad methodologies. I report on the significance of
discounting individual variation in this key ecological process when informing conservation and

3

management. My thesis examines both intrinsic and extrinsic factors on ungulate behavior and in
this way, provides a novel computational and analytical framework for assessing movement and
resource selection using Bayesian statistical methods.

4

REFERENCES

5

REFERENCES

1) Baraniuk, R. G. 2011. More is less: signal processing and the data deluge. Science, 331
(6018):717 – 719.
2) Borer, E. T., Seabloom, E. W., Jones, M. B., Schildhauer, M. 2009. Some simple guidelines
for effective data management. Bulletin of the Ecological Society of America, 90 (2):205 –
214.
3) Cagnacci, F., Boitani, L., Powell, R. A., Boyce, M. S. 2010. Animal ecology meets GPSbased radiotelemetry: a perfect storm of opportunities and challenges. Philosophical
Transactions of the Royal Society: Biological Sciences, 365 (1550):2157 – 2162.
DOI:10.1098/rstb.2010.0107.
4) Calabrese, J. M., Flemming, C. H., Gurarie, E. 2016. Ctmm: an R package for analyzing
animal relocation data as a continuous-time stochastic process. Methods in Ecology and
Evolution, 7 (9):1124 – 1132. DOI:10.1111/2041-210X.12599.
5) Clark, J. 2005. Why ecologists are becoming Bayesians. Ecology Letters, 8 (1):2 – 14.
DOI:10.1111/j.1461-0248.2004.00702.x.
6) Cooke, S. J., Hinch, S. G., Wikelski, M., Andrews, R. D., Kuchel, L. J., Wolcott, T. G.,
Butler, P. J. 2004. Biotelemetry: a mechanistic approach to ecology. Trends in Ecology and
Evolution, 19 (6):334 – 343.
7) Fleishman, E., Blockstein, D. E., Hall, J. A., Mascia, M. B., Rudd, M. A., Scott, J. M.,
Sutherland, W. J., Bartuska, A. M., Brown, A. G., Christen, C. A., Clement, J. P., Dellasala,
D., Duke, C. S., Eaton, M., Fiske, S. J., Gosnell, H., Haney, J. C., Hutchins, M., Klein, M. L.,
Marqusee, J., Noon, B. R., Nordgren, J. R., Orbuch, P. M., Powell, J., Quarles, S. P.,
Saterson, K. A., Savitt, C. C., Stein, B. A., Webster, M. S., Vedder, A. 2011. Top 40
priorities for science to inform US conservation and management policy. Bioscience, 61
(4):290 – 300.
8) Gurarie, E. 2016. What is the animal doing? Tools for exploring behavioral structure in
animal movements. Journal of Animal Ecology, 85 (1):69 – 84.
9) Hobbs, N. T., Hooten, M. B. 2015. Bayesian models: a statistical primer for ecologists.
Princeton, NJ: Princeton University Press. ISBN: 9780691159287.
10) Hussey, N. E., Kessel, S. T., Aarestrup, K., Cooke, S. J., Cowley, P. D., Fisk, A. T.,
Harcourt, R. G., Holland, K. M., Iverson, S. J., Koelk, J. F., Flemming, J. E. M., Whoriskey,
F. G. 2015. Aquatic animal telemetry: A panoramic window into the underwater world.
Science, 348 (6240):1255642-1 – 1255642-10. DOI:10.1126/science.1255642.

6

11) Johnson, D. S., London, J. M., Lea, M.-A., Durban, J. W. 2008. Continuous-time correlated
random walk model for animal telemetry data. Ecology, 89 (5):1208 – 1215.
12) Jones, M. B., Schildhauer, M. P., Reichman, O. J., Bowers, S. 2006. The new bioinformatics:
Integrating ecological data from the gene to the biosphere. Annual Review of Ecology,
Evolution, and Systematics, 37 (1):519 – 544.
13) Kays, R., Crofoot, M. C., Jetz, W., Wikelski, M. 2015. Terrestrial animal tracking as an eye
on life and planet. Science, 348 (6240):aaa2478-1 – aaa2478-9.
DOI:10.1126/science.aaa2478.
14) Keller, M., Schinel, D. S., Hargrove, W. W., Hoffman, F. W. 2008. A continental strategy for
the National Ecological Observatory Network. Frontiers in Ecology and the Environment, 6
(5):282 – 284.
15) Michener, W. K., Jones, M. B. 2012. Ecoinformatics: supporting ecology as a data-intensive
science. Trends in Ecology and Evolution, 27 (2):88 – 93.
16) Millspaugh, J. J., Marzluff, J. M. 2001. Radio-tracking and animal populations: past trends
and future needs In J. J. Millspaugh and J. M. Marzluff (Eds.), Radio tracking and animal
populations (pp. 383 – 393). Cambridge, MA: Academic Press. DOI:10.1016/B978012497781-5/50016-5.
17) Nathan, R., Getz, W. M., Revilla, E., Holyoak, M., Kadmon, R., Saltz, D., Smouse, P. E.
2008. A movement ecology paradigm for unifying organismal movement research.
Proceedings of the National Academy of Sciences of the USA, 105 (49):19052 – 19059.
DOI:10.1073/pnas.0800375105.
18) Pfeifer, M., Disney, M., Quaife, T., Marchant, R. 2012. Terrestrial ecosystems from space: A
review of earth observation products for macroecology applications. Global Ecology and
Biogeography, 21 (6):603 – 624.
19) Porter, J. H., Hanson, P. C., Lin, C.-C. 2012. Staying afloat in the sensor data deluge. Trend
in Ecology and Evolution, 27 (2):121 – 129.
20) Porter, J. H., Nagy, E., Kratz, T. K., Hanson, P., Collins, S. L., Arzberger, P. 2009. New eyes
on the world: Advanced sensors for ecology. Bioscience, 59 (5):385 – 397.
21) Purves, D., Scharlemann, J., Harfoot, M., Newbold, T., Tittensor, D. P., Hutton, J., Emmott,
S. 2013. Time to model all life on earth. Nature, 493 (7432):295 – 297.
22) Renecker, L. J., Hudson, R. J. 1986. Sesonal energy expenditures and thermoregulatory
responses of moose. Canadian Journal of Zoology, 64 (3):322 – 327.
23) Schick, R. S., Loarie, S. R., Colchero, F., Best, B. D., Boustany, A., Conde, D. A., Halpin, P.
N., Joppa, L. N., McClellan, C. M., Clark, J. S. 2008. Understanding movement data and

7

movement processes: current and emerging directions. Ecology Letters, 11 (12):1338 – 1350.
DOI:10.1111/j.1461-0248.2008.01249.x.
24) Smouse, P. E., Focardi, S., Moorcroft, P. R., Kie, J. G., Forester, J. D., Morales, J. M. 2010.
Stochastic modelling of animal movement. Philosophical Transactions of the Royal Society:
Biological Sciences, 365 (!550):2201 – 2211.
25) Steiniger, S., Hay, G. J. 2009. Free and open source geographic information tools for
landscape ecology. Ecological Informatics, 4 (4):183 – 195.
26) Urbano, F., Cagnacci, F., Calenge, C., Dettki, H., Cameron, A., Neteler, M. 2010. Wildlife
tracking data management: a new vision. Philosophical Transactions of the Royal Society:
Biological Sciences, 365 (1550): 2177 – 2185. DOI:10.1098/rstb.2010.0081.
27) Wikelski, M., Kays, R. 2017. Movebank: archive, analysis and sharing of animal movement
data. Hosted by the Max Planck Institute for Ornithology. www.movebank.org, accessed on
1 March 2017.

8

CHAPTER 1

INVESTIGATING SIGNATURES OF HEAT STRESS IN WILD-LIVING MOOSE

Abstract
Heat stress from rising ambient temperatures associated with global climate change
threatens the persistence of numerous species of wildlife as heat-induced behavioral changes can
decrease fitness and alter population dynamics. Moose (Alces alces) are purported to be heatsensitive mammals that become stressed when ambient temperatures rise above specific
temperature thresholds. However, these temperature thresholds were established in a study which
linked physiological changes in thermoregulation to ambient temperatures in a small study group
(n = 2) of captive-reared moose. Temperature thresholds for organisms vulnerable to climate
change should be rigorously evaluated from wild-living animals whenever possible. Here, I
examined the effects of ambient temperature on movements of wild-living moose to investigate
potential temperature thresholds beyond which movement is altered. Analyzing animal
movement in this manner provides a natural way to relate environmental conditions to behavior.
I fit both linear and additive mixed models to daily movement metrics (movement rate and mean
turning angle) calculated from GPS locations collected on 152 moose tracked in central Norway
during 2006-2010. I modeled these movement metrics as a function of temperature recorded by
GPS collars during different periods of the day, as well as sex, mean elevation, and month. I
found little evidence of a threshold response in moose movement to temperature during either
period, but note responses to high within-day temperature differences. I documented differences
between crepuscular and daytime movement responses to temperature, suggesting that moose

9

exhibit a compensatory strategy for mitigating the effects of high temperatures by moving more
during the crepuscular and night hours when temperatures are high. These results suggest that
onset of problematic heat stress in wild-living moose may not occur at the previously proposed
temperature thresholds, or that heat stress mitigation strategies are not well-correlated with the
movements considered. Additionally, within-day variability in ambient temperature may be an
important factor to consider when determining the onset of heat stress in moose.

10

1.1.

Introduction

Warming temperatures associated with global climate change have been causally linked
to changes in the behavior, reproduction, distribution, and abundance of a variety of species
(Harley et al., 2006; Humphries et al., 2004; Parmesan and Yohe, 2003; Root et al., 2003;
Thomas et al., 2004). Moose (Alces alces) are temperate-zone obligates that are well-adapted to
the cold and can tolerate temperatures as low as -32 °C (Karns, 2007; Renecker and Hudson,
1986), but unseasonably warm temperatures purportedly induce physiological heat stress
(Renecker and Hudson, 1986; Schwartz and Renecker, 2007). In response to heat stress, moose
are expected to engage in thermoregulatory behaviors such as the selection of habitat providing
thermal refuge, and compensatory activity schedules (Broders et al., 2012; Dussault et al., 2004).
This heat sensitivity supports the notion that a warming climate is linked with recent (within the
last three decades) declines in moose populations at the southern extent of the species’ range
globally (Dou et al., 2013; Lenarz et al., 2009; Murray et al., 2006). Other possible drivers
implicated in these declines include pathogens and underestimated effects of wolf (Canis lupus)
predation (Lankester, 2010; Mech and Fieburg, 2014; Murray et al., 2006). While temperature is
likely a factor limiting the southern geographic range of this species, and continued trends in
climate change may affect moose survivability in certain populations, the precise effect of
warming temperatures on population viability remains poorly understood (Lenarz et al., 2010;
van Beest et al., 2012). For example, a viable population of moose with relatively ideal
demographic characteristics resides in southern Ontario despite exposure to apparently
unfavorable climatic conditions (Murray et al., 2012). Identifying the mechanism(s) by which

11

rising ambient temperature affects individual animal fitness is crucial to assessing the risks that
changing climate patterns may pose to populations.
Studying animal movement is a useful way to examine climate-related impacts for many
species (McKellar et al., 2005; Nathan and Giuggioli, 2013; Patterson et al., 2009; Schick et al.,
2008). In the case of ungulates, temperature and precipitation have been linked to movements of
white-tailed deer (Odocoileus virginianus; Brinkman et al., 2005) and wild boar (Sus scrofa;
Thurfjell et al., 2014). In addition, fine-scale movements of migratory caribou (Rangifer
tarandus caribou) in the Arctic can be affected by availability of ice (Leblond et al., 2016), and
the interaction of plant phenology and climate can affect large-scale movements of red deer
(Cervus elaphus; Pettorelli). Exploring the movement-temperature relationship in moose
provides an ideal opportunity to investigate potential fitness impacts of rising temperatures, for
two reasons: 1) there are physiological temperature thresholds that have been proposed for
moose (McCan et al., 2013; Renecker and Hudson, 1986); 2) there is support for altered behavior
in response to ambient temperatures near these thresholds (Melin et al., 2014; Street et al., 2015;
van Beest et al., 2012, 2013).
It has been suggested that moose are stressed by heat when ambient temperatures rise
above 14° C (57° F) in summer and above -5° C (23° F) in winter (Renecker and Hudson, 1986,
1990). These seasonal thresholds however, were derived from a study of 2 captive-reared cow
moose kept in a provisioned fenced enclosure (Renecker and Hudson, 1986). Evidence of heat
stress in this study was based on physiological responses, namely increased respiratory rate and
open-mouth panting. The authors postulated that when moose are stressed by heat they reduce
both activity and food intake (Renecker and Hudson, 1986). Despite the narrow scope of this
assessment, these thresholds have been referenced scores of times since 1986. For example, > 15

12

papers published between 2010 and 2015 in 11 different journals studying 10 moose populations
worldwide have either used these thresholds to develop models or cited them in discussing
climate-related impacts on moose (see Table 1.1). However, because the exact mechanism by
which rising ambient temperatures could have population-level consequences for moose still
remains unclear, the utility of the thresholds presented by Renecker and Hudson (1986) has been
questioned (Lowe et al., 2010; McGraw et al., 2012).
One possible mechanistic explanation is that when ambient temperature rises above some
seasonal thresholds, moose movement would be compromised to the point where foraging effort
would be reduced. If these effects were sustained because of heat stress, weight loss would be
incurred which, in turn, would deteriorate body condition potentially leading to a reduction in
fecundity and fitness. Despite support of this explanation as a conceptual model in the literature,
the first connection (altered movement and activity beyond these thresholds) has received only
marginal support (DelGuidice et al., 2011; Sand, 1996; Testa and Adams, 1998). For instance,
while there is support for temperature-induced changes in moose movement behavior (resource
selection and habitat use) at temperatures near the proposed thresholds (Broders et al., 2012;
Melin et al., 2014; van Beest et al., 2012, 2013), another study found no difference in such
behavior of moose (e.g., habitat use) with respect to these temperature thresholds (Lowe et al.,
2010). More specific changes in wild-living moose movement behavior observed during warm
periods include reduced activity, decreased distance traveled, and selection for habitat providing
thermal refugia (Demarchi and Bunnell, 1995; Dussault et al., 2004; Street et al., 2015). A
compensatory effect has also been suggested whereby, during warm periods, moose shift the
bulk of their activity to night-time hours when conditions are comparatively cooler (Dussault et
al., 2004). That being said, no mechanistic link has connected variation in movement to

13

population-level consequences given the thresholds suggested by Renecker and Hudson (1986).
Therefore, while the Renecker and Hudson thresholds have been often applied by the research
community as benchmarks for understanding moose response to heat stress (Renecker and
Hudson, 1986; McCan et al., 2013; Table 1.1), their adequacy for describing population-level
processes in free-living moose is blurred by inconsistencies in behavioral responses to ambient
temperature and possible compensatory behaviors during warm periods.
Here, I analyzed a large set of movement data from a population of moose in central
Norway to answer a fundamental research question: are there thresholds in ambient temperature
above which wild-living moose movement changes significantly? Specifically, I searched for
changes in daytime and crepuscular movements in relation to temperature during these periods. I
also investigated the possibility of compensatory movement schedules by relating both daytime
temperature and within-day temperature difference to crepuscular movements. Analyzing animal
movement in this manner provides a framework to assess the mechanistic role that climate
change phenomena may have for species population trajectories. I discuss the implications of this
research for understanding how warming temperatures might affect moose fitness, and I examine
the utility of my approach for detecting climate-induced changes in moose behavior from
movement paths and associated environmental data.
1.2.

Methods
1.2.1. Study area
Through a partnership with the Norwegian Institute for Nature Research (NINA), I

obtained data on moose that were tracked across central Norway (64°32’ N, 12°15 E), in NordTrøndelag, Sør-Trøndelag, and southern parts of Nordland counties from 2006 to 2010 (Fig. 1.1).
The study area ranged from coastal regions to alpine zones containing mountain, boreal forest,

14

and cultivated land biomes (Van Moorter et al., 2013). The majority of the study area was
forested, consisting primarily of coniferous stands and to a lesser extent deciduous forests
(Bjørneraas et al., 2010). Cultivated land was particularly common at lower elevations, along the
fjords on the west coast (Moen, 1999).
1.2.2. GPS-tracking and temperature
During February-March or November 2006, and February-March 2007 and 2008, moose
were darted from a helicopter to immobilize them, and a sedative administered to facilitate
attachment of GPS collars (Bunnefeld et al., 2011). 171 moose (145 adults, 26 calves) were fit
with GPS/Global System for Mobile communications (GSM) collars, 7 of which were Tellus
GPS collars (Followit AB/Televilt, Lindesburg, Sweden) and the remaining 164 were GPS
PLUS/GPS PRO Light collars (VECTRONIC Aerospace GmbH, Berlin, Germany). 50 moose
(37 F, 13 M) were collared in 2006, 87 moose (63 F, 24 M) were collared in 2007, and 31 moose
(17 F, 14 M) were collared in 2008. Locations were recorded during 2009 and 2010 as well,
tracking 67 moose (54 F, 13 M) in 2009 and 25 moose (22 F, 3 M) in 2010. The collars had an
expected battery life of 3 years while attempting a locational fix every 2 hours (Bjørneraas et al.,
2011; Bunnefeld et al., 2011). These GPS collars recorded ambient temperature from an onboard sensor at each successful location fix. Although these temperatures might be biased by
factors such as prevailing conditions or placement on the body, the offset tends to be consistent
and collar temperatures provide a reliable index of ambient air temperature (Ericsson et al.,
2015).
1.2.3. Movement Metric Calculation
Using the data acquired from these efforts, I calculated moose movements from 214,024
GPS locations (107,049 daytime fixes, 106,975 crepuscular fixes) of 152 moose in central

15

Norway recorded during the summer months (July and August) of 2006-2010. I centered my
assessment on summer because this is the time of year when behavior-altering heat stress in
moose is suggested to be most obvious (Broders et al., 2012; Dussault et al., 2004), when
migratory individuals in this region are likely to be occupying summer home ranges, and to
minimize effects of parturition and rearing that occur during the early summer (Bjørneraas et al.,
2011; Rolandsen et al., 2016; Solberg et al., 2007). After screening the GPS data for erroneous
fixes following established procedures, I subdivided the data so as to assess variation in moose
movement response to temperature between daytime and crepuscular periods (Bjørneraas et al.,
2010). Crepuscular periods of each day were defined as the time from 8 p.m. to 6 a.m. the
following day, so as to encapsulate the movement path between both onset of civil twilight at
night and offset of civil twilight in the morning during July and August. This approach allowed
me to assess the effects of temperature on moose movements when activity/foraging effort is
typically at the highest, and it is consistent with daily activity patterns of other large herbivores
in temperate regions (Cederlund, 1989; Ensing et al., 2014). I defined the daytime period of each
day as the hours between 8 a.m. and 6 p.m., given the remainder of the recorded locations fell
within these hours. I calculated two movement metrics during both the daytime and crepuscular
periods of each summer day: i) movement rate - calculated as the distance travelled (the sum of
the inter-location distances) divided by the length of the period (10 hours), and ii) mean turning
angle - where the turning angle of a location is the angle formed between the line extending
beyond that location in the direction of travel from the previous point, and the line segment
between the current location and the next. I considered angles to be constrained by 180°, such
that a clockwise angle of 270° would be equivalent to a counter-clockwise angle of 90° in the
opposite direction.

16

1.2.4. Mixed-effects modeling
I modeled movement rates and mean turning angles as a function of temperature in a
linear mixed model (LMM) framework for movements during each period. In models of daytime
movements, I used the maximum collar temperature recorded (Tday) during the respective
daytime period. To directly assess crepuscular movements in relation to ambient temperature, I
used maximum collar temperature recorded during the respective crepuscular period (Tcrep) as the
temperature covariate. I used the maximum recorded temperatures because I was interested in
moose response to specific thresholds rather than average conditions; however, a preliminary
analysis using mean daily temperatures estimated by local weather stations did not differ
qualitatively from the results presented here. Additionally, to test the possibility of compensatory
movements during warm periods, I modeled crepuscular movements in two other cases, each
utilizing a different temperature covariate: Tday of the prior daytime period, and ΔT, the
difference between maximum collar temperatures recorded during the daytime prior and the
current period (ΔT = Tday – Tcrep). I included individual moose (id) and year within id as random
intercepts. Year was modeled as a random effect within id because these factors were only
partially crossed, meaning the variation in movements between years is confounded by the
differences in movements between the different groups of moose tracked each year (Pinheiro and
Bates, 2000). The remaining fixed effects for each model were sex; elevation, calculated by
taking the average of the digital elevation model estimates (resolution 25 m, altitude 1 m) at each
location recorded during the corresponding period; and month. Thus, the global LMM was of the
form:
Yijk = β0 + βtempx1ijk + βelevx2ijk + βmox3ijk + βsexx4ijk + u1j + u2jk + ξijk

17

(1)

where Yijk is the ith observation of any of the three movement metrics calculated for the jth
individual in year k; β-parameters correspond to population mean response and effects of
temperature, elevation, month (base July), and sex (base female); u1-2 correspond to normallydistributed random effects of moose id j and year k within id j; and Îľijk is a temporallyautocorrelated error term. I included both sex and month in the model because of differences
observed in preliminary analysis of within-day variation in movement rate (Fig. 1.2), and
because males may regularly travel greater distances than females (Bjørneraas et al., 2012; Van
Moorter et al., 2013). Elevation was included in the models because altitude is one of the main
environmental gradients in this region, and preliminary analyses revealed a large variation across
the study population relative to variation within individuals (Bakkestuen et al., 2008). I tested
temporal autocorrelation structures, prior to fitting LMMs, by fitting the global linear mixed
models for each movement metric with an autoregressive (AR) moving-average (MA) error
structure. I separately fit all 8 combinations of orders 0-2 for both the AR and MA components,
as well as AR(0) MA(3) and AR(3) MA(0). Models were compared using AIC, ranking models
according to ascending AIC value. I conducted these tests separately for both crepuscular and
daytime movements, and selected the AR-MA structure having the lowest AIC value for each.
To ensure reliable estimates of collar temperature and elevation were available for a given
period, I only included crepuscular and daytime periods having at least 5 (out of 6 possible)
recorded GPS locations. Although other factors influence fine-scale movements of moose I
focused on examining the relationship between moose movements and temperature while
minimizing model complexity (Leblond et al., 2010).
In addition to the linear mixed models, I also wanted to test for the potential of nonlinear
relationships in these data. Specifically, I hypothesized that the relationship between movement

18

metrics and temperature measures could be curvilinear. Therefore, I evaluated non-linear
patterns by developing an additive mixed model (AMM) framework. Additive models are
especially useful for this analysis because they would help reveal a specific temperature
threshold at which movements significantly change, as might be expected if movement is
affected at the Renecker and Hudson thresholds (i.e., 14o and 20o C in the summer (Fewster et
al., 2000; Large et al., 2013; Renecker and Hudson, 1986). I tested the same fixed effects as
above (sex and month), with thin plate regression splines placed on the continuous terms
(temperature and elevation; Wood, 2003). The same random effects structure used for the LMMs
was employed in the AMMs, as well as the best ARMA structure. These models took the form:
Yijk = β0 + ftemp(x1ijk) + felev(x2ijk) + βmox3ijk + βsexx4ijk + u1j + u2jk + ξijk

(2)

where ftemp and felev are unknown smooth functions of temperature and elevation to be estimated
from the data, and all other terms are the same as in (1). The best approximating models were
selected by the lowest AIC score (Burnham and Anderson, 2003). This resulted in a total of 16
LMMs and 12 AMMs (additive models without the continuous covariates are equivalent to linear
models, thus they were not re-fit in the AMM analysis) for each response, giving a total of 28
models under evaluation for the movement metrics in each period. I conducted all statistical
analyses in the R statistical computing environment using packages nlme and mgcv (Pinherio et
al., 2015; R Core Team, 2015; Wood, 2011).
1.2.5. Probabilistic Movement Metrics
I also developed a companion modeling effort using a probabilistic movement estimator
to quantify the moose movement metrics. I performed these additional steps to assess sensitivity
of results to differences in sample size and periods of time where GPS locations were closer in
space, as these factors could influence the bias of movement metric estimates. I predicted these

19

movement paths from continuous-time correlated random walk (CTCRW) models fit separately
for each moose in each summer using the R package crawl (Johnson et al., 2008). The CTCRW
is a state-space model based on Brownian motion, where the estimated GPS location at time =
𝑡 + 1 is a function of the current movement state. This includes the current velocity, and
parameters for velocity autocorrelation and velocity variation, which influence both the predicted
location and the associated Brownian error (Johnson et al., 2008). I then used these models to
predict GPS locations at the same two hour intervals the GPS fixes were attempted and, using
these predicted locations, calculated movement metrics for the same exact periods used before
(i.e., those having 5 or 6 locations obtained on the particular moose during the respective time of
day). With these metrics in hand, I repeated the analysis outlined above exactly as described and
evaluated the similarity in the predictions with tests of correlation. My goal here was thus to
evaluate whether probabilistic movement paths that accounted for autocorrelation and variation
in velocity revealed different patterns from those developed using the straight-line calculations
performed on the raw GPS locations.
1.3.

Results
After screening the GPS data and omitting periods having less than 5 successful GPS

fixes, 17,531 daytime periods and 16,799 crepuscular periods each consisting of 5 or 6 GPS fixes
remained for calculating movement metrics, represented by 152 moose. 95% of the maximum
daytime (Tday) and maximum crepuscular collar temperatures (Tcrep) recorded were between 17
and 41 °C and 13 and 29 °C, respectively. Temperature differences between daytime and
respective crepuscular periods (ΔT = Tday – Tcrep) fell between 0 and 17 °C for 95% of the time.
The collars consistently recorded temperatures near and above the proposed thresholds, during
both daytime and crepuscular hours (Table 1.2). The movement rates calculated in both periods

20

were positively skewed, which was improved by log-transformation prior to all analyses.
Although bounded by 0 and 180°, mean turning angles were approximately normally distributed,
with 97.5% of data falling within 55.3 and 151.2°. Variation in both movements during periods
of similar temperatures was large, although movement rates were characterized by more outliers
than mean turning angles during both daytime and crepuscular periods.
The best autocorrelation structure for crepuscular movement rate was ARMA(2,1) as
determined by both the lowest AIC values, and all other models scored best with the ARMA(1,1)
structure, including all models of daytime movement metrics (Table 1.3). The best
approximating AMMs and LMMs for both daytime and crepuscular movement rate as
determined by lowest AIC included all four terms, with the exception that month was excluded
from the top LMM for daytime movement rate (Tables 1.4 and 1.5). The best approximating
AMMs and LMMs of turning angle during the crepuscular period were of the same form,
including temperature for each estimate tested (Tday, Tcrep, and ΔT), as well as sex and month.
However, temperature was not included in the best approximating AMM for daytime turning
angle, although the global model was within ΔAIC = 1 (Table 1.5). I report the mean turning
angle predictions of the global model in this case (Fig. 1.5), as the main focus here is on the
effect of temperature, and a ΔAIC =1 suggests weak evidence to prefer one over the other
(Burnham and Anderson, 2003). The runner-up LMM for each model of mean crepuscular
turning angle was the global model, and was within ΔAIC = 2 for each estimate of temperature
tested (Table 1.4).
The top AMM for crepuscular movement rate revealed a positive relationship with
increasing Tcrep and increasing Tday, indicating that moose moved greater distances at night when
it is warm during both day and night (Figs. 1.4b and 1.4c). Further, the linear effect of ΔT on

21

crepuscular movements evidenced by the top LMMs was much less than the effects of Tcrep and
Tday on crepuscular movement rate (Table 1.7). This corresponded to a decreasing relationship
between Tday and daytime movement rates, which was also reflected in the opposite signs of the
temperature effects in the LMMs between periods (Tables 1.6 and 1.7). Interestingly, this
decreasing relationship between Tday and daytime movement rate leveled off near ~ 30 °C, while
the increasing relationship between temperatures (both Tday and Tcrep) and crepuscular
movements did not exhibit a tapering off effect (Fig. 1.4a). This suggests rate of movement
during the day may reach a limit even as temperatures during this period continue to rise
relatively high. These findings support the notion that moose move greater distances at night,
even as Tcrep increases.
A distinct temperature threshold was not evident in any of the best approximating
models, although the best approximating AMM for mean turning angle during the daytime
period, as well as that for crepuscular movement rate in relation to ΔT, both indicated a
curvilinear relationship with increasing estimates. There is a notable concave downward
relationship predicted for movement rate as ΔT approaches the extreme observed values, which
suggests the possibility of a threshold in within-day temperature variability on crepuscular
movements (Fig. 1.4d). That mean daytime turning angles began to increase beyond Tday = 30
°C, co-occurring with the leveling off of daytime movement rate, indicates that although moose
may have stopped moving less with increasing temperature, their movements became less
directed (Figs. 1.4a and 1.5a). Much of the variation in each of the movement metrics went
unexplained by the best approximating models, as evidenced by the substantial amount of
residual variation relative to the variation of the random effects (Tables 1.3 and 1.4). The effect
sizes of the temperature terms in the LMMs were small relative to the scale of the response

22

variables, and effects of the other predictors (Tables 1.3 and 1.4). Conducting the same analyses
using movement metrics calculated from CTCRW-derived GPS locations did not reveal any
differences in either LMMs or AMMs; responses were highly correlated at r = 0.949 and r =
0.993 for mean turning angle and movement rate, respectively. Subsequent predictions from
models fit to movement metrics calculated on predicted and raw GPS locations ways were nearly
identical, as evidenced by high correlations between predicted values for the two movement
modeling approaches (r > 0.995 for Figs. 1.6 and 1.7).
1.4

Discussion
I found little evidence of summer temperature thresholds beyond which moose movement

significantly changed. This result comes after I carried out a rigorous assessment of the
movement of 152 wild-living moose, evaluated across vast spatio-temporal extents and over two
different periods (daytime and crepuscular), using two separate modeling frameworks (linear and
additive) for each of two techniques for quantifying movement (raw GPS locations and predicted
locations based on CTCRW models). I did not find any concordance with the temperature
thresholds proposed by Renecker and Hudson (1986), despite repeated exposure to higher
temperatures in both periods. Rather, the effect of maximum daytime and crepuscular
temperature (Tday and Tcrep) on crepuscular movement metrics as estimated in the AMMs
indicated linear relationships between movements and temperature, and of the same direction for
each metric (Figs. 1.4b and 1.4c; 1.5b and 1.5c). These results do not dismiss the possibility that
a temperature threshold in moose movement exists. Rather, they suggest that at a 12-hour scale,
such thresholds do not scale-up to the population-level movement metrics I employed, or are not
applicable to this population of wild-living moose during the summer months. Indeed, this study
is limited by a restriction to the summer months. Additional research should quantify this
relationship during times where altered movement could have more obvious fitness implications
23

(e.g., mating season). The broad-scale nature of this approach also limits the ability to document
behavioral responses to temperature at fine scales (e.g., within-hour movements in a small
habitat patch) and does not necessarily rule out behavioral responses at such scales. Yet the lack
of any threshold response resulting from exposure to previously proposed critical temperatures is
a key finding in light of the evidence that moose alter their behavior (habitat selection and use)
near these same temperatures (Broders et al., 2012; McCan et al., 2013; Melin et al., 2014;
Renecker and Hudson, 1986; Street et al., 2015; van Beest et al., 2012, 2013). The predicted
relationship of crepuscular movement rate with ΔT is also interesting, as it provides evidence
that within-day fluctuations in ambient temperature could be a potential indicator for threshold
responses in movements at the daily scale. Thus, my study supports the notion that ‘thresholds’
may be better referred to as a general range in which ambient temperature could start to trigger
summer heat stress in moose (McGraw et al., 2012; Renecker and Hudson, 1990). Furthermore,
considering variability in temperature in addition to ambient temperature levels may allow for
more precise prediction of the onset of this heat stress.
The stronger positive effects of Tday and Tcrep on crepuscular movement rates that I
documented relative to the effect of ΔT provides evidence that moose move more at night during
warm periods in general (Figs. 1.4b and 1.4c; Table 1.4), which is consistent with previous
findings of Dussault et al., (2004). In that assessment, the authors identified a positive effect of
increasing temperature on nighttime activity, where activity was estimated from motion sensors
affixed to the GPS collars. While they did not find a negative effect of ambient temperature on
daytime moose activity, these results suggest that this may be the case if there exists a positive
relationship between moose activity and movement rates. This is a reasonable assumption, as a
recent study found decreased daytime activity levels of moose during warm periods (Street et al.,

24

2015). These results suggest that moose may exhibit compensatory movement or activity
schedules, moving more at night and less during the day as temperatures increases. That being
said, the size of this effect was much less than the effects of either temperature estimate on
crepuscular movement rate (Tables 1.6 and 1.7). This is likely due to the leveling-off
phenomenon displayed in predictions of daytime movement rates as Tday approached 30 °C (Fig.
1.4a), where the movement rate remained steady at about 48.9 m/hr (SD Âą 57.6 m/hr) thereafter.
Coupled with the concave upward trend of Tday with daytime turning angle near this temperature
range (Fig. 1.5a), these results suggest that moose may start engaging in different behavior at this
point, possibly avoiding harassment by insects or searching for thermal refuge (Renecker and
Hudson, 1990). The differences revealed between the AMMs of mean turning angle for each
period paint a similar picture. While a negative linear relationship was found between
crepuscular turning angles and Tday, Tcrep, and ΔT (Table 1.6), the nonlinear relationship during
the day indicates moose may not exhibit the same behavioral responses to warm temperatures
during the day as they do during the crepuscular and night hours (Figs. 1.4a-1.4d).
The challenge of achieving a mechanistic understanding for the demographic
consequences of high ambient temperatures is further complicated by how intrinsic factors such
as body mass and maturity may affect the movement-temperature relationship. Because I
assessed coarse population-level relationships between movements and temperature, I only
included sex as a fixed effect. This assumes calves do not move any differently than adults, and
that sex captured some of the variation in movement caused by body mass. However, the large
random effect variation suggests that including these other factors in analyses could strengthen
the model fits. Regional variability in other interacting factors, such as disease prevalence,
predation, and changes in both forage quality and selection are also important for informing a

25

mechanistic understanding. If the range of moose is limited by operative temperatures in
summer, assessing the influence of these other biotic factors would be necessary. That mature
conifer stands may be selected during warm periods, but are typically lacking in abundance of
high-quality forage compared to younger forests, suggests a possible trade-off between
thermoregulation and forage availability (Bjørneraas et al., 2011; van Beest et al., 2012). But
careful evaluation of the implications of this general concept for individual and population-level
moose fitness is necessary and currently lacking. In North America, parasites such as the winter
tick (Dermacentor albipictus) might affect fitness, and increased transmission of the meningeal
worm (Parelaphostrongylus tenuis) to moose has been implicated as an overlooked cause of
population declines (Lankester, 2010; Murray et al., 2006). Elevated infection rates have been
attributed to a warming climate, especially during winter, as this has been associated with greater
abundance of the primary mammalian host (white-tailed deer; Odocoileus virginianus) in moose
range (DelGuidice et al., 2002; Whitlaw and Lankester, 1994). Additional examination of the
multifaceted ways in which external biotic factors and ambient temperature regimes can interact
to affect moose populations is necessary if underlying mechanisms do exist.
This study highlights the complex nature of the relationship between abiotic stimuli and
animal movement at broad scales. Studying the effects of abiotic conditions on physiological
stress of captive animals (e.g. Renecker and Hudson, 1986; McCan et al., 2013) provides a
helpful, but perhaps limited, reference point from which to understand the behavior of wildliving animals. Analysis of animal movement data can be illustrative of the physiological and
behavioral states of free-ranging individuals as they interact with their environment (Gurarie et
al., 2016). For instance, increasingly high resolution data can facilitate models that account for
different “modes” or “states” of movement, which depict the interaction of an individual’s

26

internal state and the environment (Morales et al., 2004; Nathan et al., 2008; Patterson et al.,
2009; Schick et al., 2008). Employing flexible frameworks to account for the effects of
environment on these different characterizations of movement might grant an enhanced
understanding of the relationship between physiology and abiotic conditions in wild-living
moose. However, the ability to precisely characterize the abiotic conditions experienced in the
wild will continue to be challenging, particularly where compared to the access that captive
studies afford. For instance, the reliability of collar temperature measurements may diminish as
habitat selection varies with temperature (e.g. selection of open water and marshes; Street et al.,
2015). Such behavioral influence on collar estimates of ambient temperature could possibly
buffer the effects of a threshold response at the temperature extremes. Therefore, there are pluses
and minuses to studies that occur along the captive-to-wild continuum. This point emphasizes the
more philosophical position that choosing the most appropriate modeling techniques from such a
wide, often esoteric range of possible methods can prove to be restrictive for many ecologists
(Gurarie et al., 2016; Schick et al., 2008).
1.5.

Conclusions
Declines in moose population abundance have been attributed to a variety of factors

beyond rising ambient temperatures, but the southern edge of the moose range is still predicted to
shift northward as the climate continues to warm and heat stress becomes more likely (Lenarz et
al., 2010; Murray et al., 2012). However, the thermal conditions under which heat stress in
moose might negatively affect fecundity after accounting for these factors remains unknown, and
the physiological temperature thresholds for moose as inferred from a small sample of captive
animals (Renecker and Hudson, 1986; Table 1.1) should be applied with caution. Determining
temperature thresholds for wild-living moose movement would enhance our ability to predict the

27

effects of abiotic conditions on body condition, and thus survival. My study took an important
step towards that end by identifying a relatively smooth and continuous relationship between
temperature and moose movements. Additional work could build on this study by exploring this
relationship at different scales, in other seasons, or in the context of other potentially important
biotic factors (e.g., habitat heterogeneity).

28

Acknowledgements
I am thankful for a relationship between the RECaP lab and NINA, without which I would not
have been able to conduct this analysis. I acknowledge the effort put forth by personnel at NINA
to conduct a large scale tracking operation that resulted in the high-quality GPS data used herein,
and I am grateful to the Norwegian Environment Agency and the Research Council of Norway
for funding this aspect of the research.

29

APPENDIX

30

APPENDIX
Table 1.1. Review of articles published between 2010-2015 that have directly referenced the
temperature thresholds suggested by Renecker and Hudson (1986) to develop models
or discuss climate related impacts on moose.
Citing article

Journal

Context

Lenarz et al., 2010 [16]

Journal of Wildlife
Management

Lowe et al., 2010 [32]

Canadian Journal of Zoology

Bjørneraas et al., 2011 [42]

Wildlife Biology

Rolandsen et al., 2011 [67]

Ecosphere

Rempel 2011 [68]

Ecological Modeling

Broders et al., 2012 [10]

Alces

Murray et al., 2012 [18]

Canadian Journal of Zoology

van Beest et al., 2012 [17]

Animal Behaviour

Dou et al., 2013 [13]

Ecological Research

McCan et al., 2013 [27]

Canadian Journal of Zoology

van Beest et al., 2013 [28]

PLoS ONE

Melin et al., 2014 [29]

Global Change Biology

“These upper critical
temperatures,” p. 1013
“Thus, when temperatures
exceeded critical thresholds,”
p.1033
“…thresholds that are
regularly exceeded,” p. 51
“…respiratory rate, and
metabolic rate when ambient
temperatures rise above -5°C
in winter and 14°C in
summer,” p. 9
“In summer, thermal stress
begins at about 14 °C and is
fairly severe at 20 °C,” p.
3360
“…relatively low, upper
critical temperature limit,” p.
54
“…upper critical temperature
limits can cause,” p. 431
“Upper critical temperature
thresholds for moose under
captive conditions,” p. 724
“…critical temperature of the
moose,” p. 630
“…indicate heat-stress
thresholds for moose at 14
and 20 °C,” p. 893
“…locations were classified
by temperature in relation to
seasonal thermoregulation
thresholds,” p. 3
“The thresholds for thermal
stress in moose, as suggested
by,” p. 1116

31

Table 1.1 (cont’d).
Citing article

Journal

Olson et al., 2014 [69]

Alces

Street et al., 2015 [30]

Journal of Wildlife
Management

Monteith et al., 2015 [70]

Oecologia

32

Context
“…estimated upper critical
temperatures (Tuc) of moose
as,” p 105
“…particularly at
temperatures exceeding the
upper thresholds,” p. 3
“…and have the lowest upper
critical temperature of any
northern ungulate,” p. 1145

Table 1.2. Counts of temperatures recorded by GPS collars of GPS-tracked moose in central
Norway from July through August of each year during 2006-2010. Counts during
crepuscular period (20:00 – 06:00 next day) above, counts of temperatures recorded
during daytime period (08:00 – 18:00) middle, differences of temperature (ΔT)
between crepuscular and daytime periods below.
Year
2006
2007
2008
2009
2010
Year
2006
2007
2008
2009
2010

< 15
64
458
473
197
60

Tcrep (°C)
[15:17) [17 - 19) [19 - 21) [21 - 23) [23 - 25) [25 - 27) [27 - 29) ≥ 29
123
306
514
546
537
301
99
27
808
1112
1312
875
653
307
164 464
980
1190
979
546
432
271
243 196
332
556
649
365
248
59
35
67
86
141
146
91
57
12
11
47

Tday (°C)
< 21 [21 - 23) [23 - 25) [25 - 27) [27 - 29) [29 - 31) [31 - 33) [33 - 35) ≥ 35
94
94
94
94
94
94
94
94
94
591
591
591
591
591
591
591
591 591
681
681
681
681
681
681
681
681 681
290
290
290
290
290
290
290
290 290
63
63
63
63
63
63
63
63
63

Year
2006
2007
2008
2009
2010

<1
86
259
206
105
30

[1 - 3)
283
768
614
352
76

[3 - 5)
103
311
336
156
36

[5 - 7)
65
185
207
84
20

ΔT (°C)
[7 - 9) [9 - 11) [11 - 13) [13 - 15) ≥ 15
97
506
489
342 213
227
1385
1422
945 564
383
1148
1155
789 550
95
651
544
351 221
31
136
127
91
55

33

Table 1.3. Comparison of ARMA error structure for linear mixed models of crepuscular (left
half) and daytime (right half) movements of moose, sorted by ascending AIC score.
Model
df
AIC
Response = ln(Movement Rate)
ARMA(2,1) 11 33803.96
ARMA(2,2) 12 33805.50
ARMA(1,1) 10 33881.73
ARMA(3,0) 11 33956.68
ARMA(1,2) 11 33992.39
ARMA(0,3) 11 33998.61
ARMA(2,0) 10 33999.84
ARMA(1,0)
9 34040.00
ARMA(0,2) 10 34032.73
ARMA(0,1)
9 34086.27
ARMA(0,0)* 8 34491.44

Model

df

AIC

ARMA(1,1)
ARMA(2,1)
ARMA(2,2)
ARMA(3,0)
ARMA(0,3)
ARMA(2,0)
ARMA(1,2)
ARMA(0,2)
ARMA(1,0)
ARMA(0,1)
ARMA(0,0)*

10
11
12
11
11
10
11
10
9
9
8

34916.56
34917.42
34920.34
34931.68
34937.96
34944.88
34948.26
34953.78
34989.14
34996.7
35065.97

Response = Straightness Index
ARMA(1,1) 10 -529.156
ARMA(2,1) 11 -528.701
ARMA(2,2) 12 -526.248
ARMA(3,0) 11 -507.374
ARMA(0,3) 11 -505.425
ARMA(2,0) 10 -494.310
ARMA(0,2) 10 -492.866
ARMA(1,2) 11 -491.816
ARMA(1,0)
9 -485.026
ARMA(0,1)
9 -484.701
ARMA(0,0)* 8 -480.498

ARMA(1,1)
ARMA(2,0)
ARMA(0,2)
ARMA(2,1)
ARMA(3,0)
ARMA(0,3)
ARMA(1,0)
ARMA(0,1)
ARMA(1,2)
ARMA(2,2)
ARMA(0,0)*

10
10
10
11
11
11
9
9
11
12
8

-2063.72
-2062.65
-2062.47
-2061.73
-2061.5
-2061.46
-2061.22
-2060.92
-2060.65
-2059.92
-2052.76

Response = Mean turning angle
ARMA(1,1) 10 -20749.8
ARMA(2,1) 11 -20748.1
ARMA(2,2) 12 -20746.6
ARMA(3,0) 11 -20690.2
ARMA(0,3) 11 -20687.5
ARMA(2,0) 10 -20686.0
ARMA(0,2) 10 -20684.2
ARMA(1,2) 11 -20683.4
ARMA(1,0)
9 -20673.8
ARMA(0,1)
9 -20672.7
ARMA(0,0)* 8 -20655.9

ARMA(1,1)
ARMA(2,1)
ARMA(2,2)
ARMA(3,0)
ARMA(0,3)
ARMA(1,0)
ARMA(0,1)
ARMA(2,0)
ARMA(0,2)
ARMA(1,2)
ARMA(0,0)*

10
11
12
11
11
9
9
10
10
11
8

-23009.4
-23007.4
-23006.3
-22987.4
-22986.6
-22982.2
-22982.1
-22981.5
-22981.2
-22979.5
-22976.1

* denotes global model with no autocorrelation structure

34

Table 1.4. Comparison of best approximating LMMs. The top three linear mixed models of
movement rate (Rate) and mean turning angle (Turn) response variables during
crepuscular and daytime movement metrics, sorted by increasing AIC score.
Response Model

df AIC
Daytime movements
+ Elevation + Sex
+ Elevation + Sex + Month
+ Sex
+ Elevation + Sex + Month
Elevation
+ Sex + Month
Elevation
+ Sex

11
12
10
10
9
8

36261.37
36270.40
36270.55
167026.1
167027.1
167043.6

Rate
Rate
Rate
Turn
Turn
Turn

Tday
Tday
Tday
Tday

Rate
Rate
Rate
Turn
Turn
Turn

Tcrep
Tcrep
Tcrep
Tcrep
Tcrep
Tcrep

Crepuscular movements
+ Elevation + Sex + Month
+ Sex + Month
+ Elevation + Sex
+ Sex + Month
+ Elevation + Sex + Month
+ Sex

11
10
10
10
11
9

36727.61
36729.52
36767.34
158862.0
158863.3
158868.4

Rate
Rate
Rate
Turn
Turn
Turn

Tday
Tday
Tday
Tday
Tday
Tday

+ Elevation + Sex
+ Sex
+ Elevation + Sex
+ Sex
+ Elevation + Sex
+ Sex

+ Month
+ Month

11
10
10
10
11
9

35914.53
35924.94
35942.18
154935.6
154937.6
154948.4

Rate
Rate
Rate
Turn
Turn
Turn

ΔT
ΔT
ΔT
ΔT
ΔT

+ Elevation + Sex
+ Sex
+ Elevation + Sex
+ Sex
+ Elevation + Sex
Sex

+ Month
+ Month

11
10
10
10
11
9

36204.18
36208.70
36209.52
154990.4
154992.3
154993.7

+ Month
+ Month

+ Month
+ Month
+ Month

35

Table 1.5. Comparison of best approximating AMMs. The top three additive mixed models of
movement rate (Rate) and mean turning angle (Turn) response variables during
crepuscular and daytime movement metrics, sorted by increasing AIC score.
Response Model

df AIC
Daytime movements
+ s(Elevation) + Sex + Month
+ s(Elevation) + Sex
+ Sex
s(Elevation)
+ Sex + Month
+ s(Elevation) + Sex + Month
s(Elevation)
+ Sex

14
13
11
10
12
9

36206.89
36220.19
36241.85
167029.1
167030.1
167045.1

Rate
Rate
Rate
Turn
Turn
Turn

Crepuscular movements
s(Tcrep) + s(Elevation) + Sex + Month
s(Tcrep)
+ Sex + Month
s(Tcrep) + s(Elevation) + Sex
s(Tcrep)
+ Sex + Month
s(Tcrep) + s(Elevation) + Sex + Month
s(Tcrep)
+ Sex

13
11
12
11
13
10

36693.45
36731.52
36731.62
158864
158867.3
158870.4

Rate
Rate
Rate
Turn
Turn
Turn

s(Tday)
s(Tday)
s(Tday)
s(Tday)
s(Tday)
s(Tday)

+ s(Elevation)
+ s(Elevation)

13
12
11
11
13
10

35872.33
35898.85
35926.94
154937.6
154941.6
154950.4

Rate
Rate
Rate
Turn
Turn
Turn

s(ΔT)
s(ΔT)
s(ΔT)
s(ΔT)
s(ΔT)

+ s(Elevation)
+ s(Elevation)

13
12
11
11
13
11

36164.86
36170.36
36206.41
154992.4
154996.3
154997.7

Rate
Rate
Rate
Turn
Turn
Turn

s(Tday)
s(Tday)
s(Tday)
s(Tday)

+ s(Elevation)

+ s(Elevation)
s(Elevation)

+ Sex
+ Sex
+ Sex
+ Sex
+ Sex
+ Sex

+ Month

+ Sex
+ Sex
+ Sex
+ Sex
+ Sex
+ Sex

+ Month

+ Month
+ Month
+ Month

+ Month
+ Month
+ Month
+ Month

36

Table 1.6. Best approximating LMMs, daytime hours. Parameter estimates for best
approximating linear mixed models of daytime movements, as determined by AIC.
Max. daytime collar temperature (Tday) was the temperature covariate used in these
models.
Fixed effect

Estimate SE

Random effect

Std. Dev.

Movement Rate
Intercept
3.714
Temp (Tday)
-9.3E-3
Elevation
2.4E-4
Sex
0.424

0.041 Moose ID
9.6E-4 Year: Moose ID
2.4E-4
0.049

Intercept
Intercept
Residual

0.222
5.2E-4
0.696

Mean Turning Angle
Intercept
92.58
Temp (Tday)
-0.020
Elevation
-8.8E-3
Sex
-6.903

1.236 Moose ID
0.038 Year: Moose ID
1.8E-3
0.886

Intercept
Intercept
Residual

3.163
0.035
28.36

37

Table 1.7. Best approximating LMMs, crepuscular hours. Parameter estimates for best
approximating linear mixed models of crepuscular movements, as determined by AIC.
The temperature covariate used in each model is indicated in parentheses. Tcrep = max.
crepuscular collar temperature, Tday = max. daytime collar temperature, and ΔT = Tday
– Tcrep.
Fixed effect
Movement Rate
Intercept
Temp (Tcrep)
Elevation
Sex
Month

Estimate SE

Random effect

Std. Dev.

3.880
0.033
-1.6E-4
0.577
0.119

0.056 Moose ID
1.7E-3 Year: Moose ID
7.9E-5
0.064
0.018

Intercept
Intercept
Residual

0.261
7.0E-4
0.760

Mean Turning Angle
Intercept
111.61
Temp (Tcrep)
-0.485
Sex
-5.774
Month
1.249

1.254 Moose ID
0.053 Year: Moose ID
0.693
0.451

Intercept
Intercept
Residual

1.828
1.299
24.36

0.050 Moose ID
1.0E-3 Year: Moose ID
7.9E-5
0.063
0.017

Intercept
Intercept
Residual

0.259
6.5E-4
0.759

1.073 Moose ID
0.035 Year: Moose ID
0.680
0.442

Intercept
Intercept
Residual

1.707
1.261
24.39

0.041 Moose ID
1.3E-3 Year: Moose ID
8.0E-5
0.061
0.018

Intercept
Intercept
Residual

0.244
9.1E-4
0.765

0.489 Moose ID
0.045 Year: Moose ID
0.668
0.433

Intercept
Intercept
Residual

1.562
1.491
24.42

Movement Rate
Intercept
Temp (Tday)
Elevation
Sex
Month

4.069
0.020
-3.0E-4
0.570
0.098

Mean Turning Angle
Intercept
108.4
Temp (Tday)
-0.269
Sex
-5.620
Month
1.740
Movement Rate
Intercept
4.550
Temp (ΔT)
0.011
Elevation
-2.1E-4
Sex
0.567
Month
0.049
Mean Turning Angle
Intercept
Temp (ΔT)
Sex
Month

101.3
-0.105
-5.606
2.584

38

Figure 1.1. Summer locations of GPS-tracked moose in Central Norway during 2006-2010.
Norway is shaded light blue, other countries colored olive.

39

Figure 1.2. Mean ln-transformed speeds (Âą SE) of moose movement in Central Norway, during
the months of July and August of each year from 2006-2010 (n = 152 moose). Speeds
were estimated from continuous-time correlated random walk movement models fit to
each summer track using the package crawl in R.

40

Figure 1.3. Boxplots of daytime (a) and crepuscular (b) movement metrics calculated from GPS
data on moose tracked during July-August of 2006-2010. Movement metrics are
binned by maximum collar temperature recorded during the respective period for
which the movement metric was calculated.

41

Figure 1.4. Predictions based on best approximating AMMs (Âą 95% CI shaded grey) for
movement rates of male moose at median elevation during daytime (a) and
crepuscular periods (b-d). Daytime movement rates predicted as a function of daytime
temperature, and crepuscular movement rates predicted as a function of crepuscular
temperature (b), daytime temperature (c), and the difference between daytime and
crepuscular temperature (d).

42

Figure 1.5. Predictions based on best approximating AMMs (Âą 95% CI shaded grey) for mean
turning angles of male moose at median elevation during daytime (a) and crepuscular
periods (b-d). Daytime turning angles predicted as a function of daytime temperature,
and crepuscular turn angles predicted as a function of crepuscular temperature (b),
daytime temperature (c), and the difference between daytime and crepuscular
temperature (d).

43

Figure 1.6. Predictions based on best approximating AMMs of movement rates calculated from
probabilistic movement paths. Predictions for daytime (a) and crepuscular (b-d)
movement rates of male moose. Daytime movement rates predicted as a function of
daytime temperature, and crepuscular movement rates predicted as a function of
crepuscular temperature (b), daytime temperature (c), and the difference between
daytime and crepuscular temperature (d).

44

Figure 1.7. Predictions based on best approximating AMMs of mean turning angles calculated
from probabilistic movement paths. Predictions for daytime (a) and crepuscular (b-d)
movement rates of male moose. Daytime movement rates predicted as a function of
daytime temperature, and crepuscular movement rates predicted as a function of
crepuscular temperature (b), daytime temperature (c), and the difference between
daytime and crepuscular temperature (d).

45

REFERENCES

46

REFERENCES

1) Bakkestuen, V., Erikstad, L., Halvorsen, R. 2008. Step-less models for regional
environmental variation in Norway. Journal of Biogeography, 35:1906 – 1922.
2) Bjørneraas, K., Herfindal, I., Solberg, E. J., SÌther, B.-E., van Moorter, B., Rolandsen, C. M.
Habitat quality influences population distribution, individual space use, and functional
responses in habitat selection by a large herbivore. Oecologia, 168:231 – 243.
3) Bjørneraas, K., Solberg, E. J., Herfindal, I., Van Moorter, B., Rolandsen, C. M., Tremblay,
J.-P., Skarpe, C., SĂŚther, B.-E., Eriksen, R., Astrup, R. 2011. Moose Alces alces habitat use
at multiple temporal scales in a human-altered landscape. Wildlife Biology, 17:44 – 54.
4) Bjørneraas, K., van Moorter, B., Rolandsen, C. M., Herfindal, I. 2010. Screening global
positioning system location data for errors using animal movement characteristics. Journal of
Wildlife Management, 74:1361 – 1366.
5) Broders, H. G., Coombs, A. B., McCarron, J R. 2012. Ectothermic responses of moose (Alces
alces) to thermoregulatory stress on mainland Nova Scotia. Alces, 48:53 – 61.
6) Brinkman, T. J., Deperno, C. S., Jenks, J. A., Haroldson, B. S., Osborn, R. G. 2005.
Movement of female white-tailed deer: effects of climate and intensive row-crop agriculture.
Journal of Wildlife Management, 69:1099 – 1111.
7) Bunnefeld, N., BĂśrger, L., van Moorter, B., Rolandsen, C. M., Dettki, H., Solberg, E. J.,
Ericsson, G. 2011. A model-driven approach to quantify migration patterns: individual,
regional and yearly differences. Journal of Animal Ecology, 80:466 – 476.
8) Burnham, K. P., Anderson, D. R. 2003. Model selection and multimodal inference: a
practical information-theoretic approach. New York, NY: Springer-Verlag.
9) Cederlund, G. 1989. Activity patterns in moose and roe deer in a north boreal forest.
Holarctic Ecology, 12:39 – 45.
10) DelGiudice, G. D., Riggs, M. R., Joly, P., Pan, W. 2002. Winter severity, survival, and
cause-specific mortality of female white-tailed deer in north-central Minnesota. Journal of
Wildlife Management, 66:698 - 717.
11) DelGiudice, G. D., Sampson, B. A., Lenarz, M. S., Schrage, M. W., Edwards, A. J. 2011.
Winter body condition of moose (Alces alces) in a declining population in northeastern
Minnesota. Journal of Wildlife Diseases, 47:30 – 40.
12) Demarchi, M. W., Bunnell, F. L. 1995. Forest cover selection and activity of cow moose in
summer. Acta Theriologica, 40:23 – 36.

47

13) Dou, H., Jiang, G., Stott. P., Piao, R. 2013. Climate change impacts population dynamics and
distribution shift of moose (Alces alces) in Heilongjiang Province of China. Ecological
Research, 28:625 – 632.
14) Dussault, C., Ouellet, J.-P., Courtois, R., Huot, J., Breton, L., Larochelle, J. 2004.
Behavioural responses of moose to thermal conditions in the boreal forest. Ecoscience,
11:321 – 328.
15) Ensing, E. P., Ciuti, S., de Wijs, F. A. L. M., Lentferink, D., ten Hoedt, A., Boyce, M. S.,
Hut, R. A. 2014. GPS based daily activity patterns in European red deer and North American
elk (Cervus elaphus): indication for a weak circadian clock in ungulates. PLoS ONE, 9,
e1D6997.
16) Ericsson, G., Dettki, H., Neumann, W., Arnemo, J. M., Singh, N. J. 2015. Offset between
GPS collar-recorded temperature in moose and ambient weather station data. European
Journal of Wildlife Research, 61:919 – 922.
17) Fewster, R. M., Buckland, S. T., Siriwardena, G. M., Baillie, S. R., Wilson, J. D. 2000.
Analysis of population trends for farmland birds using generalized additive models. Ecology,
81:1970 – 1984.
18) Gurarie, E. 2016. What is the animal doing? Tools for exploring behavioral structure in
animal movements. Journal of Animal Ecology, 85 (1):69 – 84.
19) Gurarie, E., Bracis, C., Delgado, M., Meckley, T. D., Kojola, I., Wagner, C. M. 2016. What
is the animal doing? Tools for exploring behavioral structure in animal movements. Journal
of Animal Ecology, 85:69 – 84.
20) Harley, C. D. G., Hughes, A. R., Hultgren, K. M., Miner, B. J., Sorte, C. J. B., Thornber, C.
S., Rodriguez, L. F., Tomanek, L., Williams, S. L. 2006. The impacts of climate change in
coastal marine systems. Ecology Letters, 9:228 – 241.
21) Humphries, M. M., Umbanhowar, J., McCann, K. S. 2004. Bioenergetic prediction of climate
change impacts on northern mammals. Integrative and Comparative Biology, 44:152 – 162.
22) Johnson, D. S., London, J. M., Lea, M.-A., Durban, J. W. 2008. Continuous-time correlated
random walk model for animal telemetry data. Ecology, 89:1208 – 1215.
23) Karns, P. D. 2007. Population distribution, density and trends In: A. W. Franzmann and C. C.
Schwartz (Eds.), Ecology and Management of the North American Moose (Second edition;
pp. 125 – 139). Washington, D.C.: Smithsonian Institution.
24) Lankester, M. W. 2010. Understanding the impact of meningeal worm, Parelaphostrongylus
tenuis, on moose populations. Alces, 46:53 – 70.

48

25) Large, S. I., Fay, G., Friedland, K. D., Link, J. S. 2013. Defining trends and thresholds in
responses of ecological indicators to fishing and environmental pressures. ICES Journal of
Marine Science, 70:755 – 767.
26) Leblond, M., Dussault, C., Ouellet, J.-P. 2010. What drives fine-scale movements of large
herbivores? A case study using moose. Ecography, 33:1102 – 1112.
27) Leblond, M., St-Laurent, M. H., Côté, S. D. 2016. Caribou, water, and ice – fine-scale
movements of a migratory arctic ungulate in the context of climate change. Movement
Ecology, 4:14. DOI:10.1186/s40462-016-0079-4.
28) Lenarz, M. S., Fieberg, J., Schrage, M. W., Edwards, A. J. 2010. Living on the edge: viability
of moose in Northeastern Minnesota. Journal of Wildlife Management, 74:1013 – 1023.
29) Lenarz, M. S., Nelson, M. E., Schrage, M. W., Edwards, A. J. 2009. Temperature mediated
moose survival in Northeastern Minnesota. Journal of Wildlife Management, 73:503 – 510.
30) Lowe, S. J., Patterson, B. R., Schaefer, J. A. 2010. Lack of behavioral responses of moose
(Alces alces) to high ambient temperatures near the southern periphery of their range.
Canadian Journal of Zoology, 88:1032 – 1041.
31) McCann, N. P., Moen, R. A., Harris, T. R.2013. Warm-season heat stress in moose (Alces
alces). Canadian Journal of Zoology, 91:893 – 898.
32) McKellar, A. E., Reudink, M. W., Marra, P. P., Ratcliffe, L. M., Wilson, S. 2005. Climate
and density influence annual survival and movement in a migratory songbird. Ecology and
Evolution, 5:5892 – 5904.
33) Melin, M., Matala, J., Mehätalo, L., Tiilikainen, R., Tikkanen, O.-P., Maltamo, M., Pusenius,
J., Packalen, P. 2014. Moose (Alces alces) reacts to high summer temperatures by utilizing
thermal shelters in boreal forests – an analysis based on airborne laser scanning of the canopy
structure at moose locations. Global Change Biology, 20:1115 – 1125.
34) Mech, D. L., Fieburg, J. 2014. Re-evaluating the northeastern Minnesota moose decline and
the role of wolves. Journal of Wildlife Management. 78:1143 – 1150.
35) McGraw, A. M., Moen, R. A., Overland, L. 2012. Effective temperature differences among
cover types in northeast Minnesota. Alces, 48:45 – 52.
36) Moen A. 1999. National Atlas of Norway: Vegetation. Hønefoss, Norway: Norwegian
Mapping Authority.
37) Monteith, K., Klaver, R. W., Hersey, K. R., Holland, A. A., Thomas, T. P., Kauffman, M. J.
2015. Effects of climate and plant phenology on recruitment of moose at the southern extent
of their range. Oecologia, 178:1137 – 1148.

49

38) Morales, J. M., Haydon, D. T., Frair, J., Holsinger, K. E., Fryxell, J. M. 2004. Extracting
more out of relocation data: building movement models as mixtures of random walks.
Ecology, 85:2436 – 2445.
39) Murray, D. L., Cox, E. W., Ballard, W. B., Whitlaw, H. A., Lenarz, M. S., Custer, T. W.,
Barnett, T., Fuller, T. K. 2006. Pathogens, nutritional deficiency, and climate influences on a
declining moose population. Wildlife Monographs, 166:1 – 30.
40) Murray, D. L., Hussey, K. F., Finnegan, L. A., Lowe, S. J., Price, G. N., Benson, J.,
Loveless, K. M., Middel, K. R., Mills, K., Potter, D., Silver, A., Fortin, M.-J., Patterson, B.
R., Wilson, P. J. 2012. Assessment of the status and viability of a population of moose (Alces
alces) at its southern range limit in Ontario. Canadian Journal of Zoology, 90:422 – 434.
41) Nathan, R., Getz, W. M., Revilla, E., Holyoak, M., Kadmon, R., Saltz, D., Smouse, P. 2008.
A movement ecology paradigm for unifying organismal movement research. Proceedings of
the National Academy of Sciences, 105:19052 – 19059.
42) Nathan, R., Giuggioli, L. 2013. A milestone for movement ecology research. Movement
Ecology, 1:1. DOI:10.1186/2051-3933-1-1.
43) Olson, B., Windels, S. K., Fulton, M., Moen, R. 2014. Fine-scale temperature patterns in the
southern boreal forest: implications for the cold-adapted moose. Alces, 50:105 – 120.
44) Parmesan, C., Yohe, G. 2003. A globally coherent fingerprint of climate change impacts
across natural systems. Nature, 421:37 – 42.
45) Patterson, T. A., Basson, M., Bravington, M. V., Gunn, J. S. 2009. Classifying movement
behavior in relation to environmental conditions using hidden Markov models. Journal of
Animal Ecology, 78:1113 – 1123.
46) Pettorelli, N., Mysterud, A., Yoccoz, N. G., Langvatn, R., Stenseth, N. C. 2005. Importance
of climatological downscaling and plant phenology for red deer in heterogeneous landscapes.
Proceedings of the Royal Society: Biological Sciences, 272:2357 – 2364.
47) Renecker, L. A., Hudson, R. J. 1986. Seasonal energy expenditures and thermoregulatory
responses of moose. Canadian Journal of Zoology, 64:332 – 327.
48) Renecker, L. A., Hudson, R. J. 1990. Behavioral and thermoregulatory responses of moose to
high ambient temperatures and insect harassment in aspen-dominated forests. Alces, 26:66 –
72.
49) Pinheiro, J. C., Bates, D. M. 2000. Mixed-Effects Models in S and S-PLUS: Statistics and
Computing. New York, NY: Springer-Verlag.
50) Pinheiro, J. C., Bates, D. M., Debroy, S., Sarkar, D., R Core Team. 2015. nlme: Linear and
Nonlinear Mixed Effects Models. R package version 3.1-121.

50

51) R Core Team. 2015. R: A language and environment for statistical computing. R Foundation
for Statistical Computing, Vienna, Austria. URL https://www.R-project.org/.
52) Rempel, R. 2011. Effects of climate change on moose populations: Exploring the response
horizon through biometric and systems models. Ecological Modelling, 222:3355 – 3365.
53) Rolandsen, C. M., Solberg, E. J., Herfindal, I., Van Moorter, B., SĂŚther, B.-E. 2011. Largescale spatiotemporal variation in road mortality of moose: Is it all about population density?
Ecosphere, DOI:10.1890/ES11-00196.1.
54) Rolandsen, C. M., Solberg, E. J., SÌther, B.-E., Van Moorter, B., Herfindal, I., Bjørneraas,
K. 2016. On fitness and partial migration in a large herbivore – migratory moose have higher
reproductive performance than residents. Oikos, DOI:10.1111/oik.02996.
55) Root, T. L., Price, J. T., Hall, K. R., Schneider, S. H., Rosensweig, C., Pounds, J. A. 2003.
Fingerprints of global warming on wild animals and plants. Nature, 421:57 – 60.
56) Sand, H. 1996. Life history patterns in female moose (Alces alces): the relationship between
age, body size, fecundity and environmental conditions. Oecologia, 106:212 – 220.
57) Schick, R. S., Loarie, S. R., Colchero, F., Best, B. D., Boustany, A., Conde, D. A., Halpin. P.
N., Joppa, L. N., McClellan, C. M., Clark, J. S. 2008. Understanding movement data and
movement processes: current and emerging directions. Ecology Letters, 11:1338 – 1350.
58) Schwartz, C. C., Renecker, L. A. 2007. Nutrition and energetics In: A. W. Franzmann and C.
C. Schwartz (Eds.), Ecology and Management of the North American Moose (2nd edition; pp.
441 – 478). Washington, D.C.: Smithsonian Institution.
59) Solberg, E. J., Heim, M., Grøtan, V., SÌther, B.-E., Garel, M. 2007. Annual variation in
maternal age and calving date generate cohort effects in moose (Alces alces) body mass.
Oecologia, 154:259 – 271.
60) Street, G. M., Rodgers, A. R., Fryxell, J. M. 2015. Mid-day temperature variation influences
seasonal habitat selection by moose. Journal of Wildlife Management, 79:505 – 512.
61) Testa, J. W., Adams, G. P. 1998. Body condition and adjustments to reproductive effort in
female moose (Alces alces). Journal of Mammalogy, 79:1345 – 1354.
62) Thomas, C. D., Cameron., A., Green, R. E., Bakkenes, M., Beaumont, L. J., Collingam, Y.
C., Erasmus, B. F. N., de Siqueira, M. F., Grainger, A., Hannah, L., Hughes, L., Huntley, B.,
van Jaarsveld, A. S., Midgely, G. F., Miles, L., Ortega-Huerta, M. A., Peterson, A. T.,
Philips, O. L., Williams, S. E. 2004. Extinction risk from climate change. Nature, 427:145 –
148.

51

63) Thurfjell, H., Spong, G., Ericsson, G. 2014. Effects of weather, season, and daylight on
female wild boar movement. Acta Theriologica, 59:467 – 472.
64) van Beest, F. M., Milner, J. M. 2013. Behavioral responses to thermal conditions affect
seasonal mass change in a heat-sensitive northern ungulate. PLoS ONE,
DOI:10.1371/journal.pone.0065972.
65) van Beest, F. M., Van Moorter, B., Milner, J. M. 2012. Temperature-mediated habitat use
and selection by a heat-sensitive northern ungulate. Animal Behavior, 84:723 – 735.
66) Van Moorter, B., Bunnefeld, N., Panzacchi, M., Rolandsen, C. M., Solberg, E. J., SĂŚther, B.E. 2013. Understanding scales of movement: animals ride waves and ripples of
environmental change. Journal of Animal Ecology, 4:770 – 780.
67) Van Moorter, B., Visscher, D. R., Jerde, C. L., Frair, J. L., Merrill, E. H. 2010. Identifying
movement states from location data using cluster analysis. Journal of Wildlife Management,
74:588 – 594.
68) Whitlaw, H. A., Lankester, M. W. 1994. The co-occurrence of moose, white-tailed deer and
Parelaphostrongylus tenuis in Ontario. Canadian Journal of Zoology, 72:819 – 825.
69) Wood, S. N. 2003. Thin plate regression splines. Journal of the Royal Statistical Society (B),
65:95 – 114.
70) Wood, S. N. 2011. Fast stable restricted maximum likelihood and marginal likelihood
estimation of semiparametric generalized linear models. Journal of the Royal Statistical
Society (B), 73:3 – 36.

52

CHAPTER 2

EVALUATING THE CONSEQUENCES ON ECOLOGICAL INFERENCE OF
AGGREGATING WILDLIFE TELEMETRY DATA WHEN ESTIMATING RESOURCE
SELECTION FUNCTIONS

Abstract
Telemetry data are commonly used to assess animal-habitat relationships via the
quantification resource selection functions (RSFs). As the output of these models are often
designed to inform conservation or management practice, inference is typically desired at the
population-level. Thereby, it has become common practice to aggregate data from all telemetered
animals prior to fitting an RSF. To account for individual variation in selection among
telemetered animals, these models typically include a random effect by animal id. These
approaches are valuable for quantifying broad scale selection, but when the focal population may
be comprised of various intrinsic categories (e.g. sex or age class) or clustered spatially (e.g. two
sub-populations occupying different areas of the landscape), information relating to individual
animal decision-making may be obscured by the act of aggregating data. Here I investigated
individual variation in resource selection among a population of reintroduced elk (Cervus
elaphus) in the Missouri Ozarks. I modeled elk location data, collected from Global Positioning
System (GPS) collars, using a Bayesian discrete choice RSF model. I fit an aggregate-level
model, according to prevailing practice, and then batch-processed models at the level of each
individual elk. I compared the outputs of the top aggregate- and individual-level models via
examination of three metrics; 1) the number and composition of parameters in top models, 2) the
estimates of parameters in global models, and 3) the predicted relative probabilities of use. Via

53

this comparison I discovered substantial variation in all three metrics. Furthermore, I could not
detect any conformity at the individual level by age, sex, or year of reintroduction and telemetrycollaring. My work demonstrates that important ecological variation is lost when resource
selection analyses are aggregated at the population level. I discuss the implications of this
analysis for management and conservation practice and present some guiding principles for
developing RSFs at the individual level.

54

2.1.

Introduction
Research on animal-habitat relationships is a cornerstone of ecological inquiry providing

insights across both theoretical and applied dimensions (Johnson 1980; Morris 2003). Evaluating
animal decision-making in relation to habitat characteristics has implications for optimal
foraging, predator-prey interactions, survivorship, reproduction, life history, and,
correspondingly, population-level processes (MacArthur and Pianka 1966; Charnov 1976;
Rosenzweig 1981, 1991; Morris 2003). Applied assessments of animal-habitat relationships are
predicated on the ability to track animal movement over time. While a variety of data collection
methods inform these efforts, much of habitat selection research rely upon telemetry data
(Craighead and Craighead 1971; Craighead et al., 1971; Kenward 2001; Thomas and Taylor
2006; Montgomery and Roloff 2013). Telemetry technology has expanded dramatically in the
last 50 years with subsequent growth in the methods necessary to model these data. Sophisticated
and complex models have been developed to rigorously quantify animal movement, selection,
habitat suitability, as well as population abundance and performance (Hirzel and Lay 2008;
Gaillard et al., 2010; Thurfjell et al., 2014; Avgar et al., 2016; Boyce et al., 2016).
Examples of the formative techniques for quantifying habitat selection include
compositional analysis (Aebischer et al., 1993), Chi-square goodness-of-fit tests (Neu et al.,
1974; Byers et al., 1984), ranking methods (Arthur et al., 1996), and logistic regression
(Thomasma et al., 1991). Many of these models operate in a similar fashion in that used habitat
units (e.g., at locations detected via telemetry) are compared to habitat units which is deemed to
be either unused or available (i.e., locations not detected via telemetry). The statistical
comparison of used to unused/available was unified under the broader resource selection
function (RSF) framework (Boyce and McDonald 1999; Manly et al., 2002). Herein, the RSF

55

models the use versus unused/available data as a function of environmental variables, or
resources, where the output is a value proportional to the probability that location is used given
the resources present (Manly et al., 2002). In other words, an RSF is always proportional to an
RSPF (resource probability selection function) that quantifies the actual probability of using
landscape units. The ability to accurately estimate an RSPF, which can be desirable in a number
of scenarios, is determined by the extent to which non-use habitat units can be reliably quantified
(Manly et al., 2002). Given that estimation of availability is typically much more feasible than
determination of non-use in wildlife studies, the use/available framework tends to be most
common (Johnson 1980; Thomas and Taylor 1990, 2006). In this way, RSFs can provide insight
into the combinations of resources that are necessary to sustain wildlife populations (Manly et
al., 2002).While there are important sampling elements to consider including the ways in which
available or unused habitat units are measured, the RSF framework has become a widely-used
approach among wildlife ecologists to model animal-habitat relationships (Erickson et al., 2001;
Keating and Cherry, 2004; Arts et al., 2008; Fieburg et al., 2010; Montgomery and Roloff,
2013).
Given that management and conservation efforts are typically developed at the
population-level, researchers typically deploy telemetry collars/tags on a number of animal
subjects and then aggregate the resultant telemetry data at the population-level to fit one RSF
(Thomas and Taylor 2006). Importantly, aggregating telemetry data across animal subjects
obscures any individual variation that may exist in the ways in which these animals respond to
the environment. This can bias inference, particularly when relocation data are unbalanced across
individual animals, as is often the case (Gillies et al., 2006; Aarts et al., 2008). To account for
individual variation, researchers may fit models where animal id is included as a random effect.

56

This random effects component (fit as either a random intercept or random slope) relaxes the
assumption of independence among data points, potentially mediating issues with unbalanced
sampling efforts by allowing parameters in an RSF to vary according to an aggregate-level
probability distribution (also called population-level or top-level; Mysterud and Ims, 1998;
Gillies et al., 2006; Thomas, Johnson and Griffith 2006; Aarts et al., 2008; Hebblewhite and
Merrill 2008; Duchesne et al., 2010; Hooten et al., 2016). However, adjustments to the slope or
intercept of a model may be insufficient to account for the true variation inherent to the input
data. At risk here is the potential to misidentify resource selection which can confound inference
and problematize prevailing management or conservation philosophy. As Marzluff et al. (2004)
note, individuality in wildlife resource selection can be extreme, warranting an examination of
the consequences of that individuality.
Here I investigated the consequences of aggregating animals, where variation is assumed
to be constrained statistically, in RSF analyses. Modeling elk (Cervus elaphus) resource selection
in a reintroduced population in southern Missouri, I compared aggregate-level to individual
animal RSFs. Using a hierarchical discrete choice model to evaluate these RSFs, I treated
individuals as random effects at the aggregate level according to conventional practice. I then fit
the same model at the level of each individual elk to compare variation in selection tendencies
using multivariate techniques. I based this comparison on three areas of inference common to
resource selection analyses: 1) the number and composition of parameters in top models, 2) the
estimated RSF parameters of global models, and 3) the predicted relative probabilities of use. I
discuss the ramifications of variation across these three metrics for management and
conversation practice. Finally, I provide guidance on the framing of future analyses of RSFs and
habitat selection research more broadly.

57

2.2.

Methods
2.2.1. Resource selection data
Elk (n = 88) were captured in Kentucky for relocation to the Missouri Ozarks in June of

2011 (14 males, 19 females), 2012 (8 males, 26 females), and 2013 (3 males, 36 F). After a
period of acclimation in a fenced enclosure, these elk were introduced onto Peck Ranch
Conservation Area (Fremont, MO), a 9,327 ha plot of land managed by the Missouri Department
of Conservation. Both of these areas are contained within the broader elk restoration zone, an
896 km2 study area in southeast Missouri delimited by the Missouri Department of Conservation
because of the region’s high density of abundant public lands and high predicted habitat
suitability for elk (MDC 2010). Prior to release, all elk ≥ 1 year old were fitted with a GPS collar
(RASSL custom 3D cell collar, North Star Science and Technology, LLC, King George, VA, or
G2110E Iridium/GPS series model, Advanced Telemetry Systems, Insanti, MN). These collars
were set to a 2.5 hour fix attempt schedule, storing locations internally at this frequency and
uploading information every 5 hours (aside from two collars, set at 2 and 4 hours, respectively).
Capture and handling protocols were approved by the University of Missouri Animal Care and
Use Committee (Protocol 6909). Environmental variables were measured at the used locations
(as determined by the GPS telemetry system) and available locations, which were defined using
the radius of available habitat method (Durner et al., 2009). This involved delimiting a circle
around each used location having radius equal to 𝑐(𝑎 + 2𝑏), where a, b, and c represent the
mean hourly movement rate, the standard deviation of the movement rate, and the number of
hours between locations, respectively. Locations were then randomly sampled from within these
circles to determine the corresponding available locations for the subsequent used locations,
hereafter referred to as a “choice set”.

58

2.2.2. Environmental variables
I used a geographic database of environmental variables developed by Smith (2015) to
describe the study area. These included 11 different variables each depicted as rasters at a
resolution of 30 m. This database included percent tree canopy cover (2011 US Forest Service
National Land Cover Database (www. mrlc.gov/nlcd11_data.php, accessed 8 Jan 2015), number
of years since prescribed burn, aspect (degrees), slope (percent), distance to wooded edge (m),
the interspersion and juxtaposition index (IJI; Griffith et al., 2000), road density (km paved or
gravel road/km2 within 95 km2 circle), distance to paved road (m), distance to closed two-track
roads (m), and distance to public gravel roads (m; for more information see Smith (2015)).
Finally, I also considered habitat type, the only categorical variable considered, which had eight
categories: warm-season grassland, cool-season grassland, shrubland, woodland, savannah,
forest, glade, and forage opening.
2.2.3. RSF modeling
I fit discrete choice RSFs because they can accommodate availability data that is defined
separately for each location, which can provide a more realistic approximation of animal choice
(Cooper and Millspaugh, 1999, 2001; McCracken et al., 1998). I selected a Bayesian framework
given the difficulties of likelihood-based approaches and the flexibility of Bayesian methods
when fitting random-effects logistic regression models (Browne and Draper 2005; Gelman et al.,
2014a). For every used location I developed 5 available locations to define each discrete choice
set. I developed the RSFs at the aggregate- and individual-levels. The individual-level model is
defined as follows, where the probability of an individual elk choosing alternative 𝑙 from a set of
𝐶 feasible habitat units at unit 𝑖 is given by

59

𝑃𝑖𝑙 = ∑𝐶

𝑒 𝜷𝑿𝑖𝑙

𝑐=1 𝑒

(1)

𝜷𝑿𝑖𝑐

where
𝜷𝑿𝑖𝑙 = 𝛽1 𝑋𝑖𝑙1 + 𝛽2 𝑋𝑖𝑙2 + ⋯ + 𝛽𝑘 𝑋𝑖𝑙𝑘

(2)

is the utility of unit 𝑙 to the individual being considered, consisting of k slope parameters
measured on each used and available unit. The concept of utility derives from economic theory
and is analogous to satisfaction; eqn. (1) is valid under the assumption that the individual in
question always chooses the resource unit having the greatest utility (Cooper and Millspaugh
1999). I slightly extended eqns. (1) and (2) to develop the aggregate-level model, where I make
explicit the probability of individual 𝑗 choosing alternative 𝑙 from a set of 𝐶 feasible alternatives
to unit 𝑖, defined as

𝑃𝑖𝑗𝑙 =

𝜷 𝑿
𝑒 𝒋 𝑖𝑗𝑙

(3)

𝜷𝒋 𝑿𝑖𝑗𝑐
∑𝐶
𝑐=1 𝑒

with utility function now defined as
𝜷𝑗 𝑿𝑖𝑗𝑙 = 𝛽𝑗1 𝑋𝑖𝑗𝑙1 + 𝛽𝑗2 𝑋𝑖𝑙𝑗2 + ⋯ + 𝛽𝑗𝑘 𝑋𝑖𝑗𝑙𝑘

(4)

and 𝛽𝑗𝑘 ~𝑁𝑜𝑟𝑚𝑎𝑙(𝜇𝑘 , 𝜎𝑘2 )
where the individual parameters for selection of environmental variable 𝑘 (the 𝛽𝑗𝑘 for all 𝑗 =
1, 2, … , 88 individuals), are assumed to be normally-distributed random effects following some
population distribution. The parameters of the population distributions for each environmental
variable 𝑘 (𝜇𝑘 and 𝜎𝑘2 ) are referred to as “hyperparameters”(Hobbs and Hooten 2015). Thus, the
aggregate-level model was a hierarchical random slopes model, where inference can be made on

60

the central tendency of selection for each environmental variable 𝑘, or the “mean
hyperparameter” 𝜇𝑘 , as well as variation among individuals, or the “standard deviation
hyperparameter” 𝜎𝑘 . I conducted aggregate-level model estimation using the following
uninformative priors for all hyperparameters:
𝜇𝑘 ~𝑁𝑜𝑟𝑚𝑎𝑙(0, 10)

(5)

𝜎𝑘 ~𝑈𝑛𝑖𝑓𝑜𝑟𝑚(0, 10)
Prior to model fitting, I examined evident collinearity among the environmental variables
on a per-individual basis, and excluded redundant environmental variables until all pairwise
correlations were |r| ≤ 0.6. I then constructed a global model at the aggregate and individual
levels consisting of all remaining environmental variables such that this condition was satisfied.
I fit these models in the Bayesian package Stan (Stan Development Team 2016a) using R
and RStan as an interface (R Core Team 2016; Stan Development Team 2016b). For each
individual model, I used four chains of 1000 draws each, with a burn-in period of 200. This was
likely more simulation than necessary, given the simple structure of the models (i.e. nonhierarchical) and in most cases 1000 iterations is more than enough to reach convergence to the
posterior distribution using Stan (Vehtari et al., 2016). However, I used this approach to ensure
that any individual differences that may arise in the comparisons (next section) would not be
attributable to Monte Carlo error. Additionally, to fit all models, I used the high-performance
computing cluster developed and maintained by the Institute for Cyber-Enabled Research at
Michigan State University. With the large number of processors and excess RAM, I was able to
fit the models in parallel and remotely. I assessed convergence of all models by ensuring that for
all parameters the potential scale reduction factor, 𝑅̂ , was below 1.1 and the effective sample

61

size, 𝑛̂𝑒𝑓𝑓 , was greater than 100 (Gelman et al., 2013). I assessed goodness of model fits using
posterior predictive checks, whereby I computed the probability that a test statistic, 𝑇, calculated
on new data simulated from the model, 𝒚∗ , is more extreme than 𝑇 calculated on observed data, y
(Gelman et al., 2013; Hobbs and Hooten 2015). I used the chi-square test statistic to conduct
these checks for the aggregate-level global model and all individual-level global models as
follows
Pr(𝑇(𝒚∗ , 𝜽) ≥ 𝑇(𝒚, 𝜽)|𝒚)
𝑇(𝒚, 𝜽) = ∑𝑖

(6)

(𝑦𝑖 −𝑝𝑖 )2
𝑝𝑖

where 𝜽 represents the parameters for the fitted model, and 𝑝𝑖 is the probability associated with
the 𝑖 th choice. The first expression in (6) returns a Bayesian “𝑝-value”, 𝑃𝐵 , which I used to
diagnose lack of model fit (Gelman et al., 2013; Hobbs and Hooten 2015).
2.2.4. Comparison between aggregate- and individual-level RSFs
I assessed patterns of individual variation in resource selection among the individual and
aggregate levels by comparing: 1) the parameters included in top models, 2) parameter estimates
of all global models, and 3) predictions of the relative probability of use expressed across the
study area. To make comparisons based on the model selection approach, I fit all subsets of the
global model for each individual following examination for collinearity, carrying out estimation
as described above to perform what was essentially a “dredge” (Wiens et al., 2008; Cade 2015). I
ranked models using WAIC, a fully Bayesian estimate of out-of-sample predictive ability,
determining the top model for each individual to be that with the lowest WAIC score (Watanabe
2010; Vehtari et al., 2015; Gelman et al., 2014b). I evaluated the inclusion of parameters in the

62

top and runner-up models by calculating “inclusion rates” in the top 5% of models for each
individual. To make this calculation, I ordered all models by increasing WAIC score and, using
only the models contained in the top 5% of this ordering, divided the number of models in which
a parameter was included by the total number of models in this list. I used WAIC to rank models
because it is both fully Bayesian and computationally efficient (Vehtari et al., 2015). I also
performed this selection procedure for the aggregate-level model, where sub-models retained the
same hierarchical format used in the global model, but varied with respect to the RSF parameters
and associated hyperparameters.
I compared individual variation in parameter estimates to the estimated aggregate-level
distribution on a per-parameter basis using a graphical approach. As an aggregate-level estimate
of the individual variation in parameter estimates, I added the lower and upper limits of the mean
hyperparameter CIs to their corresponding limits for the standard deviation hyperparameters for
(𝑙)

(𝑙)

(𝑢)

(𝑢)

each environmental variable, giving the following interval [𝜇̂ 𝑘 − 𝜎̂𝑘 , 𝜇̂ 𝑘 + 𝜎̂𝑘 ] for all k
parameters, where l and u represent the lower and upper limits for the 95% CIs, respectively.
This calculates conservative intervals within which one standard deviation of the individual
random effects (~ 66% of theoretical individuals’ parameters) is expected to be contained under
the aggregate-level model assumptions. I refer to this metric as the 95% “random-effect CI”, as it
characterizes the random-effect distributions with uncertainty around the hyperparameter
estimates incorporated. Additionally, I investigated patterns of estimated selection tendencies
(i.e., the composites of individual parameter estimates) among individuals using principal
components analysis on the RSF parameter point-estimates. Ordination in parameter space has
been used previously to investigate variation in selection and movement patterns among
individuals (Hanks et al., 2011; Pape and LĂśffler 2015). I only used the individual-level

63

estimates for which uncertainty was comparable among individuals, and I compare among
intrinsic factors and year of GPS-collaring and release, hereafter referred to as “cohort”.
For the final comparison, I produced predictive maps of the relative probabilities of use
for each top model across the entire study area. I evaluated the logistic of the estimated RSFs,
exp(𝑅𝑆𝐹) /(1 + exp(𝑅𝑆𝐹), at each 30m resolution cell within the study area, returning the
predicted relative probability of use of that unit under the estimated model(s).
2.3.

Results
The 88 elk in this study were GPS-tracked between June 1, 2011 and September 15, 2014

(35 elk released in 2011 (22 female, 13 male), 24 elk released in 2012 (17 female, 7 male), and
29 elk released in 2013 (26 female, 3 male)). At the aggregate level, I fit the global discrete
choice RSF model and sub-models using 141,197 choice sets consisting of 1 used location and 5
available locations. I detected high collinearity (0.61 ≤ |𝑟| ≤ 0.84.) between the distance to gravel
road and road density variables for all individuals, and I excluded the former variable from
consideration given that two other variables (distance to paved road and distance to two-track
road) quantified proximity to roads. Thus, an identical set of sub-models was fit to all
individuals. Individual-level models were estimated based on 95 to 4865 choice sets depending
on the elk, and more than 50% of the models were fit with between 783 and 2071 choice sets. I
achieved convergence of all models fit (𝑅̂ < 1.1 for all parameters at the individual level and
hyperparameters at the aggregate level). The Bayesian p-value test for the global aggregate-level
model indicated no lack of fit, i.e. could have reasonably generated the observed data (𝑃𝐵 =
0.35, Hobbs and Hooten 2015). Bayesian p-values for all individual-level global models fell
between 0.16 and 0.58, with an average value of 0.33.

64

2.3.1. Number and composition of parameters
The global model at the aggregate level had the highest predictive ability of all submodels considered (see Table 2.1), while the top individual-level models had an average of 14.4
parameters with a range of 7 to 17 (Fig. 2.1a). The order of importance of parameters for
predictive ability was also mixed, although habitat and slope were in all but two and five of top
models, respectively. Conversely, aspect was only included as a predictor in the top model for
44% of the elk (39 of the 88 top models). The exclusion of any one parameter from top models
did not appear to be linked with patterns of exclusion for other parameters. The remaining
parameters were included in top models for between 62% and 86% of the elk (between 54 and 76
of the 88 top models, respectively). Habitat and slope also had the highest inclusion rates, having
a value of 1.0 for the majority of elk (i.e., included in the top model and all of the 5% runner-up
models; Fig. 2.1b) I observed low inclusion rates (≤ 0.2) of the distance to two-track roads
parameter for four elk and the distance to paved roads parameter for three elk, as well as IJI for
one elk (Fig. 2.1b).
2.3.2. Parameter estimates
At the aggregate level, thirteen of the seventeen 𝜇𝑘 estimates under the global model had
95% credible intervals (CIs) which did not overlap zero. These included estimates for all 𝜇𝑘 for
continuous environmental variables as well as all categorical habitat types except forest,
shrubland, glade, and warm-season grassland. However, the estimated 95% random-effect CIs
contained zero in only two cases: slope and forage opening. Aggregate-level point estimates of
individual variation (the 𝜎𝑘 estimates) varied from a minimum of 0.04 for aspect, to 1.71 for
distance to two-track. The latter was large relative to estimates of 𝜎𝑘 for the remaining
environmental variables, as all other estimates fell below 0.67 (Fig. 2). There was a high degree
65

of variation in the individual-level estimates of selection, and point estimates of the habitat type
parameters were exceptionally variable, extending far beyond random-effect CIs (Fig. 2).
However, uncertainty around these point estimates was high for most individuals, as evidenced
by the large CI’s around these estimates (Fig. 2.2). Conversely, variability among individuals in
selection for distance to two-track roads appeared much smaller than was estimated at the
aggregate-level, as virtually all point estimates fell within the random-effect CIs (gray boxes,
Fig. 2.2). Given the relatively large uncertainty around individual-level habitat type parameter
estimates, I based the principal components analysis (PCA) on parameters for only the
continuous environmental variables. Ordination by PCA revealed that some individuals were
much more different from one another in terms of combined selection estimates of the nine
resource covariates, yet none of the factors examined accounted for the observed extreme
selection tendencies (Fig. 2.3a). A small degree of clustering was evident based on cohort (Fig.
2.3b).
2.3.3. Predicted relative probabilities
Predictive mapping of the relative probability of use throughout the elk restoration zone
under the global models was vastly different among individuals, and when compared to the
aggregate-level predictions (Fig. 2.4a). Of the randomly-selected individuals I used to inform
these plots, less than half displayed patterns of selection similar to the estimated population-level
predictions. The remaining individuals varied with respect to both general areas of relative use
probability and uniformity of this use. For example, elk, at an aggregate level, were predicted to
select habitat units found in the southeastern portion of the elk restoration zone where the
predicted probability of use was very low for many of the individual elk (Fig. 2.4b). Yet even

66

among these individuals, areas of high predicted relative use varied along a gradient of uniform
patterns to being much more clustered, varying in scales (Fig. 2.4b).
2.4.

Discussion
Under the current paradigm for evaluating wildlife resource selection, telemetry data of

collared/tagged individuals are typically aggregated prior to model-fitting (Thomas and Taylor,
2006). Aggregating data among individuals is done in the interest of making inferences about
animal populations that can provide, among other things, information that is relevant to
management and conservation. My study revealed substantial differences between models fit at
the aggregate and individual animal levels across three metrics including the number and
composition of parameters in top models, the estimated parameters among global models, and
predicted relative probabilities of use. This analysis demonstrates that variation in resource
selection among individual animals may not be well accounted for by simply including a random
effect in aggregate models. While results speak to variation among one population of elk in this
ecosystem in Missouri, I find this variation to particularly interesting given that my study animal
is highly gregarious and often expected to respond similarly to the environment. These results
support existing research which demonstrates the importance of considering individual animal
differences in habitat selection, even when the species of interest exhibits sociality (Gillingham
and Parker, 2008; Putman and Flueck, 2011; Pape and LĂśffler, 2015). Further, I could not detect
any real conformity in elk-habitat relationships by age or sex class (Figs. 2.3a and 2.3b).
At the aggregate-level, the global model received the most support based on ranking by
WAIC, yet the inclusion of parameters among individual top models was highly variable. The
assumption that observations (i.e. choices between available habitat units) follow some
population process with parametric form may have influenced which variables are considered
67

important to selection by all individuals. The fact that the global aggregate-level model achieved
the best ranking by WAIC could also be affected by the large sample size, given the tendency for
information criteria-based ranking methods to favor more complex models with large sample
sizes (Piironen and Vehtari, 2016). The individual variation revealed by our model-ranking
analysis suggests that a particular resource may influence decision-making for certain elk more
than others. Differential selection may arise in cases where there is a functional response in
resource selection, or differences in resource selection/avoidance among individuals (Mysterud
and Ims 1998; Haydon et al., 2008). It is difficult to distinguish the effect of overall available
habitat from the effect of individual differences in selection, but analysis at the individual level
affords the opportunity to explore why individual animals may apparently select habitat units
differentially based on available resources. For example, a comparison of RSF models for
individual moose made by Gillingham and Parker (2008) revealed that inclusion of a particular
selection parameter in the top model varied with the availability in an individual’s home range.
Investigating why individuals may differ in this way can provide information that would be
obscured if telemetry data were aggregated and only one model was fit.
Variation in parameter estimates at the individual level was not commensurate with the
variation identified by the σk estimates (representing random effect of animal ID) of the
aggregate-level model. For example, with respect to distance to two-track roads, variation among
individual-level parameter estimates was much lower than what would be suggested by the
aggregate-level model (Fig. 2.2). This observation reflects a potential issue with standardizing
variables when evaluating RSFs at the aggregate level. It is common practice to standardize
predictor variables to facilitate comparison of effects when evaluating RSFs (e.g. Thomas et al.,
2006). In the present study, the aggregate-level hyperparameter estimates were all based on the

68

standardized environmental data at used and available habitat units for 88 elk, while our
individual-level parameter estimates were based on data standardized for each individual
separately. This should provide a more accurate depiction of the comparative effects between
environmental variables on an individual basis; in the case of distance to two-track roads, it
suggests that elk respond much more consistently to this environmental variable than would be
depicted by the aggregate-level model, where it has a small, negative effect on selection of a
resource unit, as opposed to an unknown effect that could vary in both magnitude and direction
(Fig. 2.2). Differential use of habitat units based on habitat type may, in part, account for the
relatively large uncertainty of corresponding individual-level parameter estimates (Fig. 2.2). If
avoidance of specific habitat types by individuals occurs, low use or non-used of these habitat
types can bias parameter estimates, which is a known issue with categorical covariates in such
models (problem of separation; Menard 2002). Modelling at the individual level allows further
investigation into relationships between sampling and RSF parameter estimates (Gillingham and
Parker, 2008).
Spatial variation in the predicted relative probability of use indicated broad-scale
differences in selection strategies between aggregate- and individual-level models (Figs. 2.4a and
2.4b). Because these predictive maps do not reflect parameter uncertainty, it is difficult to
quantify the extent that this variation may be attributed to error as opposed to real differences in
habitat selection strategies among individual elk. It may be that individual animals are selecting
for specific environmental features, and the configuration of those features should be expected to
be heterogeneously distributed (Aarts et al., 2008; Paton and Matthiopolous, 2015). Additionally,
given that predictive mapping with RSFs is a very useful tool for management of ungulate
species, assessing individual-based predictions could offer unique advantages (Johnson et al.,

69

2004). For example, predicting relative probability of landscape use (or the absolute probabilities
of use when using an RSPF) for individuals could help inform predictions of space use as
populations increase following a reintroduction. Predictive mapping based on individual RSF
models may present greater potential opportunities for conservation and management of solitary
animals, or animals of conservation concern whereby efforts may be focused on managing
habitat to meet the needs specific individuals.
With recent developments in technologies that can produce highly resolute movement
and remotely-sensed data, we have an unprecedented ability to model habitat use by individual
animals that can even reveal insights into specific behavioral and ecological mechanisms (Kays
et al., 2015). The idea that consistent individual differences among animals (i.e., personality)
have important evolutionary and ecological consequences has been amassing a large amount of
support for many processes, including habitat selection (RĂŠale et al., 2010; Leclerc et al., 2016;
Merrick et al., 2017; Spiegel et al., 2017). Indeed, the potential to derive fitness implications for
populations based on differences in movement among individuals has been demonstrated for elk
(Haydon et al., 2008; Morales et al. 2010). In this study, aggregate models missed important
inferences for conservation and management, including: the importance of certain resource
covariates in predicting selection, such as aspect; the inflated effects of certain resource
covariates on selection, such as distance to paved roads; that more extreme individuality in
combined selection habits could not be attributed to intrinsic factors; and that selection of
resource units across the greater landscape may not be as consistent. These results are consistent
with an emerging body of research demonstrating potentially substantial variation in habitat
selection strategies within populations of highly gregarious species (e.g., elk and reindeer
Rangifer tarandus; Sawyer et al., 2007; Pape and LĂśffler, 2015)

70

My study demonstrates that modelling resource selection at the individual level can
reveal important variation in selection strategies within populations of animals that could have
been missed by aggregate-level models. That being said, I did identify some consistencies
between the aggregate and individual levels. For example, the consistent inclusion of slope and
habitat type in top individual models compliments the result that slope and forage openings had
the greatest effects on selection at the aggregate level. However, apparent inter-individual
differences in selection as determined by RSF modelling can only be investigated by fitting
individual-level models, because the selection tendencies of individuals is a prerequisite that is
not estimable when fitting typical aggregate (i.e. hierarchical) models. Other benefits to
modelling resource selection at the individual level were beyond the scope of our study but
provide opportunities for future study. Quantifying variation in the movement process at the
individual-animal level could be particularly relevant. There have been a number of individualbased movement models developed in recent years, and methods that simultaneously estimate
resource selection and movement have been introduced (Thurfjell et al., 2014; Avgar et al.,
2016). Employing models of animal movement can be useful for defining more biologically
feasible unused or available habitat units.
In conclusion, via this assessment we do not suggest that researchers should abandon
aggregate-level RSF modelling with random-effects to quantify individual variation. In this
study, relatively strong effects of slope and forage opening on selection of a habitat unit were
found at both the aggregate and individual levels. However, I caution against discounting
individual variation in habitat selection. Recent developments in this area highlight the common
consideration of individual differences in studies of behavioral or population ecology, a
viewpoint which should be adopted in applied wildlife conservation research (Merrick and

71

Koprowski, 2017). Thus, researchers and managers should fit models at the individual level so
that they can appreciate the consequences of that variation on management practice. When using
RSFs to quantify resource selection, I reiterate the importance of multi-scale efforts that have
been advocated by others, and assert that both individual- and aggregate-level models can
provide valuable inference for management (Gillingham and Parker, 2008; Pape and LĂśffler,
2015). The methods that we used here provide a framework that other ecologists and managers
can use to assess the role of individual variation in their study system. These approaches can be
useful for evaluating the utility of population inference from RSFs which will vary among
systems and research objectives (Hanks et al., 2011). Computational improvements to
hierarchical methods for estimating movement and resource selection are also being developed
(Hooten et al., 2016). For researchers and managers limited by computational resources or
statistical expertise, fitting individual models provides a better alternative, or at least a
supplement to, aggregating across individuals for population-level inference in resource selection
analyses. In studies of resource selection, approaches that do not consider individual-level
selection are liable to giving superficial weight to certain resources and/or habitats, while
underestimating the importance of others and potentially missing interesting ecological
relationships. I therefore recommend that whenever feasible, practitioners model habitat
selection at the individual level if aggregate approaches will be used to inform management
decisions.

72

Acknowledgements
I express my deepest gratitude to T. Smith, B. Keller, and J. J. Millspaugh for preparing this
dataset and granting me access to conduct the analysis contained herein. I thank R. A.
Montgomery and J. J. Millspaugh for invaluable guidance during this study. I humbly thank the
Missouri Department of Conservation and the United States Fish and Wildlife Service for
permission to use this data as well, as it was collected as part of a project funded by the MDC
and a USFWS Wildlife Restoration grant.

73

APPENDIX

74

APPENDIX

Table 2.1. WAIC values for top 10 aggregate-level resource selection function models. +
indicates an environmental variable included in the model, spaces left blank if not
included. Environmental variable abbreviations as follows: D2t = distance to twotrack road; Dpav = Distance to paved road (m); Slope = slope (degrees); Aspect =
aspect (degrees); Rd dens = road density (km/km2); D edge = Distance to nearest
wooded edge; % Can = percent forested canopy cover; Rx burn = years since
prescribed fire; Habitat = indicator variables for 8 cover types.
D2t
+
+
+
+
+

Dpav
+
+
+

+

+
+
+

+

+

Slope
+
+
+
+
+
+
+
+

Aspect Rd dens
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+

D edge
+
+
+
+
+
+
+
+

% Can
+
+
+
+
+
+
+
+

75

Rx burn
+
+
+
+
+
+
+
+
+

Habitat
+
+
+
+
+
+
+
+
+

WAIC
321786.1
323171.8
323894.1
324310.7
324410.7
324608.6
324759.4
324865.3
326220.7

Figure 2.1. Inclusion of parameters at the individual level in top models (a panel) and top 5% of
models (b panel) as ranked by WAIC for elk reintroduced into the Missouri Ozarks
(2011-2014). Green boxes (panel a) indicate the selection parameters in line with the
row was included in the top model as ranked by WAIC. Shade of purple (panel b) is
proportional to the number of times the corresponding parameter was included in the
top 5% of sub-models, where the darkest shade corresponds to all models and lightest
shade indicates that parameter was not included in any of the top models.
a

b

76

Figure 2.2. Aggregate-level random effects distributions (boxes) and individual-level RSF
parameter estimates (blue dots) with 95% Bayesian credible intervals (CIs, in
translucent blue) for reintroduced elk in the Missouri Ozarks (2011-2014). The
middle horizontal bars within the boxes are the point estimates of the mean
hyperparameters of the random effects distributions. The light gray horizontal bars
within the boxes are point estimates for the standard deviation hyperparmeters, and
the ends of the boxes represent 95% Bayesian credible intervals on mean + standard
deviation point estimates.

77

Figure 2.3. Principal components analysis (PCA) ordination of the 88 individual elk reintroduced
into the Missouri Ozarks (2011-2014) based on individual-level parameter estimates
for selection of 9 resource covariates. Grouped by intrinsic factors (panel a) and by
GPS-collaring and release year (cohort; panel b).

78

Figure 2.4. Predicted relative probabilities of use of the elk restoration zone based on and
aggregate-level RSF (panel a) and individual-level RSFs for 20 randomly selected elk
(panel b) reintroduced to the Missouri Ozarks (2011-2014).
a

b

79

REFERENCES

80

REFERENCES

1) Aarts, G., MacKenzie, M., McConnell, B., Fedak, M., Matthiopolous, J. 2008. Estimating
space-use and habitat preference from wildlife telemetry data. Ecography, 31(1):140 – 160.
2) Aebischer, N. J., Robertson, P. A., Kenward, R. E. 1993. Compositional analysis of habitat
use from animal radio-tracking data. Ecology, 74(5):1313 – 1325.
3) Arthur, S. M., Manly, B. F. J., McDonald, L. L., Garner, G. W. 1996. Assessing habitat
selection when availability changes. Ecology, 77(1):215 – 227.
4) Avgar, T., Potts, J. R., Lewis, M. A., Boyce, M. S. 2016. Integrated step selection analysis:
bridging the gap between resource selection and animal movement. Methods in Ecology and
Evolution, 7:619 – 630.
5) Boyce, M. S., Johnson, C. J., Merrill, E. H., Nielsen, S. H., Solberg, E. J., van Moorter, B.
2016. Can habitat selection predict abundance? Journal of Animal Ecology, 85(1):11 – 20.
6) Boyce, M. S., McDonald, L. L. 1999. Relating populations to habitats using resource
selection functions. Trend in Ecology and Evolution, 14(7):268 – 272.
7) Boyce, M. S., Vernier, P R., Nielsen, S. E., Schmiegelow, F. K. A. 2002. Evaluating resource
selection functions. Ecological Modelling, 157:281 – 300.
8) Browne, W. J., Draper, D. 2006. A comparison of Bayesian and likelihood-based methods
for fitting multilevel models. Bayesian Analysis, 1(3):473 – 514.
9) Buskirk, S. W., Millspaugh, J. J. 2006. Metrics for studies of resource selection. Journal of
Wildlife Management, 70(2):358 – 366.
10) Byers, C. R., Steinhorst, R. K., Krausman, P. R. 1984. Clarification of a technique for
analysis of utilization-availability data. Journal of Wildlife Management, 48(3):1050 – 1053.
11) Cade, B. 2015. Model averaging and muddled multimodal inference. Ecology, 96(9):2370 –
2382.
12) Charnov, E. L. 1976. Optimal foraging, the marginal value theorem. Theoretical Population
Biology, 9:129 – 136.
13) Cooper, A. B., Millspaugh, J. J. 1999. The application of discrete choice models to wildlife
resource selection studies. Ecology, 80(2):566 – 575.
14) Cooper, A. B., Millspaugh, J. J. 2001. Accounting for variation in resource availability and
animal behavior in resource selection studies In J. J. Millspaugh and J. M. Marzluff (Eds.),
Radio Tracking and Animal Populations (pp. 4 – 15). San Diego, CA: Academic Press.

81

15) Craighead, J. J., Craighead, F., C. 1971. Grizzly bear – man relationships in Yellowstone
National Park. Bioscience, 21(16):845 – 857.
16) Craighead, J. J., Craighead, F. C., Varney, J. R., Cote, C. E. 1971. Satellite monitoring of
black bear. Bioscience, 21(24):1206 – 1212.
17) Duchesnse, T., Fortin, D., Courbin, N. 2010. Mixed conditional logistic regression for habitat
selection studies. Journal of Animal Ecology, 79(3):548 – 555.
18) Fieburg, J., Matthiopolous, J., Hebblewhite, M., Boyce, M. S., Frair, J. L. 2010. Correlation
and studies of habitat selection: problem, red herring or opportunity? Philosophical
Transactions of the Royal Society: Biological Sciences, 365(1550):2233 – 2244.
19) Gaillard, J.-M., Hebblewhite, M., Loison, A., Fuller, M., Powell, R., Basille, M., van
Moorter, B. 2010. Habitat-performance relationships: finding the right metric at a given
spatial scale. Philosophical Transactions of the Royal Society: Biological Science,
365(1550), 2255-2265.
20) Gelman, A., Carlin, J. B., Stern, H. S., Dunson, D. B., Vehtari, A., Rubin, D. B. 2014a.
Bayesian Data Analysis (Third Edition). Boca Raton, Florida: CRC Press.
21) Gelman, A., Hwang, J., Vehtari, A. 2014b. Understanding predictive information criteria for
Bayesian models. Statistical Computing, 24:997 – 1016.
22) Gillies, C. S., Hebblewhite, M., Nielsen, S. E., Krawchuk, M. A., Aldridge, C. L., Frair, J. L.,
Saher, D. J., Stevens, C. E., and Jerde, C. L. 2006. Application of random effects to the study
of resource selection by animals. Journal of Animal Ecology, 75(4):887 – 898.
23) Gillingham, M. P., Parker, K. L. 2008. The importance of individual variation in defining
habitat selection by moose in northern British Columbia. Alces, 44:7 – 20.
24) Griffith, J. A., Martinko, E. A., Price, K. P. 2000. Landscape structure analysis of Kansas at
three scales. Landscape and Urban Planning, 52(1):45 – 61.
25) Hanks, E. M., Hooten, M. B., Johnson, D. S., Sterling, J. T. 2011. Velocity-based movement
modeling for individual and population level inference. PLoS ONE, 6(8):e22795.
26) Hebblewhite, M., and Merrill, E. 2008. Modelling wildlife-human relationships for social
species with mixed-effects resource selection models. Journal of Applied Ecology, 45(3):834
– 844.
27) Hobbs, N. T., Hooten, M. B. 2015. Bayesian models: a statistical primer for ecologists.
Princeton, NJ: Princeton University Press.

82

28) Hooten, M. B., Buderman, F. E., Brost, B. M., Hanks, E. M., Ivan, J. S. 2016. Hierarchical
animal movement models for population-level inference. Environmetrics, 27(6):322 – 333.
29) Hirzel, A. H., Le Lay, G. 2008. Habitat suitability modelling and niche theory. Journal of
Applied Ecology, 45(5):1372 – 1381.
30) Johnson, D. H. 1980. The comparison of usage and availability measurements for evaluating
resource preference. Ecology, 61(1):65 – 71.
31) Johnson, D. J., Thomas, D. L., Ver Hoef, J. M., Christ, A. 2008. A general framework for the
analysis of animal resource selection from telemetry data. Biometrics, 64(3):986 – 976.
32) Kays, R. Crofoot, M. C., Jetz, W., Wikelski, M. 2015. Terrestrial animal tracking as eye on
life and planet. Science, 348(6240):1222 – 1232.
33) Keating, K. A., Cherry, S. 2004. Use and interpretation of logistic regression in habitatselection studies. Journal of Wildlife Management, 68(4):774 – 789.
34) Kenward, R. E. 2001. Historical and Practical Perspectives In J. J. Millspaugh and J. M.
Marzluff (Eds.), Radio Tracking and Animal Populations (pp. 243 – 273). San Diego, CA:
Academic Press.
35) Merrick, M. J., Koprowski, J. L. 2017. Should we consider individual behavior differences in
applied wildlife conservation studies? Biological Conservation, 209:34 – 44.
36) Leclerc, M., Vander Wal, E., Zedrosser, A., Swenson, J. E., Kindberg, J., Pelltier, F. 2015.
Quantifying consistent individual differences in habitat selection. Oecologia, 180:697 – 705.
37) MacArthur, R. H., Pianka, E. R. 1966. On optimal use of a patchy environment. The
American Naturalist, 100(916):603 – 609.
38) Manly, B. F., McDonald, L. L., Thomas, D. L., McDonald, T. L., Erickson, W. P. 2002.
Resource Selection by Animals: Statistical Design and Analysis for Field Studies (Second
Edition). Dordrecht, Netherlands: Kluwer Academic Publishers.
39) Marzluff, J. M., Millspaugh, J. J., Hurvitz, P., Handcock, M. S. 2004. Relating resources to a
probabilistic measure of space use: forest fragments and Stellar’s Jays. Ecology, 85(5):1411
– 1427.
40) McCracken, M. L., Manly, B. F. J., and Vander Heyden, M. 1998. The use of discrete-choice
models for evaluating resource selection. The Journal of Agricultural, Biological, and
Environmental Statistics, 3(3):268 – 279.
41) Missouri Department of Conservation (MDC). 2010. Elk restoration in Missouri. Missouri
Department of Conservation, Jefferson City, USA.

83

42) Montgomery, R. A., and Roloff, G. J. 2013. Habitat Selection In S. A. Levin (Ed.),
Encyclopedia of Biodiversity (Second edition, Volume 4; pp. 59 – 69). Waltham, MA:
Academic Press.
43) Morris, D. W. 2003. Toward an ecological synthesis: a case for habitat selection. Oecologia,
136:1 – 13.
44) Mysterud, A., Ims, R. A. 1998. Functional responses in habitat use: availability influences
relative use in trade-off situations. Ecology 79(4):1435 – 1441.
45) Neu, C. W., Byers, C. R., Peek, J. M. 1974. A technique for analysis of utilizationavailability data. Journal of Wildlife Management, 38(3):541 – 545.
46) Paton, R. S., Matthiopolous, J. 2015. Defining the scale of habitat availability for models of
habitat selection. Ecology, 97(5):1113 – 1122.
47) Pape, R., LĂśffler, J. 2015. Ecological dynamics in habitat selection of reindeer: an interplay
of spatial scale, time, and individual animal’s choice. Polar Biology, 38:1891 – 1903.
48) Putman, R., Flueck, W. T. 2011. Intraspecific variation in biology and ecology of deer:
magnitude and causation. Animal Production Science, 51(4):277 – 291.
49) R Core Team. 2016. R: A language and environment for statistical computing. R Foundation
for Statistical Computing, Vienna, Austria. URL https://www.R-project.org/.
50) Ramsey, F. L., Usner, D. 2003. Persistence and heterogeneity in habitat selection studies
using radio telemetry. Biometrics, 3(2):332 – 340.
51) RĂŠale, D., Dingemanse, N. J., Kazem, A. J. N., Wright, J. 2010. Evolutionary and ecological
approaches to the study of personality. Philosophical Transactions of the Royal Society B.,
365:3937 – 3946.
52) Rosensweig, M. L. 1981. A theory of habitat selection. Ecology, 62(2):327 – 335.
53) Rosenzweig, M. L. 1991. Habitat selection and population interactions: the search for
mechanism. The American Naturalist, 137(Supplement: Habitat Selection), S5 – S28.
54) Thomas, D. L., Johnson, D., Griffith, B. 2006. A Bayesian random effects discrete-choice
model for resource selection: population-level selection inference. Journal of Wildlife
Management, 70(2):404 – 412.
55) Sawyer, H., Nielsen, R. M., Lindzey, F., Keith, L., Powell, J. H., Abraham, A. A. 2007.
Habitat selection of rocky mountain elk in a nonforested environment. Journal of Wildlife
Management, 71(3):868 – 874.

84

56) Sawyer, H., Nielsen, R. M., Lindzey, F., McDonald, L. L. 2006. Winter habitat selection of
mule deer before and during development of a natural gas field. Journal of Wildlife
Management, 70(2):396 – 403.
57) Smith, T. N. 2015. Broad-scale resource selection and food habits of a recently reintroduced
elk population in Missouri (Master’s thesis). Retrieved from
https://mospace.umsystem.edu/xmlui/bitstream/handle/10355/49138/research.pdf?sequence=
2andisAllowed=y
58) Spiegel, O., Leu, S. T., Bull, C. M., Sih, A. 2017. What’s your move? Movement as a link
between personality and spatial dynamics in animal populations. Ecology Letters, 20:3 – 18.
59) Stan Development Team 2016a. The Stan C++ Library, Version 2.14.0. http://mc-stan.org.
60) Stan Development Team 2016b. Rstan: the R interface to Stan, R package version 2.14.1.
http://mc-stan.org.
61) Thomas, D. L., Taylor, E. J. 1990. Study designs and tests for comparing resource use and
availability. Journal of Wildlife Management, 54(2):322 – 330.
62) Thomas, D. L., Taylor, E. J. 2006. Study designs ad tests for comparing resource use and
availability II. Journal of Wildlife Management, 70(2):324 – 336.
63) Thomasma, L. E., Drummer, T. D., Peterson, R. O. 1991. Testing the habitat suitability index
model for the fisher. Wildlife Society Bulletin, 19(3):291 – 297.
64) Thurfjell, H., Ciuti, S., Boyce, M. S. 2014. Applications of step-selection functions in
ecology and conservation. Movement Ecology, 2(1):4.
65) Vehtari, A., Gelman, A., Gabry, J. 2016. Practical Bayesian model evaluation using leaveone-out cross-validation and WAIC. Journal of Statistics and Computing,
DOI:10.1007/s11222-016-9696-4.
66) Watanabe, S. 2010. Asymptotic equivalence of Bayes cross validation and widely applicable
information criterion in singular learning theory. Journal of Machine Learning Research,
11:3571 – 3594.
67) Wiens, T. S., Dale, B. C., Boyce, M. S., and Kershaw, G. P. 2008. Three way k-fold crossvalidation of resource selection functions. Ecological Modelling, 212(3-4):244 – 255.

85

CONCLUSION

Conservation and management of ungulate species will require a mechanistic
understanding of the ways in which extrinsic and intrinsic factors change movement and habitat
selection. The ability to obtain location data on these animals at relatively fine temporal
resolutions for long periods of time suggests that the movement ecology paradigm is well-suited
for discovering these mechanistic links. When coupled with on-board data-loggers (e.g.,
temperature) and data from remote sensing platforms or autonomous sensor networks, there is
opportunity to study the influence of both extrinsic characteristics of an individual’s
surroundings, and of intrinsic properties. My study took important steps towards identifying the
underlying mechanisms while developing quantitative techniques to model these data with
increased accuracy. I envision the contributions of my thesis to be two-fold. First, in chapter 1
non-linear modeling revealed little evidence of heat stress thresholds in movement of moose in
central Norway. Using a population-level model, I was able to demonstrate that physiological
thresholds of moose established about 30 years ago may not affect moose movement as much as
had previously been thought. The results further verified the importance of studying free-living
animals whenever possible. Second, in chapter 2 I evaluated the pitfalls of the practice of
analyzing the process of resource selection at the aggregate level. Even when accounting for
individual variation by incorporating mixed-effects models (i.e., random effect by individual
animal) to fit resource selection functions, this work demonstrated that potential inferences could
still be missed. This finding has important implications for wildlife management and
conservation practice.

86

The main strength of this study is to simultaneously treatment of both the aggregate- and
individual-level perspectives, while remaining in the Lagrangian framework of treating the
individual as the sample unit, in accordance with prevailing discourse in movement ecology.
Specifically, I used movement data derived from individual tracks and associated temperature
readings from collar-borne loggers to evaluate my research question in an aggregate perspective.
Additionally, I employed a Bayesian discrete choice modelling framework to estimate resource
selection, which is predicated on the notion of individual choice and was therefore highly
appropriate for estimating choice at the individual level. In this way, it makes intuitive sense to
think about individual heterogeneity in the process of resource selection that may not be
effectively captured by aggregate methods. Admittedly, a salient weakness of my study is that I
did not account for both processes (habitat selection and movement) in either of the focal
analyses. These processes are intricately linked; habitat likely plays a major role in deciding how
and when ungulates move, and vice versa (Thurfjell et al., 2014). Thus, quantifying both
simultaneously could better inform the mechanistic roles of important intrinsic and extrinsic
dimensions.
Clearly, efforts to conserve ungulate species will be greatly benefitted by the quantitative
analysis of movement and habitat selection, in addition to other processes (Nathan et al., 2013).
However, scaling up inference to the population level is paramount for making management
decisions. Development of holistic methods to do this is an area of active research, and
computationally efficient hierarchical modeling approaches show much promise (Hooten et al.,
2016; Jonsen, 2017). However, besides the fact that many recent open-source software packages
have been developed with analysis of individual movement tracks in mind (e.g., Calabrese et al.,
2016), the notion that consistent differences among individuals (i.e., personality) could have

87

implications for management is gaining attention, which raises questions about the scenarios in
which such aggregate methods make biological or ecological sense, or are even practical
(Merrick and Koprowski, 2017; Speigel et al., 2017). Nevertheless, it remains critical to combine
the power of open-source analytical tools and mechanistic models of movement and habitat
selection with data collected from distributed sensor networks, remote sensing platforms, and
telemetered individuals, as well as intrinsic biological information. This thesis is a testament to
the need for multiple analytical perspectives when nurturing a mechanistic understanding of
ungulate space use and fitness at multiple scales.

88

REFERENCES

89

REFERENCES

1) Hooten, M. B., Buderman, F. E., Brost, B. M., Hanks, E. M., Ivan, J. S. 2016. Hierarchical
animal movement models for population-level inference. Environmetrics, 27:322 – 333.
DOI:10.1002/env.2402.
2) Jonsen, I. 2017. Joint estimation over multiple individuals improves behavioural state
inference from animal movement data. Nature, 6:20625. DOI: 10.1038/srep20625.
3) Merrick, M. J., Koprowski, J. L. 2017. Should I consider individual differences in applied
wildlife conservation studies? Biological Conservation, 209:34 – 44.
4) Nathan, R., Giuggioli, L. A milestone for movement ecology research. Movement Ecology, 1
(1):1. DOI:10.1186/2051-3933-1-1.
5) Speigel, O., Leu, S. T., Bull, C. M., Sih, A. What’s your move? Movement as a link between
personality and spatial dynamics in animal populations. Ecology Letters, 20:3 – 18.
DOI:10.1111/ele.12708.
6) Thurfjell, H., Ciuti, S., Boyce, M. S. 2014. Applications of step-selection functions in
ecology and conservation. Movement Ecology, 2:4.

90