[R-sig-ME] random effect variance per treatment group in lmer
Alan Cobo-Lewis
alanc at umit.maine.edu
Fri Jul 13 18:27:24 CEST 2007
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
More information about the R-sig-mixed-models
mailing list