Mike Harwood <harwood262 <at> gmail.com> writes:
> Is the labeling/naming of levels in the documentation for the
> predict.glmmPQL function "backwards"? The documentation states
"Level
> values increase from outermost to innermost grouping, with level zero
> corresponding to the population predictions". Taking the sample in
> the documentation:
>
> fit <- glmmPQL(y ~ trt + I(week > 2), random = ~1 | ID,
> family = binomial, data = bacteria)
>
> > head(predict(fit, bacteria, level = 0, type="response"))
> [1] 0.9680779 0.9680779 0.8587270 0.8587270 0.9344832 0.9344832
> > head(predict(fit, bacteria, level = 1, type="response"))
> X01 X01 X01 X01 X02 X02
> 0.9828449 0.9828449 0.9198935 0.9198935 0.9050782 0.9050782
> > head(predict(fit, bacteria, type="response")) ## population
prediction
> X01 X01 X01 X01 X02 X02
> 0.9828449 0.9828449 0.9198935 0.9198935 0.9050782 0.9050782
>
> The returned values for level=1 and level=<unspecified> match, which
> is not what I expected based upon the documentation.
Well, the documentation says: "Defaults to the highest or innermost level
of
grouping", which is level 1 in this case -- right?
> Exponentiating
> the intercept coefficients from the fitted regression, the level=0
> values match when the random effect intercept is included
Do you mean "is NOT included" here?
0.9680779 (no random effect, below) matches the level=0 prediction above
0.9828449 (include random effect, below) matches the level=1 prediction,
which is also the default prediction, above.
>
> > 1/(1+exp(-3.412014)) ## only the fixed effect
> [1] 0.9680779
> > 1/(1+exp(-1*(3.412014+0.63614382))) ## fixed and random effect
intercepts
> [1] 0.9828449
This all matches my expectations. If your expectations still go
in the other direction, could you explain in more detail?
By the way, I recommend r-sig-mixed-models at r-project.org for
mixed model questions in general ...
Ben Bolker