[R-sig-ME] Specifying correlation structure

Gang Chen gangchen at mail.nih.gov
Tue May 15 12:27:52 CEST 2012


I have some time series data (with 100 time points) collected from one
subject. At each time point three measures were taken under three
conditions, and there is an explanatory variable for each condition.

Below is just some fake data to demonstrate the data structure:

# 100 time points per condition
tp <- 100                             # time points
trials <- seq(1, tp, 1)          # time counter
condition <- c('a', 'b', 'c')   # three conditions

# fake data
set.seed(5)
Dat <- data.frame(time=rep(trials, 3), cond=rep(condition, each=tp),
res=rnorm(3*tp))

reg1 <- c(rnorm(tp), rep(0, 2*tp))            # explanatory variable
for condition 1
reg2 <- c(rep(0, tp), rnorm(tp), rep(0, tp))   # explanatory variable
for condition 2
reg3 <- c(rep(0, 2*tp), rnorm(tp))             # explanatory variable
for condition 3

First I thought that I'd start with gls in nlme package with an ARMA
model for the correlation structure:

(fm <- gls(res ~ 1+reg1+reg2+reg3, correlation=corARMA(c(0.02, 0.03),
form=~time, p=1,q=1), data=Dat))

The above model does not work because of the following error:

Error in Initialize.corARMA(X[[1L]], ...) :
  Covariate must have unique values within groups for corARMA objects

Moreover, I would like to account for the fact that the ARMA structure
should be similar or the same across the three conditions. My
questions are:

1) how to impose a same ARMA structure across the three conditions?
2) alternatively I'd like to have something like:

(fm <- gls(res ~ 1+reg1+reg2+reg3, correlation=corARMA(c(0.02, 0.03),
form=~time:condition, p=1,q=1), data=Dat))

but this would not work either, and seems to crash R!

3) maybe construct a multivariate gls model, but how to do that? what package?

Thanks,
Gang



More information about the R-sig-mixed-models mailing list