[R] nlme cross-over and fixed nested

[Dieter Menne]

I have problem getting the concept of a nested fixed variable into the nlme
scheme. I fear the question is very stupid.  In the past I had asked this
before, and never got a reply (in other cases, the response was within
hours). I also checked the S-list, where several similar enquiries of other
people are orphaned.

We have a cross-over design, where patient are treated two weeks with
placebo, two weeks with a drug ("verum"), and a washout phase in between.
For each phase, there is ONE number (let's call it "glucose level")
describing the treatment effect. Blind treatement sequence is Seq=1 for
placebo/verum, Seq=2 for verum/placebo.

In Pinheiro/Bates, cross-over trials are mentioned on page 107, but the
parameter "Side" in the example is not exactly what we have.
We used the scheme below to analyze the data (there is a small effect of
imperfect washout in the example below), but a reviewer complained that the
analysis is not correct, because Patient should be nested in Seq. As I can
see, the results look reasonable, but I admit that nesting is nowhere
indicated.
I know how to handle the case is Seq if assumed a random effect, but I don't
see why it should be random, AND I want an estimate of the effect.

So, if the question is stupid, please help me out at least an tell my WHY.


# cross-over design
library(nlme)
NSubj<-20
# Placebo Phase
# Seq: 1 for "Placebo first phase", 2 for "Verum first phase"
glup<-data.frame(Patient=as.factor(1:NSubj),
                Seq=as.factor(rep(c(1,2),NSubj/2)),
                Treat=rep("Pla",NSubj),
                Glu=rnorm(NSubj,30,4)
   )
# Verum Phase
gluv<-data.frame(Patient=as.factor(1:NSubj),
                Seq=as.factor(rep(c(1,2),NSubj/2)),
                Treat=rep("Ver",NSubj),
                Glu=glup$Glu-rnorm(NSubj,6,2)
   )

# Simulate imperfect washout when verum was in first phase
# == Interaction of Seq and Treat
glu<-rbind(gluv,glup)
glu$Glu<-glu$Glu+as.integer((glu$Seq=="2")&(glu$Treat=="Pla"))*rnorm(NSubj,2
,3)

# Sequence as fixed effect
glu.lme<-lme(Glu~Treat*Seq,data=glu,random=~1|Patient)
print(summary(glu.lme))

# result looks perfectly valid for me
#                  Value Std.Error DF   t-value p-value
#(Intercept)   24.694617  1.602160 18 15.413324  <.0001
#TreatPla       5.106889  1.096139 18  4.658980  0.0002
#Seq2           0.466712  2.265797 18  0.205981  0.8391
#TreatPla:Seq2  3.525737  1.550174 18  2.274413  0.0354

# Sequence as random ???????
glu.lme<-lme(Glu~Treat,data=glu,random=~1|Seq/Patient)
print(summary(glu.lme))



Dieter Menne


---------------------------------------
Dr. Dieter Menne
Biomed Software
72074 Tübingen
Tel (49) (7071) 52176
FAX (49) (7071) 55 10 46
dieter.menne at menne-biomed.de
www.menne-biomed.de



-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-
r-help mailing list -- Read http://www.ci.tuwien.ac.at/~hornik/R/R-FAQ.html
Send "info", "help", or "[un]subscribe"
(in the "body", not the subject !)  To: r-help-request at stat.math.ethz.ch
_._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._