[R] GLMM with lme4 and octopus behaviour

[Ben Bolker]

Samaritan <d.harliwich <at> gmail.com> writes:
 
> [snip]

  You might want to ask follow-up questions on the R-sig-mixed-models list

> At the most basic level, I'm testing the effect of sleep deprivation on
> various behaviours (e.g. amount of time spent awake, amount of time spend
> expressing difference textures, patterns, colours etc). I take a video
> sample of each octopus every hour for 72 hours and score their behaviour in
> Jwatcher.
> 
> I'm using GLMM so that I can nest Individual as a random factor within Time,
> which I'm told will reduce the problem of making repeated measures (is this
> effectively blocking by time and by octopus?). Currently my model looks like
> this:
> 
> octopus.lmer<-lmer(awake~as.factor(Treatment)+Sex+Weight+(1|Time/Octopus))

  You say you want to use GLMM -- presumably awake is a binary variable
that you want to treat as such? If so, you need the argument
'family=binomial' in your model. (You might want to use the 'glmer'
function instead, for clarity, although in practice R takes care of
this for you.)

> 
> where "Treatment" is "0" or "1" representing sleep-deprived or
> sleep-allowed. 

  As a general practice you should probably code your data as
"deprived"/"allowed" rather than 0/1: you won't have to use as.factor
and you will automatically be able to keep track of the coding.

When I try to fit the model I get the following error
> messages:
> 
> Error: length(f1) == length(f2) is not TRUE
> In addition: Warning messages:
> 1: In Octopus:Time :
>   numerical expression has 190 elements: only the first used
> 2: In Octopus:Time :
>   numerical expression has 190 elements: only the first used
> 
> This is the point at which I become lost. What does this mean? Clearly I'm
> not doing something right, perhaps in my data preparation? So far as I can
> see the length of each of the variables is the same (although I'm no certain
> as to what f1 and f2 refer to).

  This is coming from inside lmer.  My guess is that you want
to make sure Time is a factor.

dat <- expand.grid(Time=1:5,Octopus=1:5)
dat$awake <- sample(0:1,prob=c(0.5,0.5),replace=TRUE,size=nrow(dat))

dat$Time <- factor(dat$Time)
dat$Octopus <- factor(dat$Octopus)

library(lme4)

lmer(awake~(1|Time/Octopus),data=dat)

  I got a similar error when both Time and Octopus were numeric.
When I turned one or the other but not both into a factor I got
different errors.  When they are both factors I get an error which
is related to the fact that I don't have enough data (only one
observation per block, so the block effects are confounded with
the residual variation)

lmer(awake~(1|Time)+(1|Octopus),data=dat)

which represents a crossed effect, does work.

> If anyone could offer me some kind of advice about this I would appreciate
> it very, VERY much. Both of my supervisors have no experience with R and so
> have kind of washed their hands, so I'm alone in this and your expertise
> would be a big help.
> 
> Dean

  Even though R is wonderful, if your supervisors use a different
statistical package, I would strongly recommend you stick to what 
they use so that you can get properly trained,  unless you are *very*
self-sufficient.