[R] generalized mixed linear models, glmmPQL and GLMER give very different results that both do not fit the data well...

[Bert Gunter]

I suggest that you post this on the R-sig-mixed-models list where you
are more likely to find those with bothe interest and expertise in
these matters.

-- Bert

On Thu, Mar 17, 2011 at 7:44 AM, Franssens, Samuel
<Samuel.Franssens at econ.kuleuven.be> wrote:
> Hi,
>
> I have the following type of data: 86 subjects in three independent groups (high power vs low power vs control). Each subject solves 8 reasoning problems of two kinds: conflict problems and noconflict problems. I measure accuracy in solving the reasoning problems. To summarize: binary response, 1 within subject var (TYPE), 1 between subject var (POWER).
>
> I wanted to fit the following model: for problem i, person j:
> logodds ( Y_ij ) = b_0j + b_1j TYPE_ij
> with b_0j = b_00 + b_01 POWER_j + u_0j
> and b_1j = b_10 + b_11 POWER_j
>
> I think it makes sense, but I'm not sure.
> Here are the observed cell means:
>                conflict                 noconflict
> control 0.6896552            0.9568966
> high      0.6935484            0.9677419
> low         0.8846154            0.9903846
>
> GLMER gives me:
> summary(glmer(accuracy~f_power*f_type + (1|subject), family=binomial,data=syllogisms))
> Generalized linear mixed model fit by the Laplace approximation
> Formula: accuracy ~ f_power * f_type + (1 | subject)
>   Data: syllogisms
>  AIC   BIC logLik deviance
> 406 437.7   -196      392
> Random effects:
> Groups  Name        Variance Std.Dev.
> subject (Intercept) 4.9968   2.2353
> Number of obs: 688, groups: subject, 86
>
> Fixed effects:
>                            Estimate Std. Error z value Pr(>|z|)
> (Intercept)                  1.50745    0.50507   2.985  0.00284 **
> f_powerhp                    0.13083    0.70719   0.185  0.85323
> f_powerlow                   2.04121    0.85308   2.393  0.01672 *
> f_typenoconflict             3.28715    0.64673   5.083 3.72e-07 ***
> f_powerhp:f_typenoconflict   0.21680    0.93165   0.233  0.81599
> f_powerlow:f_typenoconflict -0.01199    1.45807  -0.008  0.99344
> ---
> Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
>
> Correlation of Fixed Effects:
>            (Intr) f_pwrh f_pwrl f_typn f_pwrh:_
> f_powerhp   -0.714
> f_powerlow  -0.592  0.423
> f_typncnflc -0.185  0.132  0.109
> f_pwrhp:f_t  0.128 -0.170 -0.076 -0.694
> f_pwrlw:f_t  0.082 -0.059 -0.144 -0.444  0.308
>
> glmmPQL gives me:
> summary(glmmPQL(fixed=accuracy~f_power*f_type, random=~1|subject, family=binomial, data=syllogisms))
> iteration 1
> iteration 2
> iteration 3
> iteration 4
> iteration 5
> iteration 6
> Linear mixed-effects model fit by maximum likelihood
> Data: syllogisms
>  AIC BIC logLik
>   NA  NA     NA
>
> Random effects:
> Formula: ~1 | subject
>        (Intercept)  Residual
> StdDev:    1.817202 0.8045027
>
> Variance function:
> Structure: fixed weights
> Formula: ~invwt
> Fixed effects: accuracy ~ f_power * f_type
>                                Value Std.Error  DF  t-value p-value
> (Intercept)                 1.1403334 0.4064642 599 2.805496  0.0052
> f_powerhp                   0.0996481 0.5683296  83 0.175335  0.8612
> f_powerlow                  1.5358270 0.6486150  83 2.367856  0.0202
> f_typenoconflict            3.0096016 0.4769761 599 6.309754  0.0000
> f_powerhp:f_typenoconflict  0.1856061 0.6790046 599 0.273350  0.7847
> f_powerlow:f_typenoconflict 0.0968204 1.0318659 599 0.093830  0.9253
> Correlation:
>                            (Intr) f_pwrh f_pwrl f_typn f_pwrh:_
> f_powerhp                   -0.715
> f_powerlow                  -0.627  0.448
> f_typenoconflict            -0.194  0.138  0.121
> f_powerhp:f_typenoconflict   0.136 -0.182 -0.085 -0.702
> f_powerlow:f_typenoconflict  0.089 -0.064 -0.153 -0.462  0.325
>
> Standardized Within-Group Residuals:
>         Min           Q1          Med           Q3          Max
> -12.43735991   0.06243699   0.22966010   0.33106978   2.23942234
>
> Number of Observations: 688
> Number of Groups: 86
>
>
> Strange thing is that when you convert the estimates to probabilities, they are quite far off. For control, no conflict (intercept), the estimation from glmer is 1.5 -> 81% and for glmmPQL is 1.14 -> 75%, whereas the observed is: 68%.
>
> Am I doing something wrong?
>
> Any help is very much appreciated.
> Sam.
>
>        [[alternative HTML version deleted]]
>
> ______________________________________________
> R-help at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>



-- 
Bert Gunter
Genentech Nonclinical Biostatistics