[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