On 2/4/07, Steffen Barembruch <steffen@barembruch.de> wrote:
Hello!
Sorry to bother you again!
When fitting an LME model to my data, everything works fine. here is the
code:
lme(y ~ weeks + weeks:conditions + weeks:package, data=mydata, random=
~weeks | batch)
where y and weeks are numeric, and package, conditions and batch are
factors.
But now I played around a little for implementing an algorithm which will
use the lme function.
mydata$weeks1 <- mydata$weeks
if i now replace the effect "weeks" by "weeks1" as following
lme(y ~ weeks1 + weeks:conditions + weeks:package, data=mydata, random=
~weeks | batch)
i get an error: singularity in backsolve
furthermore
if i replace the constant term in the model by a vector of
one's i will get the same error.
The code would now look
mydata$intercept <- rep(1, dim(mydata)[1])
lme(y ~ intercept - 1 + weeks + weeks:conditions + weeks:package,
data=mydata, random= ~weeks | batch)
I don't understand that, because the design matrices look the same for all
models, don't they. And as far as I understood the maximum likelihood
estimation of these parameters, the design matrices are the only thing that
matters.
Did you check if the model matrices for the different models you tried
to fit are indeed the same? I don't think they will be.
In the S language the conversion from a formula and data set to a
model matrix goes through several stages. At one stage there is a
symbolic analysis of the formula to choose the set of contrasts used
to represent a term. In the symbolic analysis there is a big
difference between weeks + weeks:conditions and weeks1 +
weeks:conditions. In the first expression it can be determined that
the main effect for one of the factors in the interaction also occurs
as a main effect. That influences the choice of contrasts used to
represent the interaction. In the second case it seems that neither
weeks nor conditions are present as main effects so a different set of
contrasts will be used, resulting in a singularity.
Similarly, there is no way for the formula language to know that the
variable 'intercept' happens to be a constant so it blithely adds
another constant to the model matrix, resulting in yet another
singularity.
I'm not sure what extractors can be used on the S-PLUS version of an
lme fitted model so I would recommend that you compare the model
matrices for a linear model, just to see the differences.
Try
model.matrix(lm(y ~ intercept - 1 + weeks + weeks:conditions +
weeks:package, mydata))
and your other variations to see what the software actually does do,
rather than guessing at what you think it should do.
If you were using R I would suggest starting with
mm <- model.matrix(lm(y ~ intercept - 1 + weeks + weeks:conditions +
weeks:package, mydata))
str(mm)
to see the structure of that matrix, which could be rather large. I
don't know if S-PLUS has a function similar to str() to provide a
concise description of the structure of an object.
I'd appreciate any hints on what S-plus, Version 7.06 on Windows XP, is
doing.
Thanks a lot
Steffen Barembruch
|