Regresión logística
A continuación se mostrará un modelado de regresion logistica para el tema de los camarones, para poder hacer una estimación de la probabilidad de una variable binaria en función de una variable cuantitativa. Todo esto mediante un proceso que nos ayudará a determinar los resultados.
Caso de estudio: Camarones
Se analizara una empresa que produce camarón, en la cual a los camarones se le dan cierta cantidad de comida con el objetivo que lleguen a 12 gramos, lo cual se vería traducido en una buena cosecha.
Imagen representativa de camarones
Importar datos
library(pacman)
p_load("prettydoc", "DT", "xfun","readr")camarones <- read_csv("camarones.csv")## Rows: 12 Columns: 2## -- Column specification --------------------------------------------------------
## Delimiter: ","
## dbl (2): AlimentoDiario, Exito##
## 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.datatable(camarones)Datos recolectados de 12 estanques distintos
La columna “AlimentoDiario” representa los kilos de comida diarios durante la última semana de alimentación del estanque, en la columna “éxito” se pone 1 cuando los camarones tienen sus 12 gramos y la venta se puede realizar y con un 0 en el caso contrario, en esta tabla se obtuvó que 3 de los estanques tuvieron éxito.
Tabla de frecuencia de los datos
tabla <- table(camarones$Exito)
tabla##
## 0 1
## 9 3Ilustrando gráficamente
colores <- NULL
colores[camarones$Exito == 1] <- "green"
colores[camarones$Exito == 0] <- "red"
plot(camarones$AlimentoDiario, camarones$Exito, pch=21, bg=colores,
xlab="Cantidad de alimento diario", ylab="Probabilidades de éxito")
legend('bottomleft', c('Éxito', 'Sin éxito'), pch = 21, col = c('green', 'red'))
Haciendo incapie en estas observaciones, se determina que la cantidad de alimento dado al estanque es influyente. Pueden haber muchos factores externos que no permitan el crecimiento de los camarones.
Regresión logística
La variable AlimentoDiario será empleada para el modelo de regresión logistica y se determinará si el éxito es verdadero.
regresion <- glm(Exito ~ AlimentoDiario, data = camarones, family=binomial)
summary(regresion)##
## Call:
## glm(formula = Exito ~ AlimentoDiario, family = binomial, data = camarones)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -1.28965 -0.68424 -0.39705 -0.00008 2.00729
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -35.1229 25.8776 -1.357 0.175
## AlimentoDiario 0.1194 0.0901 1.325 0.185
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 13.496 on 11 degrees of freedom
## Residual deviance: 11.311 on 10 degrees of freedom
## AIC: 15.311
##
## Number of Fisher Scoring iterations: 5Predicción para valores nuevos con el modelo ajustado
Representando la funcion logistica se calcula la probabilidad de fallo estimadas. Para eso se creará un vector de alimentos entre 250 y 350 kilos.
datos_nuevos <- data.frame(AlimentoDiario = seq(250, 300, 1))
datos_nuevos## AlimentoDiario
## 1 250
## 2 251
## 3 252
## 4 253
## 5 254
## 6 255
## 7 256
## 8 257
## 9 258
## 10 259
## 11 260
## 12 261
## 13 262
## 14 263
## 15 264
## 16 265
## 17 266
## 18 267
## 19 268
## 20 269
## 21 270
## 22 271
## 23 272
## 24 273
## 25 274
## 26 275
## 27 276
## 28 277
## 29 278
## 30 279
## 31 280
## 32 281
## 33 282
## 34 283
## 35 284
## 36 285
## 37 286
## 38 287
## 39 288
## 40 289
## 41 290
## 42 291
## 43 292
## 44 293
## 45 294
## 46 295
## 47 296
## 48 297
## 49 298
## 50 299
## 51 300calculo de las nuevas probabilidades
probabilidades_nuevas <- predict(regresion, datos_nuevos, type="response")
probabilidades_nuevas## 1 2 3 4 5 6
## 0.005057235 0.005694741 0.006412092 0.007219149 0.008126955 0.009147866
## 7 8 9 10 11 12
## 0.010295693 0.011585859 0.013035568 0.014663985 0.016492425 0.018544562
## 13 14 15 16 17 18
## 0.020846632 0.023427652 0.026319639 0.029557817 0.033180821 0.037230872
## 19 20 21 22 23 24
## 0.041753922 0.046799751 0.052421995 0.058678083 0.065629069 0.073339313
## 25 26 27 28 29 30
## 0.081876000 0.091308444 0.101707168 0.113142702 0.125684103 0.139397160
## 31 32 33 34 35 36
## 0.154342290 0.170572147 0.188128986 0.207041842 0.227323661 0.248968480
## 37 38 39 40 41 42
## 0.271948870 0.296213806 0.321687180 0.348267152 0.375826503 0.404214109
## 43 44 45 46 47 48
## 0.433257569 0.462766946 0.492539476 0.522365017 0.552031926 0.581333034
## 49 50 51
## 0.610071349 0.638065168 0.665152328Por defecto esta probabilidad se calcula como: $ log p_{i}/(1-p_{i}) $
colores[camarones$Exito == 1] <- "green"
colores[camarones$Exito == 0] <- "red"
plot(camarones$AlimentoDiario, camarones$Exito, pch=21, bg=colores,
xlab="Cantidad de alimento diario", ylab="Probabilidades de éxito")
legend('bottomleft', c('Éxito', 'Sin éxito'), pch = 21, col = c('green', 'red'))
lines(datos_nuevos$AlimentoDiario, probabilidades_nuevas, col="blue", lwd=1)
Los datos pueden ser analizados en la siguiente gráfica