[R] inverse prediction and Poisson regression
Spencer Graves
spencer.graves at pdf.com
Fri Jul 25 12:29:57 CEST 2003
Dear Prof. Ripley & M. Philion:
First some commentary then questions for Prof. Ripley and M. Philion.
COMMENTARY
Prof. Ripley said, "to fit a curve of mean response vs dose,
and find the dose at which the mean response is half of that at dose 0.
That one is easy." Unfortunately, it is not obvious to me at the
moment. From "www.r-project.org" -> search -> "R site search" ->
"LD50", I found "dose.p", described on p. 193, sec. 7.2, of Venables and
Ripley (2002) Modern Applied Statistics with S, 4th ed. (Springer).
Then I cut the data set down to a size that I could easily play with,
and fit Poisson regression:
Phytopath <- data.frame(x=c(0, 0.03, 0.1),
y=c(28, 21, 11))
(fitP100 <- glm(y~log(x+0.015), data=Phytopath[rep(1:3, 100),],
family="poisson"))
Call: glm(formula = y ~ log(x + 0.015), family = "poisson", data =
Phytopath[rep(1:3, 100), ])
Coefficients:
(Intercept) log(x + 0.015)
1.6088 -0.4203
(LD50P100 <- dose.p(fitP100, p=14))
Dose SE
p = 14: -2.451018 0.04858572
To get a 95% confidence interval from this, in S-Plus 6.1, I did:
exp(LD50 + c(-2,0, 2)*LD50 at SE)-0.015
[1] 0.01762321 0.07120579 0.21279601
R 1.7.1 seemed to require a different syntax, which I couldn't parse
in my present insomniac state (3:20 AM in California).
QUESTIONS:
PROF. RIPLEY: Is this what you said was easy?
M. PHILION: Does this provide sufficient information for you to now
solve your problem?
hope this helps. spencer graves
Prof Brian Ripley wrote:
> On Fri, 25 Jul 2003, Vincent Philion wrote:
>
>
>>Hello and thank you for your interest in this problem.
>>
>>"real life data" would look like this:
>>
>>x y
>>0 28
>>0.03 21
>>0.1 11
>>0.3 15
>>1 5
>>3 4
>>10 1
>>30 0
>>100 0
>>
>>x y
>>0 30
>>0.0025 30
>>0.02 25
>>0.16 25
>>1.28 10
>>10.24 0
>>81.92 0
>>
>>X Y
>>0 35
>>0.00025 23
>>0.002 14
>>0.016 6
>>0.128 5
>>1.024 3
>>8.192 2
>>
>>X Y
>>0 43
>>0.00025 35
>>0.002 20
>>0.016 16
>>0.128 11
>>1.024 6
>>8.192 0
>>
>>Where X is dose and Y is response.
>>the relation is linear for log(response) = b log(dose) + intercept
>
>
> Is that log(*mean* response), that is a log link and exponential decay
> with dose?
>
>
>>Response for dose 0 is a "control" = Ymax. So, What I want is the dose
>>for 50% response. For instance, in example 1:
>>
>>Ymax = 28 (this is also an observation with Poisson error)
>
>
> Once you observe Ymax, Y is no longer Poisson.
>
>
>>So I want dose for response = 14 = approx. 0.3
>
>
> What exactly is Ymax? Is it the response at dose 0? The mean
response at
> dose 0? The largest response? About the only thing I can actually
> interpret is that you want to fit a curve of mean response vs dose, and
> find the dose at which the mean response is half of that at dose 0.
> That one is easy.
>
> I think you are confusing response with mean response, and we can't
> disentangle them for you.
>
More information about the R-help
mailing list