[R-sig-ME] Response variable with zero variance

[Joshua Wiley]

On Sat, Aug 11, 2012 at 6:31 PM, Ben Bolker <bbolker at gmail.com> wrote:
> Tonya Haff <tonyahaff at ...> writes:
>
>> I used the lmer package, using the following code
>>
>> FL.glmm <- lmer(RESPONSE ~ fEXPN + fPBN + fEXPN*fPBN + (1|fNID), data =
>> FL.1, family = binomial)
>
>   By the way, the fixed effect part of the formula can also
> be written as fEXPN*fPBN alone (since A*B is equivalent to A+B+A:B;
> unlike in some other stats packages, * means "crossed" rather
> than "interaction")
>
>    It may be worth double-checking this with glmmPQL:
>
>   library(MASS)
>   glmmPQL(RESPONSE~fEXPN*fPBN,random=~1|fNID,data=FL.1,
>     family="binomial")
>
> although if the data are Bernoulli (0/1) then PQL is a little
> bit questionable.
>
>>
>> Here is the summary
>>
>> Generalized linear mixed model fit by the Laplace approximation
>> Formula: RESPONSE ~ fEXPN + fPBN + fEXPN * fPBN + (1 | fNID)
>>    Data: FL.1
>>    AIC   BIC logLik deviance
>>  104.1 163.0 -35.04    70.09
>> Random effects:
>>  Groups Name        Variance Std.Dev.
>>  fNID   (Intercept) 0.047168 0.21718
>> Number of obs: 236, groups: fNID, 15
>>
>> Fixed effects:
>>                Estimate Std. Error z value Pr(>|z|)
>> (Intercept)     -0.4102     0.5304  -0.774    0.439
>> fEXPN2           1.4329     0.7880   1.818    0.069 .
>> fEXPN3          21.9715 12840.4523   0.002    0.999
>> fEXPN4          21.9715 12409.4403   0.002    0.999
>> fPBN2           21.9715 12409.4408   0.002    0.999
>> fPBN3            0.5453     0.7390   0.738    0.461
>> fPBN4           -1.4788     0.9286  -1.593    0.111
>> fEXPN2:fPBN2    -1.4171 17619.4645  -1e-04    1.000
>> fEXPN3:fPBN2   -21.9715 21994.2881  -0.001    0.999
>> fEXPN4:fPBN2   -21.9715 21493.7828  -0.001    0.999
>> fEXPN2:fPBN3    20.0155 12547.8259   0.002    0.999
>> fEXPN3:fPBN3    -0.5390 18188.0409   0.000    1.000
>> fEXPN4:fPBN3    -0.5390 17577.5297   0.000    1.000
>> fEXPN2:fPBN4   -21.1113 12447.6988  -0.002    0.999
>> fEXPN3:fPBN4   -41.6658 18265.7721  -0.002    0.998
>> fEXPN4:fPBN4   -41.6660 17647.6445  -0.002    0.998
>
>   This looks like a Hauck-Donner effect to me (you can look it
> up -- it occurs in GLM(M)s when there are strong binomial
> effects).  I'm a little worried about the exact equivalence
> of some of the parameter estimates though, and apparent
> overfitting ... if these are Bernoulli data, how many
> positive responses are there overall?  Have you tried
> nAGQ=8 ?
>
>> So here I had no F statistic at all, and the results of the fixed effects
>> don't make sense, when compared against the actual data (where the pattern
>> between PBN (playback type) and EXPN (experiment) are quite stark.
>>
>> Because this didn't make sense to me, and because I am not a whiz at R, I
>> ran the same model in SPSS 20.0, and that is where I got the huge F
>> statistic. Here is that model summary:
>>
>>  F = 1, 229, 656.019, df1 = 15, df2=220, P<0.0001
>> EXPN F = 83.83, df1=3, df2=220, P<0.0001
>> PBN F=2,2280.08, df1=3, df2=220, P<0.0001
>> EXPNxPNB F=82.96, df1=9, df2=220, P<0.0001
>
>   But: what did you do in SPSS, which doesn't fit GLMMs (as
> far as I know)?

SPSS does fit GLMMs (using PQL, I believe).  It is not a pleasant
interface, but here are examples:

http://www.ats.ucla.edu/stat/spss/code/hdp.htm

which are meant to analyze this (simulated) dataset:

http://www.ats.ucla.edu/stat/r/pages/mesimulation.htm

although I had to give up the right censored mixed effects poisson
(just ran as poisson in SPSS) and the mixed effects beta (ran as
binomial in SPSS).

>
>> >
>> > Ideally, I would like to run a GLMM on my data, which are are binomial (so
>> > I'm trying to run a model with a binomial distribution and a logit link).
>> >  I have two fixed effects, fledgling age and playback type.  I have
>> > fledging ID as a random effect.  The problem I'm running into is zero
>> > variance in my response variable - ie fledglings always respond the same
>> > way to some of the playback types.  I can run a model and it converges, but
>> > it spits out an F value of over a million. So I'm wondering (1) is it
>> > possible to use a GLMM with a zero variance response; and (2) if so, which
>> > R package is the most appropriate?
>> >
>> >
>
> _______________________________________________
> R-sig-mixed-models at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models



-- 
Joshua Wiley
Ph.D. Student, Health Psychology
Programmer Analyst II, Statistical Consulting Group
University of California, Los Angeles
https://joshuawiley.com/