Statistical Analyses:
All statistical analyses were conducted in R v.4.0.3 (R Core Team, 2020)
using packages climwin v1.2.3 (Bailey & van de Pol, 2016; van de Polet al., 2016), lme4 v1.1.26 (Bates et al., 2015), asreml-R
v4.1 (Butler, 2021), ggplot2 v.3.3.3 (Wickham, 2016), tidyverse v.1.3.1
(Wickham et al., 2019), nadiv v.2.17.1 (Wolak, 2012), lubridate
v.1.7.10 (Grolemund & Wickham, 2011), gridExtra v.2.3 (Auguie, 2017),
ggeffects v.1.1.2 (Lüdecke, 2018), lmerTest v3.1.3 (Kuznetsova et
al., 2017), and car 3.0.10 (Fox & Weisberg, 2019).
One of the key assumptions in a marmot’s life history is that they
reproduce immediately following their emergence. We wanted to
investigate this assumption by estimating the relationship between
female emergence from hibernation and pup emergence date from weaning in
a given year. To do so, we used a restricted dataset from the years
2003-2017 where we had all female emergence dates from hibernation and
pup emergence from weaning in each year. We then ran a linear mixed
model with the package lme4 (Bates et al., 2015) using pup
emergence date as our response variable. We used female emergence date
as the only fixed effect. Year and female identity were fitted as random
effects. Since we were expecting a 1:1 ratio between pup and female
emergence date, we also used a t-test to compare the slope of the female
emergence date to 1.
To determine during which phenological window temperature and snowpack
had the strongest association with pup emergence date, we used the R
package climwin (Bailey & van de Pol, 2016; van de Pol et al.,2016). For both environmental variables, we fit a model in climwin with
pup emergence date as our response variable and either average
temperature or snowpack as our independent variable. Our specified
reference day was June 1st, and the starting date of
the window varied from June 1st to November
13th (200 days before June 1st) the
previous year. We allowed any length of window from 1 to 200 days.
To investigate whether climate change was impacting the timing of
reproduction in the yellow-bellied marmot, we fitted a univariate animal
model of pup emergence date using the asreml-R package (Butler, 2021).
For this model, we restricted the
dataset to include those females of known age, with a known mass in
June, who reproduced that year. Fixed effects included mother’s age,
valley (up- or down-valley), average spring snowpack, average spring
temperature, and the mother’s mass in June. Random effects were the
year, permanent environment (see Kruuk & Hadfield, 2007), additive
genetic, and colony effects. Year was added to control for inter-annual
variation in conditions experienced. Permanent environment effect was
added to control for any inter-individual variation in pup emergence
date not due to genetic effects. Additive genetic effects were added to
estimate the amount of variation in the phenotype associated with
additive genetic variation. Colony was added to control for potential
micro-environmental differences between them. Significance of the random
effects was determined using a log-likelihood ratio test. Starting from
our full model, random effects were dropped one at a time and the
log-likelihood ratio of each model were compared.
To investigate inter-individual variation in the plasticity of the
timing of reproduction, we extended the previous animal model by adding
random slopes for average spring snowpack and temperature in two
separate models using asreml-R (Butler, 2021). We tested the
significance of the random slope terms and the random intercept by
comparing the log-likelihood of the models with and without each of
these effects. We did not include additive genetic effects in this
model.
We estimated the existence of selection (directional and/or stabilizing)
on pup emergence date. To do this, we ran generalized linear mixed
models with three different fitness proxies in lme4 (Bates et
al., 2015). For each litter, we used the total number of pups surviving
to one year old, the proportion of pups surviving to one year old, and
the total number of pups. We used the same model structure for the three
models. Fixed effects were the linear and quadratic orthogonal
polynomials of pup emergence date, mother’s age, mother’s mass in June,
average spring snowpack, average spring temperature, and valley. Random
effects were female identity, colony, and year. We used a Poisson
distribution for number of pups surviving to one year old and the number
of pups in the litter. We used a Binomial distribution for the
proportion of pups surviving to one year old.
For all models, all continuous variables fit as fixed effects were
scaled to have a variance of one and a mean of zero.
Results: