random effect variance per treatment group in lmer
Hmm, could it be a word-wrap issue? I just verified that the code works on R 2.5 with lme4 and lattice packages installed. See http://www.umaine.edu/visualperception/lme4het The model I posted assumes that the subj-to-subj variability in baseline score has higher variability in the "D" condition than in the "P" condition. You are correct that there's no random effect for drug. I don't think you *want* there to be a random effect for drug, since drug doesn't have levels sampled from some population (it's properly a fixed effect). Instead, the code I posted has the random effect of Patient with a higher variability in one level of drug than in the other level of drug. If you're looking some the drug to have a different time course for some patients than for other patients I think that'd be a random effect of Time. (To avoid overparameterization I think you'd almost surely want to treat time as a quantitative predictor if you're going to model it as having a random effect.) One of the reasons that I constructed a full example with a new sim was to get a big-enough data set together so that the lmer fit fairly obviously matched the parameters of the sim (to verify the correctness of the sim and the suggested model formula). Another reason is that the code of your sim didn't actually have different subj-to-subj variability in the "D" condition than it did in the "P" condition, so a model with such an effect was overparameterized for your original simulated data. Could this be why the lmer call I suggested failed for you? BTW, nice to see someone from the UM School of Business (I didn't notice until just now what your email address said). I worked there about 20 years ago on the Applied AI Reporter when I was an undergraduate psychology major. But, looking at the SBA web site, I recognize barely anyone, and those I do recognize probably wouldn't recognize me. alan
"Afshartous, David" <afshart at exchange.sba.miami.edu> on Friday, July 13, 2007 at 10:44 AM -0500 wrote:
Alan, Thanks for the suggestion. I noticed the following error msg for that lmer call:
( fm.het <- lmer( dv ~ rep(1:n.timepoints, n.subj.per.tx*2)*drug + ( 0
+ as.numeric(drug=="D") | Patient ) + ( 0 + as.numeric(drug=="P") | Patient ), data=dat.new ) ) Error: length(term) == 3 is not TRUE I tried a few changes to the model but the error still exists; I'll keep checking. I assume the rationale for the structure of your lmer call, where you use as.numeric as opposed to just drug above, is to insure that you do not introduce a random effect for drug into the model? Regards, Dave
-----Original Message-----
From: Alan Cobo-Lewis [mailto:alanc at umit.maine.edu]
Sent: Wednesday, July 11, 2007 6:40 PM
To: r-sig-mixed-models at r-project.org
Cc: " "Afshartous at basalt.its.maine.edu; Afshartous, David;
Andrew Robinson
Subject: Re: random effect variance per treatment group in lmer
Dave,
How about using stratifying variance on level of drug using (
0 + as.numeric(drug=="D") | Patient ) + ( 0 +
as.numeric(drug=="P") | Patient ) Here's some code (whose sim
also builds in a fixed effect of time that applies only to
the Drug condition).
set.seed(500)
n.timepoints <- 8
n.subj.per.tx <- 20
sd.d <- 5; sd.p <- 2; sd.res <- 1.3
drug <- factor(rep(c("D", "P"), each = n.timepoints, times =
n.subj.per.tx)) drug.baseline <- rep( c(0,5),
each=n.timepoints, times=n.subj.per.tx ) Patient <-
rep(1:(n.subj.per.tx*2), each = n.timepoints)
Patient.baseline <- rep( rnorm( n.subj.per.tx*2, sd=c(sd.d,
sd.p) ), each=n.timepoints ) time <- factor(paste("Time-",
rep(1:n.timepoints, n.subj.per.tx*2), sep="")) time.baseline
<- rep(1:n.timepoints,n.subj.per.tx*2)*as.numeric(drug=="D")
dv <- rnorm( n.subj.per.tx*n.timepoints*2,
mean=time.baseline+Patient.baseline+drug.baseline, sd=sd.res
) dat.new <- data.frame(time, drug, dv, Patient) xyplot(
dv~time|drug, group=Patient, type="l", data=dat.new ) # fit
model treats time as a quantitative predictor ( fm.het <-
lmer( dv ~ rep(1:n.timepoints, n.subj.per.tx*2)*drug + ( 0 +
as.numeric(drug=="D") | Patient ) + ( 0 +
as.numeric(drug=="P") | Patient ), data=dat.new ) )
HTH,
alan
From: "Afshartous, David" <afshart at exchange.sba.miami.edu>
asked:
All, I didn't receive a response to the query below sent to the general R-help mailing list so figured I'd try this mailing list.
Apologies in
advance if this is an overly simplistic question for this
list; I just
started w/ lmer after not using lme for awhile. Cheers, Dave
__________________________________________________________ All, How does one specify a model in lmer such that say the random effect for the intercept has a different variance per treatment group? Thus, in the model equation, we'd have say b_ij represent the random effect for patient j in treatment group i, with variance
depending on
i, i.e, var(b_ij) = tau_i. Didn't see this in the docs or Pinherio & Bates (section 5.2 is specific for modelling within group errors). Sample
repeated measures
code below is for a single random effect variance, where the random
effect corresponds to patient.
cheers,
dave
z <- rnorm(24, mean=0, sd=1)
time <- factor(paste("Time-", rep(1:6, 4), sep="")) Patient
<- rep(1:4,
each = 6) drug <- factor(rep(c("D", "P"), each = 6, times =
2)) ## P =
placebo, D = Drug dat.new <- data.frame(time, drug, z,
Patient) fm =
lmer(z ~ drug + time + (1 | Patient), data = dat.new )
-- Alan B. Cobo-Lewis, Ph.D. (207) 581-3840 tel Department of Psychology (207) 581-6128 fax University of Maine Orono, ME 04469-5742 alanc at maine.edu http://www.umaine.edu/visualperception
-- Alan B. Cobo-Lewis, Ph.D. (207) 581-3840 tel Department of Psychology (207) 581-6128 fax University of Maine Orono, ME 04469-5742 alanc at maine.edu http://www.umaine.edu/visualperception