[R-sig-ME] Specifying correlation structure in lme

Ben Bolker bbolker at gmail.com
Tue May 8 19:22:34 CEST 2012


Angelina Mukherjee <angelina.mukherjee88 at ...> writes:


> I'm trying to use AR1 as a correlation structure in modeling data. My
> variables, Patient and Region, are factors with Region nested in Patient.
> Expression is my Response variable.
> 
> *Patient <- c(rep('1', 10), rep('2', 12))*
> 
> *Region <- c(rep('1', 3), rep('2', 3), rep('3', 4), rep('4', 4), rep('5',
> 4), rep('6', 4))*
> library(nlme)

  Are you really using character variables?  I don't think that works
(although other aspects of your coding appear to have been mangled
in transition, so maybe that's another one).

I can replicate your error as follows:

df <- data.frame(Patient=factor(rep(1:2,c(10,12))),
                 Region=rep(1:6,c(3,3,4,4,4,4)),
                 Expression=rnorm(22))

library(nlme)
m0 <- lme(Expression ~ Patient + Region,
    random = ~1 |Patient/Region,
    correlation=corAR1(form=~1 | Patient/Region), data=df)

lme(Expression ~ Patient  + Region,
    random = ~1 | Patient/Region,
    correlation=corAR1(form=~Region | Patient/Region), data=df)

You're basically stuck here, the constraint on gls models is
as suggested by the error message -- covariate values must be
unique within grouping factors.

By the way, there are several other potential issues here (some
of which I think I've pointed out previously)

 * I don't remember whether you said you really only have two
patients; in that case treating patient as a random effect is
unlikely to work very well.
 * In any case, including a categorical predictor (such as Patient)
as both a random and a fixed effect probably leads to an
unidentifiable model (intervals(m0) indicates an non-positive-definite
variance-covariance matrix, a sign of trouble)
 * ideally I would probably prefer to use a random effect of
interaction(Patient,Region), and put correlations on Region
within Patient, but that is really beyond the scope of lme.

  Accepting a cruder model, or finding a local statistical
consultant, may be your choices ...



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