R help - Aug 2006 - Finney's fiducial confidence intervals of LD50

I am working with Probit regression (I cannot switch to logit) can anybody help
me in finding out how to obtain with R Finney's fiducial confidence
intervals for the levels of the predictor (Dose) needed to produce a proportion
of 50% of responses(LD50, ED50 etc.)?
  If the Pearson chi-square goodness-of-fit test is significant (by default), a
heterogeneity factor should be used to calculate the limits.
   
  Response<-c(0,7,26,27,0,5,13,29,0,4,11,25)
  Tot<-rep(30.5,12)
  Dose<-rep(c(10,40,160,640),3)
  probit<-glm(formula = Response/Tot~ log10(Dose), family=quasibinomial
  (link=probit))
  D50<- round(10^(dose.p(probit,cf=1:2,p=0.5)))
   
  #This is what SPSS calculates. I would like to reproduce these results with R:
  #SPSS RESULTS:
  #PRNT50= 140,83525
  #CI = 98,37857;205,34483
  #Regr.coeff= 1,91676 (S.E.=0,16765)
  #Intercept=-4,11856 (S.E.=0,36355) 
  Thanks a lot for your help.
   
  Carlos

 __________________________________________________



	[[alternative HTML version deleted]]

I don't know what Finney's fiducial confidence interval is but if your
problem is to handle the output of dose.p (from MASS), you can do as
follows:
> library(MASS)
> Response <- c(0, 7, 26, 27, 0, 5, 13, 29, 0, 4, 11, 25)
> Tot <- rep(30.5, 12)
> Dose <- rep(c(10, 40, 160, 640), 3)
> fm <- glm(Response/Tot ~ log10(Dose), family = quasibinomial(link =
probit))
> logD50 <- dose.p(fm, cf = 1:2, p = 0.5)
> D50 <- 10^c(logD50 + c(1, -1.96, 1.96) * attr(logD50, "SE"))
> names(D50) <- c("D50", "lower", "upper")
> D50
     D50    lower    upper
164.9506 103.3171 191.9777

Best,

Renaud

2006/8/21, carlos riveira <carlos.riveira at
yahoo.com>:
>   I am working with Probit regression (I cannot switch to logit) can
anybody help me in finding out how to obtain with R Finney's fiducial
confidence intervals for the levels of the predictor (Dose) needed to produce a
proportion of 50% of responses(LD50, ED50 etc.)?
>   If the Pearson chi-square goodness-of-fit test is significant (by
default), a heterogeneity factor should be used to calculate the limits.
>
>   Response<-c(0,7,26,27,0,5,13,29,0,4,11,25)
>   Tot<-rep(30.5,12)
>   Dose<-rep(c(10,40,160,640),3)
>   probit<-glm(formula = Response/Tot~ log10(Dose), family=quasibinomial
>   (link=probit))
>   D50<- round(10^(dose.p(probit,cf=1:2,p=0.5)))
>
>   #This is what SPSS calculates. I would like to reproduce these results
with R:
>   #SPSS RESULTS:
>   #PRNT50= 140,83525
>   #CI = 98,37857;205,34483
>   #Regr.coeff= 1,91676 (S.E.=0,16765)
>   #Intercept=-4,11856 (S.E.=0,36355)
>   Thanks a lot for your help.
>
>   Carlos
>
>  __________________________________________________
>
>
>
>         [[alternative HTML version deleted]]
>
> ______________________________________________
> R-help at stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide
http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>

-- 
Renaud LANCELOT
D?partement Elevage et M?decine V?t?rinaire (EMVT) du CIRAD
Directeur adjoint charg? des affaires scientifiques

CIRAD, Animal Production and Veterinary Medicine Department
Deputy director for scientific affairs

Campus international de Baillarguet
TA 30 / B (B?t. B, Bur. 214)
34398 Montpellier Cedex 5 - France
T?l   +33 (0)4 67 59 37 17
Secr. +33 (0)4 67 59 39 04
Fax   +33 (0)4 67 59 37 95

carlos riveira <carlos.riveira <at> yahoo.com> writes:
> 
> I am working with Probit regression (I cannot switch to logit) can anybody 
> help me in finding out how to  obtain with R Finney's fiducial
confidence
> intervals for the levels of the predictor (Dose) needed to
> produce a proportion of 50% of responses(LD50, ED50 etc.)? 
See the example in MASS on <dose.p>. Don't know if it is comes out the
same as
SPSS, though.

ldose <- rep(0:5, 2)
numdead <- c(1, 4, 9, 13, 18, 20, 0, 2, 6, 10, 12, 16)
sex <- factor(rep(c("M", "F"), c(6, 6)))
SF <- cbind(numdead, numalive = 20 - numdead)
budworm.lg0 <- glm(SF ~ sex + ldose - 1, family = binomial)

dose.p(budworm.lg0, cf = c(1,3), p = 1:3/4)
dose.p(update(budworm.lg0, family = binomial(link=probit)),
       cf = c(1,3), p = 1:3/4)

Dieter