I want to compare diving behaviour between penguin populations breeding at different locations using generalized linear mixed effects models (glmer). The data looks like this:
# A tibble: 12 × 7
# Groups: ID [1]
ID begdesc bout maxdep island phase diurnal
<chr> <dttm> <dbl> <dbl> <chr> <chr> <chr>
1 Nelson_Bro1_HP01 2019-01-01 19:53:11 1 52.0 Nelson Chick day
2 Nelson_Bro1_HP01 2019-01-01 20:08:18 2 5.34 Nelson Chick day
3 Nelson_Bro1_HP01 2019-01-01 20:14:39 2 52.0 Nelson Chick day
4 Nelson_Bro1_HP01 2019-01-01 20:24:46 3 64.1 Nelson Chick day
5 Nelson_Bro1_HP01 2019-01-01 20:28:44 3 75.5 Nelson Chick day
6 Nelson_Bro1_HP01 2019-01-01 20:39:44 4 68.5 Nelson Chick day
7 Nelson_Bro1_HP01 2019-01-01 20:46:58 4 62.8 Nelson Chick day
8 Nelson_Bro1_HP01 2019-01-01 20:52:19 4 62.0 Nelson Chick day
9 Nelson_Bro1_HP01 2019-01-01 20:56:43 4 60.7 Nelson Chick day
10 Nelson_Bro1_HP01 2019-01-01 20:59:56 4 62.1 Nelson Chick day
11 Nelson_Bro1_HP01 2019-01-01 21:04:05 4 62.8 Nelson Chick day
12 Nelson_Bro1_HP01 2019-01-01 21:10:05 4 51.5 Nelson Chick day
ID= Individual penguin, begdesc = date and time of dive, bout = bout number, maxdep = dive depth, Island, breeding phase (incubation or chick-rearing) and diurnal = day/night.
The response variable is maximum dive depth and the fixed effects are Island ('Kop' or 'Nel') and Diurnal (Day or Night). I include the ID of each penguin as a random effect to account for the repeated measures. The diving depths have a Gamma error distribution and each row represents a single dive. The times of dives are irregularly spaced. I used the following syntax for the generalized linear mixed models.
glmer(maxdepth ~ island + diurnal+ (1|ID), data = dive.stats, family = Gamma(link = "log")
I realized that the dives closer in time are highly autocorrelated. ACF plot without an important time variable
So now I also have a variable 'bout' that group dives that are closer together in time in the same group and dives that are separated by more than 5 minutes into a different dive bout so that the timing of dives is accounted for. I account for repeated measures of bouts within individuals as a nested random effect:
glmer(maxdepth ~ island + diurnal+ (1|ID/bout), data = dive.stats, family = Gamma(link = "log")
Once I account for bouts of dives within individuals, the ACF plot improved a lot. ACF plot with a variable that accounts for time.
But I still want to include an autocorrelation structure in the glmer model to improve the model estimates. I have looked online for how to include autocorrelation structures in a generalized linear mixed effects model, but I can't find any good answers. Has anyone been able to do this before?
Other options I explored are to use 'lme' instead, where I can add autocorrelation structures easily, but that also means that I can't specify that the data is Gamma distributed. So I will need to log transform the diving depth variable so that the errors are normally distributed, which has been discouraged before because interpretation is not easy to do with transformed data. The ACF plots look even better when I use 'lme' compared to the glmer models.
I've also tried 'glmmTMB' where I can specify the Gamma (link: log) error distribution, but I am struggling to figure out how to include the AR1 correlation structure. Can I specify 'bout' as the time factor needed for the ar1 correlation structure, since groups of dives within 5 min of each other in the same bout? Or should I try and include another time variable, e.g. hour when dive took place?
According to this vignette, if I have irregular spaced dives I should rather use the ou (Ornstein-Uhlenbeck) correlation structure, but then I need to include the locations of dives?
I would appreciate any suggestions you might have!
