[R-sig-ME] Using a gamm4 model to predict estimates in new [SEC=UNCLASSIFIED]

[Ben Bolker]

Grace, Justin <justin.grace at ...> writes:


> Thanks for your comment. Unfortunately the predict function does not
> work with lme4, and therefore nor with gamm4.
 
> The issue is what to do with random effects - and the fact that some
> random effects are to be incorporated (s-terms for additive model)
> but some are not (i-terms, or random effects at the individual
> level).

The development version of lme4 _does_ have a predict() method,
and an interface for specifying for which effects conditional
predictions should be made, and for which the population-level
predictions should be made (i.e. random effects terms should
be set to zero) [there is no option to make the predictions
_marginal_ over the random effects].

The latest version of gamm4 (0.2-1) does work with the
development version of lme4, but it looks like it will take
a little more hacking to allow predictions (especially in general)
without the grouping random effects.  This is a start, but
it might be easiest to pull out Z, X, beta, random effects
and construct X%*%beta + Z%*%u yourself, with the relevant
u terms zeroed out ...

library(gamm4)
## from ?gamm4
set.seed(0) 
dat <- gamSim(1,n=400,scale=2) ## simulate data
g <- as.factor(sample(1:20,400,replace=TRUE))
dat$f <- dat$f + model.matrix(~ g-1)%*%rnorm(20)*2
dat$y <- rpois(400,exp(dat$f/7))

b2r<-gamm4(y~s(x0)+s(x1)+s(x2)+s(x3),family=poisson,
           data=dat,random = ~ (1|g))

p0 <- predict(b2r$mer)
p1 <- predict(b2r$mer,REform=~(0+s(x0)|Xr)+(0+s(x1)|Xr.0)+(0+s(x2)|Xr.1)+
              (0+s(x3)|Xr.2))
## fails