Regresión logística
En el presente documento se realizara un modelado de regresion logistica para el tema de los camarones, para asi estimar la probabilidad de una variable binaria en función de una variable cuantitativa. Todo esto con un proceso de metodos que nos ayudaran a determinar los resultados.
¿Por qué regresión logística y no lineal?
Si una variable cualitativa con dos niveles se codifica como 1 y 0, matemáticamente es posible ajustar un modelo de regresión lineal por mínimos cuadrados y= β0+β1x . El problema de esta aproximación es que, al tratarse de una recta, para valores extremos del predictor, se obtienen valores de Y menores que 0 o mayores que 1, lo que entra en contradicción con el hecho de que las probabilidades siempre están dentro del rango [0,1].
Caso de estudio: Camarones
En este caso de estudio se analiza una empresa productora de camarones, en el cual se brinda cierta cantidad de comida con el objetivo de llegar a 12 gramos, esto se representa como una buena cosecha.
Imagen ilustrativa para el caso de estudio 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)Aquí se muestran datos de 12 estanques.
La columna “AlimentoDiario” se utiliza para representar los kilos de comida por día durante la última semana de alimentación a un determinado estanque, la columna “éxito” es utilizada para poner un 1 cuando los camarones tienen sus 12 gramos y la venta se puede realizar, y con un 0 cuando no es posible realizar la venta, solamente 3 estanques tuvieron éxito.
Tabla de frecuencia de los datos
tabla <- table(camarones$Exito)
tabla##
## 0 1
## 9 3
Mediante esta table es posible apreciar los 3 estanques de 12 en total posibles llegaron a ser considerados como exitosos, y de esta manera poder ser vendidos.
Ilustrando 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'))Luego de determinar los datos, es posible observar que la cantidad de alimento tiene influencia en el caso de estudio, y en cambio en otros casos no es así. Existen muchos factores causantes de la falta de crecimiento de los camarones, debido a que es una especie frágil.
Regresión logística
La variable AlimentoDiario se utiliza en el modelo de regresión logística, con esto se establece la hipótesis de 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: 5
Predicción para valores nuevos con el modelo ajustado
Al representar la función logística se calculan estimaciones de la probabilidad de fallo. También se emplea 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 300
Cálculo 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.665152328
Por defecto esta probabilidad se calcula de esta manera: $ 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)A través de la gráfica es posible analizar los datos.