Differences between R and SPSS in probit analysis

Previous Topic Next Topic
classic Classic list List threaded Threaded
1 message Options
Reply | Threaded
Open this post in threaded view
Report Content as Inappropriate

Differences between R and SPSS in probit analysis

Hi everybody,

I'm working on the effects of alternative larvicides on Aedes aegypti. Right now, I am doing binary mortality response with a single explanatory variable (dose) on 4 concentrations of one larvicide (+ control). Our university is fond of SPSS, and I have learned to conduct the basic probit model with it, including a natural logarithm transformation on my dosis data.
Not so long ago, I've started working with R, and through a combination of the 'glm' and 'dose.p' functions, I get the same slope and intercept, as well as LD50 calculations. Nevertheless, the standard errors and Z-scores calculated through the Probit model in SPSS comes out completely differently in R. Additionally, the 95% confidence intervals for the LD50 come out very differently between the two programs. I really don't have a clue on how I am getting the same slopes, intercepts and LD50's, but totally different SE, Z, and 95% CI. Can anybody help me so I can get the same results in R??

I'll pass you the script and hypothetical data:

dose <- c(6000, 4500, 3000, 1500, 0)
total <- c(100, 100, 100, 100, 100)
affected <- c(91, 82, 69, 49, 0)
finney71 <- data.frame(dose, total, affected)

fm1 <- glm(affected/total ~ log(dose),
 family=binomial(link = probit), data=finney71[finney71$dose != 0, ])

xp1 <- dose.p(fm1, p=c(0.5,0.9))
xp.ci <- xp1 + attr(xp1, "SE") %*% matrix(qnorm(1 - 0.05/2)*c(-1,1), nrow=1)
EAUS.Aa <- exp(cbind(xp1, attr(xp1, "SE"), xp.ci[,1], xp.ci[,2]))
dimnames(EAUS.Aa)[[2]] <- c("LD", "SE", "LCL","UCL")

So, this is the regression results I get with R:

Deviance Residuals:
       1         2         3         4  
 0.06655  -0.02814  -0.06268   0.03474  

            Estimate Std. Error z value
(Intercept)  -6.8940    10.7802  -0.640
log(dose)     0.9333     1.3441   0.694
(Intercept)    0.522
log(dose)      0.487

(Dispersion parameter for binomial family taken to be 1)

    Null deviance: 0.513878  on 3  degrees of freedom
Residual deviance: 0.010356  on 2  degrees of freedom
AIC: 6.5458

Number of Fisher Scoring iterations: 5

And the LD50 and CI transformed:

               LD       SE      LCL      UCL
p = 0.5: 1614.444 3.207876 164.3822 15855.91
p = 0.9: 6373.473 3.764879 474.1600 85669.72

These are the values I get on SPSS (just replacing the values on R output) :

            Estimate Std. Error z value
(Intercept)  -6.8940    1.082  -6.373
(dose)         2.149     0.311   6.918
And the LD50 and CI transformed:

                 LD          LCL      UCL
p = 0.5: 1614.444 1198.932 1953.120
p = 0.9: 6373.473 5145.767 9013.354

So, please if somebody can help me with this, I'd be grateful. If working with those functions won't do it, I'll use another, the one you recommend.

Thank you very much!