[R-sig-ME] simulate lme function value

Ben Bolker bbolker at gmail.com
Wed Dec 4 19:46:08 CET 2013


Andrzej Galecki <agalecki at ...> writes:

> 
> Hello Bonitta,
> 
> You may want to  try:
> 
> library(nlme)
> fm1 <- lme(distance ~ age, data = Orthodont) # random is ~ age
> 
> library(nlmeU)
> simY <- simulateY(fm1, nsim =5)
> 
> Internallly, simulateY  employs  getVarCov.lme method  function from nlme
> package. For this reason it will not work for models with multiple levels
> of nesting.
> 
> More details about simulateY function in Section 20.4 of the  Galecki,
> Burzykowski book "Linear Mixed-Effects Models Using R: A step-by-step
> approach" (2013).

  For what it's worth, if you're willing to move to lme4 (and I
know there are reasons not to do this), simulate() works more as 
one would expect from the ?simulate help page (simulate.lme was
written *before* the generic simulate method was introduced into R),
and you can simulate either from a fitted model or from a one-sided
formula (with specified data and parameters).  Of course, you can't
simulate models with R-side effects such as correlation or
heteroscedasticity, but apparently these can't be simulated with 
simulate.lme either:

   if (length(fit1$modelStruct) > 1) {
        stop("models with \"corStruct\" and/or
     \"varFunc\" objects not allowed")
    }


e.g.

library(lme4)
data(Orthodont,package="nlme")
fm1 <- lmer(distance ~ age+(1|age), data = Orthodont)
ss <- simulate(fm1,nsim=5)

or

## this requires the development version of lme4 ...
ss2 <- simulate(~age+(1|age),
                newdata=Orthodont,
                family=gaussian,
                newparams=list(beta=c(1,1),theta=1,sigma=1),
                nsim=5)



More information about the R-sig-mixed-models mailing list