Paquetes a utilizar
library(pacman)
p_load(rmdformats, readr, readxl, ggplot2, plotly, DT, xfun, gridExtra, leaflet, GGally, psych, corrplot, cluster)
Regresión Logística
En estadística, la regresión logística es un tipo de análisis de regresión utilizado para predecir el resultado de una variable categórica (una variable que puede adoptar un número limitado de categorías) en función de las variables independientes o predictoras. Es útil para modelar la probabilidad de un evento ocurriendo en función de otros factores. El análisis de regresión logística se enmarca en el conjunto de Modelos Lineales Generalizados (GLM por sus siglas en inglés) que usa como función de enlace la función logit. Las probabilidades que describen el posible resultado de un único ensayo se modelan como una función de variables explicativas, utilizando una función logística.
La regresión logística es usada extensamente en las ciencias médicas y sociales. Otros nombres para regresión logística usados en varias áreas de aplicación incluyen modelo logístico, modelo logit, y clasificador de máxima entropía.
Funcion sigmoide
IMAGEN
Caso de estudio: Accidente del challenger
El accidente del transbordador espacial Challenger se produjo el martes 28 de enero de 1986 a las 16:39:13 UTC, 1 cuando el transbordador espacial Challenger (misión STS-51-L) se desintegró 73 segundos tras el lanzamiento, 2 provocando la muerte de los siete miembros de la tripulación —Francis “Dick” Scobee, Michael J. Smith, Ronald McNair, Ellison Onizuka, Gregory Jarvis, Judith Resnik y Christa McAuliffe.3 La nave se desintegró sobre el océano Atlántico, frente a la costa del centro de Florida (Estados Unidos) a las 11:38 EST (16:38 UTC) Ha sido calificado como el accidente más grave en la conquista del espacio.
La desintegración del vehículo entero comenzó después de que una junta tórica de su cohete acelerador sólido (SRB) derecho fallara durante el despegue. El fallo de la junta tórica causó la apertura de una brecha, permitiendo que el gas caliente presurizado del interior del motor del cohete sólido saliera al exterior y contactara con la estructura adyacente de conexión con el SRB y el tanque externo de combustible. Esto provocó la separación de la conexión posterior del SRB derecho y el fallo estructural del depósito externo. Las fuerzas aerodinámicas destruyeron rápidamente el orbitador.
. También ignoraron las advertencias de los ingenieros sobre los peligros en el lanzamiento provocados por las frías temperaturas de aquella mañana y no habían informado adecuadamente a sus superiores de estas preocupaciones. La Comisión Rogers hizo nueve recomendaciones a la NASA que debía poner en práctica antes de continuar con los vuelos de transbordadores
IMAGEN
En 1986, el transbordador espacial Challenger tuvo un accidente catastrófico debido a un incendio en una de las piezas de sus propulsores. Era la vez 25 en que se lanzaba un transbordador espacial. En todas las ocasiones anteriores se habían inspeccionado los propulsores de las naves, y en algunas de ellas se habían encontrando defectos.
El fichero challenger contiene 23 observaciones de las siguientes variables: defecto, que toma los valores 1 y 0 en función de si se encontraron defectos o no en los propulsores; y temp, la temperatura (en grados Fahrenheit) en el momento del lanzamiento. 1. Cargar datos y contamos las frecuencias de casos sin y con defectos:
datos <- read_csv("challenger.csv")
## Rows: 23 Columns: 2
## -- Column specification --------------------------------------------------------
## Delimiter: ","
## dbl (2): temp, defecto
##
## i Use `spec()` to retrieve the full column specification for this data.
## i Specify the column types or set `show_col_types = FALSE` to quiet this message.
datos
## # A tibble: 23 x 2
## temp defecto
## <dbl> <dbl>
## 1 66 0
## 2 70 1
## 3 69 0
## 4 68 0
## 5 67 0
## 6 72 0
## 7 73 0
## 8 70 0
## 9 57 1
## 10 63 1
## # ... with 13 more rows
Observe que los datos están etiquetados en éxito= 1 y fracaso= 0 variable categórica “Y” donde si se encontraron defectos=1 o no=0 en los propulsores
Tabla de frecuencia
tabla <- table(datos$defecto)
tabla
##
## 0 1
## 16 7
Según la tabla anterior se presentaron 7 fallas de 23 ensayos o inspecciones y 16 de 23 no presentaron fallas en las turbinas
Representación gráfica
colores <- NULL
colores[datos$defecto == 0] <- "green"
colores[datos$defecto == 1] <- "red"
plot(datos$temp, datos$defecto, pch = 21, bg = colores,
xlab = "Temperatura", ylab = "Probabilidad de defectos")
legend("bottomleft", c("No defecto", "Si defecto"), pch = 21, col = c("green",
"red"))
Parece razonable, a la vista de los datos, pensar que la temperatura puede influir en la probabilidad de que los propulsores tengan defectos o No. Luego la idea es ajustar por medio de la regresión logistica,un modelo de regresión logística para estudiar la posible relación entre Y=1 dado diferentes valores de temperatura: P(Y=1|X). Para ajustar el modelo se usa el comando glm (para modelos lineales generalizados) indicando que la respuesta es binomial mediante el argumento family:
regresion <- glm(defecto ~ temp, data= datos, family = binomial)
summary(regresion)
##
## Call:
## glm(formula = defecto ~ temp, family = binomial, data = datos)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -1.00342 -0.58426 -0.18732 0.08337 2.00587
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 23.4033 11.8316 1.978 0.0479 *
## temp -0.3610 0.1755 -2.057 0.0397 *
## ---
## 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: 14.377 on 21 degrees of freedom
## AIC: 18.377
##
## Number of Fisher Scoring iterations: 7
La formula del modelo de regresion logistica es:
\[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 Yi−p̂i . En la salida anterior estas cantidades se denominan deviance residuals.
Analisis del modelo
En el modelo de regresión logística:
\[P(Y=1|X)=\dfrac{e^{15.0429-0.2322x}}{1+e^{15.0429-0.2322x}}\] En la tabla coeficientes, la temperatura es significativa con pvalor= 0.0320 el cual es menor que el valor de significancia de 0.05. Lo anterior nos permite concluir que la hipoesis nula se rechaza
\[H_0:\beta_{1}=0\] \[H_1:\beta_{1} \neq 0\]
Con respecto al intercepto: También es sinificativo dado que p valor= 0.0415< 0.05. Lo que significa que se rechaza la hipotesis nula, que afirma el intercepto es cero. Más aun, significa que bajísimos, lo que demuestra que ambas variables son importantes para explicar la variable dependiente (probabilidad de éxito).
Predicción para valores nuevos con el modelo ajustado
Notas importantes: Tenga en cuenta que si ha usado en la fórmula nombres de variables que implican el data.frame que las contiene, el nuevo data.frame debe llamarse igual y con los mismos nombres de variables. Es decir, si la fórmula es ay a x, para predecir debe usar el mismo nombre, a, para el data.frame, y por supuesto los mismos nombres de variables, aunque los valores sean distintos. Si usa sólo el nombre de las variables en la fórmula debe construir el data.frame nuevo con las columnas cuyos nombres representen a dichas variables.
##Ejemplo
Para representar gráficamente la función logística estimada, calculamos las probabilidades de fallo estimadas (usando el comando predict) para un vector adecuado de nuevas temperaturas (entre 50 y 85 grados): 1. Primero generamos una nueva secuencia de datos, el cula llamaremos datos nuevos
Para representar gráficamente la función logística estimada, calculamos las probabilidades de fallo estimadas (usando el comando predict) para un vector adecuado de nuevas temperaturas (entre 50 y 85 grados): 1. Primero generamos una nueva secuencia de datos, el cula llamaremos datos nuevos
datos_nuevos <- data.frame(temp = seq(50, 85, 0.1))
datos_nuevos
## temp
## 1 50.0
## 2 50.1
## 3 50.2
## 4 50.3
## 5 50.4
## 6 50.5
## 7 50.6
## 8 50.7
## 9 50.8
## 10 50.9
## 11 51.0
## 12 51.1
## 13 51.2
## 14 51.3
## 15 51.4
## 16 51.5
## 17 51.6
## 18 51.7
## 19 51.8
## 20 51.9
## 21 52.0
## 22 52.1
## 23 52.2
## 24 52.3
## 25 52.4
## 26 52.5
## 27 52.6
## 28 52.7
## 29 52.8
## 30 52.9
## 31 53.0
## 32 53.1
## 33 53.2
## 34 53.3
## 35 53.4
## 36 53.5
## 37 53.6
## 38 53.7
## 39 53.8
## 40 53.9
## 41 54.0
## 42 54.1
## 43 54.2
## 44 54.3
## 45 54.4
## 46 54.5
## 47 54.6
## 48 54.7
## 49 54.8
## 50 54.9
## 51 55.0
## 52 55.1
## 53 55.2
## 54 55.3
## 55 55.4
## 56 55.5
## 57 55.6
## 58 55.7
## 59 55.8
## 60 55.9
## 61 56.0
## 62 56.1
## 63 56.2
## 64 56.3
## 65 56.4
## 66 56.5
## 67 56.6
## 68 56.7
## 69 56.8
## 70 56.9
## 71 57.0
## 72 57.1
## 73 57.2
## 74 57.3
## 75 57.4
## 76 57.5
## 77 57.6
## 78 57.7
## 79 57.8
## 80 57.9
## 81 58.0
## 82 58.1
## 83 58.2
## 84 58.3
## 85 58.4
## 86 58.5
## 87 58.6
## 88 58.7
## 89 58.8
## 90 58.9
## 91 59.0
## 92 59.1
## 93 59.2
## 94 59.3
## 95 59.4
## 96 59.5
## 97 59.6
## 98 59.7
## 99 59.8
## 100 59.9
## 101 60.0
## 102 60.1
## 103 60.2
## 104 60.3
## 105 60.4
## 106 60.5
## 107 60.6
## 108 60.7
## 109 60.8
## 110 60.9
## 111 61.0
## 112 61.1
## 113 61.2
## 114 61.3
## 115 61.4
## 116 61.5
## 117 61.6
## 118 61.7
## 119 61.8
## 120 61.9
## 121 62.0
## 122 62.1
## 123 62.2
## 124 62.3
## 125 62.4
## 126 62.5
## 127 62.6
## 128 62.7
## 129 62.8
## 130 62.9
## 131 63.0
## 132 63.1
## 133 63.2
## 134 63.3
## 135 63.4
## 136 63.5
## 137 63.6
## 138 63.7
## 139 63.8
## 140 63.9
## 141 64.0
## 142 64.1
## 143 64.2
## 144 64.3
## 145 64.4
## 146 64.5
## 147 64.6
## 148 64.7
## 149 64.8
## 150 64.9
## 151 65.0
## 152 65.1
## 153 65.2
## 154 65.3
## 155 65.4
## 156 65.5
## 157 65.6
## 158 65.7
## 159 65.8
## 160 65.9
## 161 66.0
## 162 66.1
## 163 66.2
## 164 66.3
## 165 66.4
## 166 66.5
## 167 66.6
## 168 66.7
## 169 66.8
## 170 66.9
## 171 67.0
## 172 67.1
## 173 67.2
## 174 67.3
## 175 67.4
## 176 67.5
## 177 67.6
## 178 67.7
## 179 67.8
## 180 67.9
## 181 68.0
## 182 68.1
## 183 68.2
## 184 68.3
## 185 68.4
## 186 68.5
## 187 68.6
## 188 68.7
## 189 68.8
## 190 68.9
## 191 69.0
## 192 69.1
## 193 69.2
## 194 69.3
## 195 69.4
## 196 69.5
## 197 69.6
## 198 69.7
## 199 69.8
## 200 69.9
## 201 70.0
## 202 70.1
## 203 70.2
## 204 70.3
## 205 70.4
## 206 70.5
## 207 70.6
## 208 70.7
## 209 70.8
## 210 70.9
## 211 71.0
## 212 71.1
## 213 71.2
## 214 71.3
## 215 71.4
## 216 71.5
## 217 71.6
## 218 71.7
## 219 71.8
## 220 71.9
## 221 72.0
## 222 72.1
## 223 72.2
## 224 72.3
## 225 72.4
## 226 72.5
## 227 72.6
## 228 72.7
## 229 72.8
## 230 72.9
## 231 73.0
## 232 73.1
## 233 73.2
## 234 73.3
## 235 73.4
## 236 73.5
## 237 73.6
## 238 73.7
## 239 73.8
## 240 73.9
## 241 74.0
## 242 74.1
## 243 74.2
## 244 74.3
## 245 74.4
## 246 74.5
## 247 74.6
## 248 74.7
## 249 74.8
## 250 74.9
## 251 75.0
## 252 75.1
## 253 75.2
## 254 75.3
## 255 75.4
## 256 75.5
## 257 75.6
## 258 75.7
## 259 75.8
## 260 75.9
## 261 76.0
## 262 76.1
## 263 76.2
## 264 76.3
## 265 76.4
## 266 76.5
## 267 76.6
## 268 76.7
## 269 76.8
## 270 76.9
## 271 77.0
## 272 77.1
## 273 77.2
## 274 77.3
## 275 77.4
## 276 77.5
## 277 77.6
## 278 77.7
## 279 77.8
## 280 77.9
## 281 78.0
## 282 78.1
## 283 78.2
## 284 78.3
## 285 78.4
## 286 78.5
## 287 78.6
## 288 78.7
## 289 78.8
## 290 78.9
## 291 79.0
## 292 79.1
## 293 79.2
## 294 79.3
## 295 79.4
## 296 79.5
## 297 79.6
## 298 79.7
## 299 79.8
## 300 79.9
## 301 80.0
## 302 80.1
## 303 80.2
## 304 80.3
## 305 80.4
## 306 80.5
## 307 80.6
## 308 80.7
## 309 80.8
## 310 80.9
## 311 81.0
## 312 81.1
## 313 81.2
## 314 81.3
## 315 81.4
## 316 81.5
## 317 81.6
## 318 81.7
## 319 81.8
## 320 81.9
## 321 82.0
## 322 82.1
## 323 82.2
## 324 82.3
## 325 82.4
## 326 82.5
## 327 82.6
## 328 82.7
## 329 82.8
## 330 82.9
## 331 83.0
## 332 83.1
## 333 83.2
## 334 83.3
## 335 83.4
## 336 83.5
## 337 83.6
## 338 83.7
## 339 83.8
## 340 83.9
## 341 84.0
## 342 84.1
## 343 84.2
## 344 84.3
## 345 84.4
## 346 84.5
## 347 84.6
## 348 84.7
## 349 84.8
## 350 84.9
## 351 85.0
Calculo de las nuevas probabilidades
probabilidades_nuevas <- predict(regresion, datos_nuevos, type = "response")
probabilidades_nuevas
## 1 2 3 4 5 6
## 0.9952848912 0.9951124030 0.9949336370 0.9947483670 0.9945563591 0.9943573708
## 7 8 9 10 11 12
## 0.9941511515 0.9939374415 0.9937159721 0.9934864653 0.9932486334 0.9930021785
## 13 14 15 16 17 18
## 0.9927467927 0.9924821571 0.9922079419 0.9919238061 0.9916293968 0.9913243489
## 19 20 21 22 23 24
## 0.9910082850 0.9906808149 0.9903415348 0.9899900274 0.9896258614 0.9892485908
## 25 26 27 28 29 30
## 0.9888577547 0.9884528767 0.9880334646 0.9875990099 0.9871489871 0.9866828535
## 31 32 33 34 35 36
## 0.9862000486 0.9856999936 0.9851820908 0.9846457234 0.9840902545 0.9835150269
## 37 38 39 40 41 42
## 0.9829193625 0.9823025619 0.9816639034 0.9810026430 0.9803180134 0.9796092239
## 43 44 45 46 47 48
## 0.9788754595 0.9781158803 0.9773296213 0.9765157917 0.9756734740 0.9748017242
## 49 50 51 52 53 54
## 0.9738995707 0.9729660138 0.9720000255 0.9710005490 0.9699664980 0.9688967562
## 55 56 57 58 59 60
## 0.9677901774 0.9666455844 0.9654617693 0.9642374926 0.9629714834 0.9616624388
## 61 62 63 64 65 66
## 0.9603090239 0.9589098714 0.9574635819 0.9559687235 0.9544238319 0.9528274105
## 67 68 69 70 71 72
## 0.9511779309 0.9494738326 0.9477135236 0.9458953809 0.9440177511 0.9420789507
## 73 74 75 76 77 78
## 0.9400772670 0.9380109593 0.9358782594 0.9336773731 0.9314064813 0.9290637415
## 79 80 81 82 83 84
## 0.9266472893 0.9241552403 0.9215856921 0.9189367260 0.9162064099 0.9133928002
## 85 86 87 88 89 90
## 0.9104939449 0.9075078865 0.9044326649 0.9012663210 0.8980068999 0.8946524554
## 91 92 93 94 95 96
## 0.8912010531 0.8876507756 0.8839997264 0.8802460351 0.8763878618 0.8724234031
## 97 98 99 100 101 102
## 0.8683508970 0.8641686289 0.8598749373 0.8554682202 0.8509469413 0.8463096364
## 103 104 105 106 107 108
## 0.8415549205 0.8366814943 0.8316881515 0.8265737857 0.8213373978 0.8159781032
## 109 110 111 112 113 114
## 0.8104951392 0.8048878726 0.7991558064 0.7932985878 0.7873160149 0.7812080439
## 115 116 117 118 119 120
## 0.7749747960 0.7686165639 0.7621338184 0.7555272143 0.7487975960 0.7419460029
## 121 122 123 124 125 126
## 0.7349736743 0.7278820537 0.7206727921 0.7133477520 0.7059090091 0.6983588546
## 127 128 129 130 131 132
## 0.6906997963 0.6829345583 0.6750660816 0.6670975216 0.6590322469 0.6508738356
## 133 134 135 136 137 138
## 0.6426260719 0.6342929408 0.6258786226 0.6173874862 0.6088240816 0.6001931315
## 139 140 141 142 143 144
## 0.5914995222 0.5827482932 0.5739446266 0.5650938355 0.5562013519 0.5472727134
## 145 146 147 148 149 150
## 0.5383135500 0.5293295705 0.5203265473 0.5113103024 0.5022866922 0.4932615922
## 151 152 153 154 155 156
## 0.4842408817 0.4752304287 0.4662360744 0.4572636184 0.4483188034 0.4394073012
## 157 158 159 160 161 162
## 0.4305346978 0.4217064805 0.4129280242 0.4042045787 0.3955412574 0.3869430254
## 163 164 165 166 167 168
## 0.3784146896 0.3699608884 0.3615860837 0.3532945522 0.3450903786 0.3369774495
## 169 170 171 172 173 174
## 0.3289594480 0.3210398492 0.3132219170 0.3055087012 0.2979030359 0.2904075386
## 175 176 177 178 179 180
## 0.2830246100 0.2757564350 0.2686049838 0.2615720141 0.2546590743 0.2478675066
## 181 182 183 184 185 186
## 0.2411984510 0.2346528502 0.2282314544 0.2219348269 0.2157633498 0.2097172304
## 187 188 189 190 191 192
## 0.2037965077 0.1980010590 0.1923306072 0.1867847273 0.1813628543 0.1760642901
## 193 194 195 196 197 198
## 0.1708882111 0.1658336752 0.1608996294 0.1560849173 0.1513882858 0.1468083923
## 199 200 201 202 203 204
## 0.1423438120 0.1379930444 0.1337545198 0.1296266061 0.1256076150 0.1216958077
## 205 206 207 208 209 210
## 0.1178894010 0.1141865729 0.1105854679 0.1070842019 0.1036808677 0.1003735389
## 211 212 213 214 215 216
## 0.0971602748 0.0940391243 0.0910081301 0.0880653321 0.0852087711 0.0824364920
## 217 218 219 220 221 222
## 0.0797465464 0.0771369962 0.0746059157 0.0721513937 0.0697715366 0.0674644697
## 223 224 225 226 227 228
## 0.0652283390 0.0630613135 0.0609615862 0.0589273755 0.0569569267 0.0550485127
## 229 230 231 232 233 234
## 0.0532004352 0.0514110255 0.0496786450 0.0480016860 0.0463785720 0.0448077583
## 235 236 237 238 239 240
## 0.0432877320 0.0418170126 0.0403941518 0.0390177336 0.0376863746 0.0363987233
## 241 242 243 244 245 246
## 0.0351534608 0.0339492998 0.0327849851 0.0316592928 0.0305710301 0.0295190352
## 247 248 249 250 251 252
## 0.0285021767 0.0275193533 0.0265694930 0.0256515534 0.0247645206 0.0239074088
## 253 254 255 256 257 258
## 0.0230792600 0.0222791434 0.0215061548 0.0207594160 0.0200380746 0.0193413030
## 259 260 261 262 263 264
## 0.0186682984 0.0180182814 0.0173904965 0.0167842108 0.0161987136 0.0156333160
## 265 266 267 268 269 270
## 0.0150873503 0.0145601694 0.0140511465 0.0135596741 0.0130851640 0.0126270465
## 271 272 273 274 275 276
## 0.0121847698 0.0117577999 0.0113456198 0.0109477289 0.0105636429 0.0101928932
## 277 278 279 280 281 282
## 0.0098350263 0.0094896034 0.0091562002 0.0088344061 0.0085238242 0.0082240705
## 283 284 285 286 287 288
## 0.0079347737 0.0076555749 0.0073861271 0.0071260948 0.0068751536 0.0066329901
## 289 290 291 292 293 294
## 0.0063993013 0.0061737946 0.0059561869 0.0057462049 0.0055435844 0.0053480702
## 295 296 297 298 299 300
## 0.0051594157 0.0049773828 0.0048017412 0.0046322689 0.0044687510 0.0043109803
## 301 302 303 304 305 306
## 0.0041587565 0.0040118861 0.0038701824 0.0037334651 0.0036015600 0.0034742989
## 307 308 309 310 311 312
## 0.0033515194 0.0032330649 0.0031187838 0.0030085301 0.0029021626 0.0027995453
## 313 314 315 316 317 318
## 0.0027005466 0.0026050395 0.0025129017 0.0024240147 0.0023382645 0.0022555409
## 319 320 321 322 323 324
## 0.0021757375 0.0020987517 0.0020244845 0.0019528401 0.0018837264 0.0018170542
## 325 326 327 328 329 330
## 0.0017527377 0.0016906939 0.0016308427 0.0015731069 0.0015174120 0.0014636861
## 331 332 333 334 335 336
## 0.0014118597 0.0013618659 0.0013136401 0.0012671198 0.0012222449 0.0011789574
## 337 338 339 340 341 342
## 0.0011372013 0.0010969224 0.0010580687 0.0010205898 0.0009844371 0.0009495639
## 343 344 345 346 347 348
## 0.0009159250 0.0008834766 0.0008521768 0.0008219850 0.0007928620 0.0007647701
## 349 350 351
## 0.0007376727 0.0007115348 0.0006863224
Por defecto esta probabilidad se calcula como:
\[log p_{i}/(1-p_{i})\]
Representación grafica del ajuste
plot(datos$temp, datos$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_nuevos$temp, probabilidades_nuevas, col = "blue", lwd = 2)