[R] LD50 contrasts with lmer/lme4

Dieter Menne dieter.menne at menne-biomed.de
Mon Feb 26 17:27:51 CET 2007


Dear R-list,

I have a data set from 20 pigs, each of which is tested at crossed 9 doses
(logdose -4:4) and 3 skin treatment substances when exposed to a standard
polluted environment. So there are 27 patches on each pig. The response is
irritation=yes/no.

I want to determine "equally effective 50% doses" (similar to old LD50), and
to test the treatments against each other. I am looking for something like
dose.p in MASS generalized to lmer (or glmmPQL or whatever). The direct as
output by lmer are not useful, because saying "30% irritation with A and 40%
with B at dose xx" has less meaning than giving "equivalent effective
doses".

Dieter

----- Simulated data -----
library(lme4)
animal = data.frame(ID = as.factor(1:20), da = rnorm(1:20))
treat = data.frame(treat=c('A','B','C'), treatoff=c(1,2,1.5),
       treatslope = c(0.5,0.6,0.7))

gr = expand.grid(animal=animal$ID,treat=treat$treat,logdose=c(-4:4))
gr$resp = as.integer(treat$treatoff[gr$treat]+
  treat$treatslope[gr$treat]*gr$logdose+
  animal$da[gr$animal] +  rnorm(nrow(gr),0,2) >0)

gr.lmer = lmer(resp ~ treat*logdose+(1|animal),data=gr,family=binomial)
summary(gr.lmer)

------- Output
Fixed effects:
               Estimate Std. Error z value Pr(>|z|)
(Intercept)      0.9553     0.3074    3.11   0.0019 **
treatB           0.8793     0.3313    2.65   0.0079 **
treatC           0.5516     0.3077    1.79   0.0730 .
logdose          0.3733     0.0774    4.82  1.4e-06 ***
treatB:logdose   0.3081     0.1323    2.33   0.0198 *
treatC:logdose   0.2666     0.1249    2.13   0.0328 *

----- Goal
                Value SD  p
50% logdose (A-B)  xx   xx  xx
50% logdose (A-C)  yy   yy  yy



More information about the R-help mailing list