'glmmTMB with Ornstein-Uhlenbeck autocorrelation. Different autocorrelations for different phases of the data
I promise this is going to be simple and easy to read. Or else please shoot me dead in the comments section.
I've got a dependent variable which has temporal autocorrelation, which I would like to analyse with glmmtmb. Here's the tricky part, there are times when the dependent variable is MORE autocorrelated, and times when it is LESS autocorrelated.
Lets make up some data.
library(ggplot2)
# function to get time series of autocorrelated data. high autocorrelation = low value for diV parameter
yfunct = function(diV){
y <- c(rnorm(1,0,sd=1)) # Mean=0, standard deviation=1
for(i in 1:(100-1)) { # Number of the observations is 1000
y[i+1]<-( y[i] + rnorm(1,0,sd=1))/diV #
}
y
}
# Create three groups with 6 phases. 3 high autocorrelation and 3 low autocorrelation
set.seed(11)
y = replicate (3, as.vector(replicate ( 3, c(yfunct(1), yfunct( 1.4)) ))) # high autocorrelation (divide by 1) followed by low autocorrelation (divide by 1.4)
# combine into dataframe.
df = data.frame ( y = as.vector ( y) ,
t = rep ( 1:600,3),
group = as.factor ( rep(1:3,each = 600)),
phase = rep(c("lowAC","highAC"),each=100))
# plot
ggplot ( df, aes ( x=t , y= y, color= group))+
geom_line()+
geom_vline( xintercept = seq( 0,600,100))
OK so now I'll prepare for and perform the statistical analysis. First adding a dummy independent variable, and also strip out some of the rows. The dataset I'm (really) working with has lots of missing rows. This is why I'm using Ornstein-Uhlenbeck autocorrelation. https://cran.r-project.org/web/packages/glmmTMB/vignettes/covstruct.html
My method for dealing with the autocorrelation differences between the phases is to combine group number with phase type ("high" or "low" autocorrelation), and then put this combined variable in the ou() autocorrelation function.
Is this correct? Is there something I haven't considered? Your wisdom will be hugely appreciated.
Dan
library(glmmTMB)
# set the seed
set.seed(11)
# add dummy independent variable
df$x = runif ( nrow( df))
# subsample the dataset
df = df [ sample ( 1:nrow(df), 300),]
# combine group and phase
df$comb = apply ( cbind ( df$group , df$phase) ,1, function(x)paste(x, collapse = "-"))
# numFactor the times
df$t = numFactor(df$t)
# y exists between 0 and 1 for glmmTMB
df$y = df$y + abs( min(df$y))+0.001
min(df$y)
df$y = df$y / (max(df$y)+0.001)
max(df$y)
# Model
mod1 <- glmmTMB(y~x
+ (1|group)
+ ou(t + 0 |comb)
, data=df
, beta_family(link="logit"))
# Summary
summary(mod1)
Sources
This article follows the attribution requirements of Stack Overflow and is licensed under CC BY-SA 3.0.
Source: Stack Overflow
| Solution | Source |
|---|

