I maintain the Devore5 package for R. This package provides the data sets from Jay Devore's text "Probability and Statistics for Engineering and the Sciences (5th ed)". I am having difficulty reproducing some logistic regression results from the textbook. Perhaps this is because I am not using the glm function correctly. The data from Example 13.5 (page 559 for those with a copy of the book) are "the launch temperatures and the incidence of failure for O-rings in 24 space shuttle launches prior to the Challenger disaster of January 1986." Here are the data (Temperature in degrees Fahrenheit).> library(Devore5) > data(xmp13.05) > xmp13.05Temperature Failure 1 53 Y 2 56 Y 3 57 Y 4 63 N 5 66 N 6 67 N 7 67 N 8 67 N 9 68 N 10 69 N 11 70 N 12 70 Y 13 70 N 14 70 Y 15 72 N 16 73 N 17 75 N 18 75 Y 19 76 N 20 76 N 21 78 N 22 79 N 23 80 N 24 81 N> summary(xmp13.05)Temperature Failure Min. :53.00 N:18 1st Qu.:67.00 Y: 6 Median :70.00 Mean :69.92 3rd Qu.:75.25 Max. :81.00 Devore shows the parameter estimates for a logistic regression model fit by JMP as Term Estimate Std Error ChiSquare Prob>ChiSq Intercept 10.8753321 5.7031291 3.64 0.0565 temp -0.1713202 0.0834419 4.22 0.0401 I have not been able to reproduce these results in R. My results are > fm1 <- glm(Failure ~ Temperature, data = xmp13.05, family = "binomial") > summary(fm1) Call: glm(formula = Failure ~ Temperature, family = "binomial", data = xmp13.05) Deviance Residuals: Min 1Q Median 3Q Max -1.12489 -0.72335 -0.40139 -0.04056 2.22441 Coefficients: Estimate Std. Error z value Pr(>|z|) (Intercept) 11.74641 6.01916 1.952 0.0510 . Temperature -0.18843 0.08906 -2.116 0.0344 * --- Signif. codes: 0 `***' 0.001 `**' 0.01 `*' 0.05 `.' 0.1 ` ' 1 (Dispersion parameter for binomial family taken to be 1) Null deviance: 26.992 on 23 degrees of freedom Residual deviance: 20.482 on 22 degrees of freedom AIC: 24.482 Number of Fisher Scoring iterations: 4 Am I doing something wrong? -.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.- r-help mailing list -- Read http://www.ci.tuwien.ac.at/~hornik/R/R-FAQ.html Send "info", "help", or "[un]subscribe" (in the "body", not the subject !) To: r-help-request at stat.math.ethz.ch _._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._
Doug, It looks to me like a convergence issue. Here is the output from "another, not dissimilar" facility for data analysis and graphics (which I just happened to be using at the time): ---> dug <- read.table("dug.txt") > dm <- glm(Failure ~ Temperature, binomial, dug, trace=T)GLM linear loop 1: deviance = 21.387 GLM linear loop 2: deviance = 20.52 GLM linear loop 3: deviance = 20.482 GLM linear loop 4: deviance = 20.482> summary(dm)$coefValue Std. Error t value (Intercept) 11.74622 6.003206 1.9567 Temperature -0.18843 0.088811 -2.1217 --- If you look back at your R output it is very close to yours but with differences at about the fourth sig. fig. Now slacken off the convergence criterion to something pretty lax: ---> dm <- glm(Failure ~ Temperature, binomial, dug, trace=T, eps=0.1)GLM linear loop 1: deviance = 21.387 GLM linear loop 2: deviance = 20.52> summary(dm)$coefValue Std. Error t value (Intercept) 10.79758 4.880483 2.2124 Temperature -0.17356 0.071172 -2.4386 --- This is not altogether dissimilar to the result reported in Devore. Now couple that with different the starting values and I think you have a fair suggestion of what's going on. No, there's nothing wrong with the way you are using glm, but the important points I would make are 1. Always put trace = T in glm model fitting. I really don't know why it is not the default. 2. If in any doubt at all, see what difference it makes to strengthen the convergence criterion. This is surprisingly often necessary, especially in "another, not dissimilar" working environment where they don't have Brian Ripley working on their core team.... Regards, Bill. -- Bill Venables, Statistician, CMIS Environmetrics Project CSIRO Marine Labs, PO Box 120, Cleveland, Qld, AUSTRALIA. 4163 Tel: +61 7 3826 7251 Email: Bill.Venables at cmis.csiro.au Fax: +61 7 3826 7304 http://www.cmis.csiro.au/bill.venables/ -.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.- r-help mailing list -- Read http://www.ci.tuwien.ac.at/~hornik/R/R-FAQ.html Send "info", "help", or "[un]subscribe" (in the "body", not the subject !) To: r-help-request at stat.math.ethz.ch _._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._