[R-sig-ME] Mixed effect logistic regression help

[Douglas Bates]

On Tue, Mar 29, 2011 at 4:19 AM, CJ Griffiths
<Christine.Griffiths at bristol.ac.uk> wrote:
> Dear all,
>
> I want to specify a mixed conditional logistic regression to model
> microhabitat selection, but am unsure whether my dataframe and model are
> correct. I want to compare parameters such as wind and temperature at the
> location of the animal (1) to a random observation, where the animal was
> absent (0). Each Response (1/0) is thus paired by the variable Micro. To
> account for this pairing, I specified the random effect as 1|Micro.
> However, I repeatedly sampled 11 animals (Ind). Random effect = 1|Ind
>
> Is this the correct way to specify my random effects?
> Have I labelled the “Ind” correctly? Originally I had each
> Micro (pair of observations) labelled, i.e. where there was no animal I
> still identified the “Ind” which that sample was collected
> for.
>
> My data looks like this:
>
> 'data.frame':   100 obs. of  7 variables:
> $ Ind: Factor w/ 11 levels "a","b","c","d",..: 11 2 11 5 11 7 11 9 11 10 ...
>  $ Response    : int  0 1 0 1 0 1 0 1 0 1 ...
>  $ Micro   : int  1 1 7 7 8 8 9 9 10 10 ...
>  $ Slope   : int  12 0 26 2 24 2 23 4 30 2 ...
>  $ Wind    : int  4 4 3 3 2 1 2 2 3 2 ...
>  $ Temp    : num  24.2 24.1 25.9 25.3 26.6 ...
>  $ Cover   : int  0 60 2 25 5 70 20 90 30 95 ...
>
> head(dataset)
>   Ind Response Micro Slope Wind Temp Cover
> 1   na    0       1    12    4   24.2     0
> 2    b    1       1     0    4   24.1    60
> 3   na    0       7    26    3   25.9     2
> 4    e    1       7     2    3   25.3    25
> 5   na    0       8    24    2   26.6     5
> 6    g    1       8     2    1   29.2    70
>
> MICRO<-as.factor(Micro)
> y=cbind(Response,1-Response)
>
> lmer(y~Temp+Wind+Slope+Cover+(1|MICRO)+(1|Ind),data=dataset,family=binomial,REML=0)
>
> Generalized linear mixed model fit by the Laplace approximation
> Formula: y ~ Temp + Wind + Slope + Cover + (1 | MICRO) + (1 | Ind)
>   Data: dataset
>   AIC   BIC logLik deviance
>  17.39 35.63 -1.696    3.391
> Random effects:
>  Groups   Name        Variance   Std.Dev.
>  MICRO    (Intercept) 5.6244e-12 2.3716e-06
>  Ind      (Intercept) 7.8395e+03 8.8541e+01
> Number of obs: 100, groups: MICRO, 50; Ind, 11
>
> Fixed effects:
>             Estimate Std. Error z value Pr(>|z|)
> (Intercept)  14.49702  194.79843   0.074    0.941
> Temp          0.03105    5.76866   0.005    0.996
> Wind         -0.20947   15.36092  -0.014    0.989
> Slope        -0.03506    2.19051  -0.016    0.987
> Cover         0.02349    0.54182   0.043    0.965
>
> Correlation of Fixed Effects:
>      (Intr) Temp   Wind   Slope
> Temp  -0.876
> Wind  -0.506  0.498
> Slope -0.143  0.071  0.086
> Cover -0.112  0.060 -0.353 -0.173
>
> According to the results from the above model, none of the variables have
> a significant influence on the probability of finding an animal. Yet, when
> I plot the data for say cover, there is a clear trend with animals being
> more readily associated with high vegetation cover. Surely this should
> result in a greater slope and a significant p-value.
>
> The occurrence of this trend makes me suspect that the random effect
> 1|MICRO is not actually resulting in the comparison of parameters for 1 to
> 0 at a particular paired site. I would be grateful for confirmation and
> advice on my model specification.

One quick comment.  The estimated standard deviation of the MICRO
random effects is essentially zero and the fact that it is not exactly
zero is a consequence of the optimizer being used.  You could refit
the model without the random effects for MICRO.

A minor point in your usage of the function - you only need to create
the two column matrix for the response in the model (the
cbind(Response,1-Response) construction) when you don't have a 0/1
encoding.  You could put Response on the left hand side of the
formula.