R help - Jun 2012 - lme random effects in additive models with interaction

Hello,
I work with a mixed model with 4 predictor variables Time, Size, Charge,
Density and Size, Charge, Density are factors, all with two levels. Hence I
want to put their interactions with Time into the model. But, I have two
data sets (Replication 1 and 2) and I want that Replication is random
effect. Here is my code:


knots <- default.knots(Time)

z <- outer(Time, knots, "-") 
z <- z * (z > 0)
z<-z^2

i.size50 <- I(Size==50)
i.chargepos <- I(Charge=="+")
i.densitylow <- I(Density==0.001)


X <- cbind( I(Time^2),Time*i.size50,Time*i.chargepos,Time*i.densitylow)
Z <- cbind( z, z*i.size50, z*i.chargepos,z*i.densitylow)
K <- length(knots)

block.ind <- list(1:K, (K+1):(2*K),(2*K+1):(3*K),(3*K+1):(4*K))
Z.block <- list()
for (i in 1:length(block.ind))
Z.block[[i]] <-
as.formula(paste("~Z[,c(",paste(block.ind[[i]],collapse=","),")]-1"))

group <- rep(1, length(Time))
model.data <- groupedData(y~X|group, data=data.frame(X, y))

fit <- lme(y~-1+X, data=model.data, random=pdBlocked(list(
pdBlocked(Z.block,pdClass="pdIdent"), pdIdent(~-1+ Replication) ))
,control=list(maxIter=1000, msMaxIter=1000, niterEM=1000))


It gives errror: "Error:  getResponseFormula(el) : "Form" must be
a two
sided formula"

Does anybody help how can I write random part? 
Thanks..


--
View this message in context:
http://r.789695.n4.nabble.com/lme-random-effects-in-additive-models-with-interaction-tp4634067.html
Sent from the R help mailing list archive at Nabble.com.

elifnurdogruoz <elifnurdogruoz <at> mynet.com> writes:
> 
> Hello,
> I work with a mixed model with 4 predictor variables Time, Size, Charge,
> Density and Size, Charge, Density are factors, all with two levels. Hence I
> want to put their interactions with Time into the model. But, I have two
> data sets (Replication 1 and 2) and I want that Replication is random
> effect. Here is my code:
  Sorry not to give you a full answer to this, but I have a few 
quick responses:
 
  * although you have given us your code, this isn't fully
reproducible (see http://tinyurl.com/reproducible-000 ) --- that
isn't always completely necessary, but it can help a lot ...

  * What is default.knots() ? I can't find it anywhere ...

  *  The rest of your code is rather complicated for understanding
at a glance -- some comments (##) about what's going on might be
helpful

  * You might want to re-post this on r-sig-mixed-models <at>
r-project.org , specifying that it was previously posted here.

  * Most specifically (but alas most negatively), you may have some
serious problems with treating a two-level factor (replication) as a
random effect: in principle and philosophically, there is nothing
wrong with supposing that your replications are randomly chosen
from a hypothetical population of possibilities, but it may be
difficult both practically and statistically to get meaningful
results with n=2 (see e.g.
https://stat.ethz.ch/pipermail/r-sig-mixed-models/2010q2/003709.html
for a discussion on this topic) ...