Analissidel accidente del transbordador espacial challenger

  • Antecedentes

poner imagen del challenger ## Regresión logistica

challenger <- read.table("http://verso.mat.uam.es/~joser.berrendero/datos/challenger.txt", header = TRUE)
table(challenger$defecto)
## 
##  0  1 
## 16  7
  • Representar los defectos vs no defestos de forma gráfica
colores <- NULL
colores[challenger$defecto==0] <- "green"
colores[challenger$defecto==1] <- "red"
plot(challenger$temp, challenger$defecto, pch=21, bg=colores, xlab = "Temperatura", ylab = "Probabilidad de defectos")
legend("bottomleft", c("No defecto", "Si defecto"), pch=21, col=c("green", "red"))

  • Para ajustar el modelo se usa el comando glm (para modelos lineales generalizados) indicando que la respuesta es binomial mediante el argumento family:
reg <- glm(defecto ~ temp, data=challenger, family = binomial)
summary(reg)
## 
## Call:
## glm(formula = defecto ~ temp, family = binomial, data = challenger)
## 
## Deviance Residuals: 
##     Min       1Q   Median       3Q      Max  
## -1.0611  -0.7613  -0.3783   0.4524   2.2175  
## 
## Coefficients:
##             Estimate Std. Error z value Pr(>|z|)  
## (Intercept)  15.0429     7.3786   2.039   0.0415 *
## temp         -0.2322     0.1082  -2.145   0.0320 *
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 28.267  on 22  degrees of freedom
## Residual deviance: 20.315  on 21  degrees of freedom
## AIC: 24.315
## 
## Number of Fisher Scoring iterations: 5

En el modelo de regresió logística la raíz de las desviaciones representa el papel de los residuos:

\[D_i = \mp \sqrt{=2 [Y_i\log \hat p_i + (1-Y_i)\log(1-\hat p_I)]},\] Donde el signo coincide con el signo de \(Y_i - \hat p_i\)

En la salida anterior estas cantidades se denominan devíance residuals. Para calcular estos pseudo-residuos, podemos ejecutar res = resid(reg).

res <- resid(reg)
res
##          1          2          3          4          5          6          7 
## -1.0611168  1.7145343 -0.7996042 -0.8817559 -0.9690847 -0.5865724 -0.5267644 
##          8          9         10         11         12         13         14 
## -0.7229433  0.5506685  1.0063470  1.7145343 -0.3018707 -0.9690847  0.3540506 
##         15         16         17         18         19         20         21 
## -0.9690847 -0.4229076 -0.7229433 -0.2143127 -0.3782680 -0.2694144  2.2175345 
##         22         23 
## -0.3782680  0.6127353
  • Pasaremos a calcular las probabilidades de fallo estimadas (usando el comando predict) para un vector adecuado de nuevas temperaturas (entre 50 y 85 grados)
datos <- data.frame(temp = seq(50,85,0.1))
probabilidades <- predict(reg, datos, type = "response")
#Gráfia
plot(challenger$temp, challenger$defecto, pch=21, bg=colores, xlab = "Temperatura", ylab = "Probabilidad de defectos")
legend("bottomleft", c("No defecto", "Si defecto"), pch=21, col=c("green", "red"))
lines(datos$temp, probabilidades, col="blue", lwd=3)