There is probably a way to do what you want, but I don't know how.
You are to be commended for providing a self-contained example citing an
interesting article from the Journal of Statistical Software and
including the modification necessary to make the S-Plus "lme" call
work
in R.
I tried several things including the following:
fit.r1 <- lmer(y~-1+X+(-1+Z|grp))
This seemed to run but gave an answer different from that from lme.
I also noted that the lmer documentation said the argument "start" was
"a list of relative precision matrices for the random effects. This has
the same form as the slot '"Omega"' in a fitted model. Only
the upper
triangle of these symmetric matrices should be stored." This suggested
I provide "start" with a list consiting of R from the QR decomposition
of Z:
> grp.dat.qr <- qr(grp.dat[2:26])
> grp.dat.R <- qr.R(grp.dat.qr)
> fit.r <- lmer(y~-1+X+(-1+Z|grp), start=list(grp=grp.dat.R))
Error in lmer(y ~ -1 + X + (-1 + Z | grp), start = list(grp = grp.dat.R)) :
Leading 1 minor of Omega[[1]] not positive definite
Some of the diagonal elements of R were negative, so I just changed
the sign and tried again:
> grp.dat.R2 <- (diag(sign(diag(grp.dat.R)))
+ %*% grp.dat.R)
> fit.r2 <- lmer(y~-1+X+(-1+Z|grp), start=list(grp=grp.dat.R2))
Error in lmer(y ~ -1 + X + (-1 + Z | grp), start = list(grp =
grp.dat.R2)) :
Leading 2 minor of Omega[[1]] not positive definite
If I were to do more with this, I'd first review all the
documentation I could find on lmer including Doug Bates, "Fitting linear
mixed models in R", R News, 5(1): 27-30, May 2005, and "Linear mixed
model implementation in lmer", July 26, 2005, distributed with lme4 and
stored on my hard drive under
"~\R\R-2.2.0\library\lme4\doc\Implementation.pdf". I might also
consult
Pinheiro and Bates (2000) Mixed-Effects Models in S and S-PLUS
(Springer), which is my primary source for mixed models generally. If I
couldn't figure it out from there, I'd then try to work through the code
line by line. Since "lmer" calls "standardGeneric",
it's not completely
obvious how to get the code. Moreover, "methods" won't get it,
because
"lmer" follows the S4 standard (if I understand correctly).
Fortunately, 'showMethods("lmer")' produced the following:
Function "lmer":
formula = "formula"
With this information, I then tried, 'getMethod("lmer",
"formula")', which gave me the desired source code. I could then
copy
it into a script file, walk through it line by line, and learn something.
hope this helps.
spencer graves
Mark Lyman wrote:
> Is there a way to specify a Z matrix using the lmer function, where the
> model is written as y = X*Beta + Z*u + e?
>
> I am trying to reproduce smoothing methods illustrated in the paper
> "Smoothing with Mixed Model Software" my Long Ngo and M.P. Wand.
> published in the /Journal of Statistical Software/ in 2004 using the
> lme4 and Matrix packages. The code and data sets used can be found at
> http://www.jstatsoft.org/v09/i01/.
>
> There original code did not work for me without slight modifications
> here is the code that I used with my modifications noted.
>
> x <- fossil$age
> y <- 100000*fossil$strontium.ratio
> knots <- seq(94,121,length=25)
> n <- length(x)
> X <- cbind(rep(1,n),x)
> Z <- outer(x,knots,"-")
> Z <- Z*(Z>0)
> # I had to create the groupedData object with one group to fit the model
> I wanted
> grp <- rep(1,n)
> grp.dat<-groupedData(y~Z|grp)
> fit <- lme(y~-1+X,random=pdIdent(~-1+Z),data=grp.dat)
>
> I would like to know how I could fit this same model using the lmer
> function. Specifically can I specify a Z matrix in the same way as I do
> above in lme?
>
> Thanks,
> Mark
>
> ______________________________________________
> R-help at 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
--
Spencer Graves, PhD
Senior Development Engineer
PDF Solutions, Inc.
333 West San Carlos Street Suite 700
San Jose, CA 95110, USA
spencer.graves at pdf.com
www.pdf.com <http://www.pdf.com>
Tel: 408-938-4420
Fax: 408-280-7915