Skip to content

lme4 and lmeSplines

4 messages · Kevin Wright, Spencer Graves, Douglas Bates +1 more

Original
Kevin Wright
# 
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
6 days later
Spencer Graves
# 
At least at one time, to use lmer after lme, you had to quit R 
before loading 'lme4' because of substantive conflicts between 'nlme' 
and 'lme4' / 'Matrix'.  That may be the source of your problem.  I can 
think of two possible ways to get around this:

	  1.  I might try quitting R and restarting after your use of 
'lmeSplines' but before 'lmer'.

	  2.  If that failed, I might try making local copies of everything I 
needed from 'lmeSplines' and using them with 'lmer'.  If that failed, 
I'd use 'debug' to walk through the code (or local copies of whatever I 
needed, possibly with name changes) until I figured it out.

	  I know this is not what you wanted to hear, but I hope it helps.
	  Spencer Graves
Kevin Wright wrote:
Douglas Bates
#
On 8/2/06, Kevin Wright <kwright68 at gmail.com> wrote:
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.
4 days later
Rod Ball
# 
The source of the problem is that lmeSplines requires a set of iid random 
effects corresponding to a transformed smoothing spline basis (`pdIdent' 
structure in nlme) i.e. pdIdent(~Zs -1),  where Zs is the transformed set of 
spline basis functions. Kevin's model does not work because each of his 
random effects has a separately estimated variance parameter, equivalent to a 
`pdDiag' structure in nlme. lmeSplines really does require nlme.

There is currently no way in lme4 to get the equivalent of a pdIdent 
structure, i.e. set of iid random effects with common variance parameter (see 
mesage from Doug Bates below), unless these correspond to a grouping factor 
as in (~1|row). There are a number of other nlme models e.g. genetics models, 
that cannot be fitted in lme4, examples are in my message to Doug Bates 
below.

Regards,
Rod Ball