Re: [R] lme4 and lmeSplines

From: Douglas Bates <bates_at_stat.wisc.edu>
Date: Thu 10 Aug 2006 - 04:16:28 EST

On 8/2/06, Kevin Wright <kwright68@gmail.com> wrote:
> I'm trying to use the lmeSplines package together with lme4.
>
> Below is (1) an example of lmeSplines together with nlme (2) an
> attempt to use lmeSplines with lme4 (3) then a comparison of the
> random effects from the two different methods.
>
> (1)
>
> require(lmeSplines)
> data(smSplineEx1)
> dat <- smSplineEx1
> dat.lo <- loess(y~time, data=dat)
> plot(dat.lo)
> dat$all <- rep(1,nrow(dat))
> times20 <- seq(1,100,length=20)
> Zt20 <- smspline(times20)
> dat$Zt20 <- approx.Z(Zt20, times20, dat$time)
> fit1.20 <- lme(y~time, data=dat, random=list(all=pdIdent(~Zt20-1)))
> # Loess model
> dat.lo <- loess(y~time, data=dat)
> plot(dat.lo)
> # Spline model
> with(dat, lines(fitted(fit1.20)~time, col="red"))
> # Save random effects for later
> ranef.nlme <- unlist(ranef(fit1.20))
>
> (2) Now an attempt to use lme4:
>
> library(lmeSplines)
> detach(package:nlme)
> library(lme4)
> data(smSplineEx1)
> # Use 20 spline in lme4
> dat <- smSplineEx1
> times20 <- seq(1,100,length=20)
> Zt20 <- smspline(times20)
> dat <- cbind(dat, approx.Z(Zt20, times20, dat$time))
> names(dat)[4:21] <- paste("Zt",names(dat)[4:21],sep="")
> dat$all <- rep(1, nrow(dat))
> fit1.20 <- lmer(y~time
> +(-1+Zt1|all)+(-1+Zt2|all)+(-1+Zt3|all)+(-1+Zt4|all)+(-1+Zt5|all)+(-1+Zt6|all)
> +(-1+Zt7|all)+(-1+Zt8|all)+(-1+Zt9|all)+(-1+Zt10|all)+(-1+Zt11|all)+(-1+Zt12|all)
> +(-1+Zt13|all)+(-1+Zt14|all)+(-1+Zt15|all)+(-1+Zt16|all)+(-1+Zt17|all)+(-1+Zt18|all),
> data=dat)
> #summary(fit1)
> # Plot the data and loess fit
> dat.lo <- loess(y~time, data=dat)
> plot(dat.lo)
> # Fitting with splines
> with(dat, lines(fitted(fit1.20)~time, col="red"))
> ranef.lme4 <- unlist(ranef(fit1.20))
>
> (3) Compare nlme lme4 random effects
>
> plot(ranef.nlme~ranef.lme4)
>
> The plot of fitted values from lme4 is visually appealing, but the
> random effects from lme4 are peculiar--three are non-zero and the rest
> are essentially zero.
>
> Any help in getting lme4 + lmeSplines working would be appreciated.
> It is not unlikely that I have the lmer syntax wrong.

It is not surprising that you get different answers from lme and lme4 because you are fitting different models. The variances of the random effects for the spline basis in the lme fit are constrained to be equal. In the lmer fit they are not constrained to be equal. It is interesting that you get all but three of the variances essentially zero. That means that there are only three active components in your spline basis, out of 20, for the fit.

I exchanged some mail off-list with Rod Ball about the definition of the random effects needed for lmeSplines and we concluded that the current capabilities in lmer are not sufficiently flexible to use it for lmeSplines. However, the sources for lmer are freely available and any enterprising programmer who would like to use the components for a more flexible model is welcome to do so.

The point is that the tools are available in lmer to represent a mixed-effects model and to evaluate the log-likelihood or restricted log-likelihood from such a model very efficiently. To optimize a model such as that being used in the lme fit one needs to go from the parameters to the model representation. This is the part that would need to be written and after that you could hook into existing code.

> Kevin Wright
>
> ______________________________________________
> R-help@stat.math.ethz.ch 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.
>



R-help@stat.math.ethz.ch 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. Received on Thu Aug 10 04:23:11 2006

Archive maintained by Robert King, hosted by the discipline of statistics at the University of Newcastle, Australia.
Archive generated by hypermail 2.1.8, at Thu 10 Aug 2006 - 06:24:27 EST.

Mailing list information is available at https://stat.ethz.ch/mailman/listinfo/r-help. Please read the posting guide before posting to the list.