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.
Kevin Wright
______________________________________________
[email protected] 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.
