Regresión logística
library(pacman)
p_load("DT","prettydoc","xfun")Regresión logística
La Regresión Logística Simple, desarrollada por David Cox en 1958, es un método de regresión que permite estimar la probabilidad de una variable cualitativa binaria en función de una variable cuantitativa. Una de las principales aplicaciones de la regresión logística es la de clasificación binaria, en el que las observaciones se clasifican en un grupo u otro dependiendo del valor que tome la variable empleada como predictor. Por ejemplo, clasificar a un individuo desconocido como hombre o mujer en función del tamaño de la mandíbula.
Es importante tener en cuenta que, aunque la regresión logística permite clasificar, se trata de un modelo de regresión que modela el logaritmo de la probabilidad de pertenecer a cada grupo. La asignación final se hace en función de las probabilidades predichas.
La existencia de una relación significativa entre una variable cualitativa con dos niveles y una variable continua se puede estudiar mediante otros test estadísticos tales como t-test o ANOVA (un ANOVA de dos grupos es equivalente al t-test). Sin embargo, la regresión logística permite además calcular la probabilidad de que la variable dependiente pertenezca a cada una de las dos categorías en función del valor que adquiera la variable independiente.
¿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 β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].
En el siguiente ejemplo se modela la probabilidad de fraude por impago (default) en función del balance de la cuenta bancaria (balance).
Caso de estudio: Estanques de camarones
En este caso de estudio se hará un analisis acerca de los datos sobre un estanque de camarones. Los datos se plasmaran en una tabla tales como el número del estanque, la cantidad de comida en kilos, y si el camarón llego a su peso ideal y fue vendido se representará con el número 1, en caso contrario con el número 0.
Estanque de camarones
Importar datos
datos <- read.csv("camarones.csv")
datatable(datos)Observe que los datos estan etiquetados en exito=1 y fracaso=0 dentro de la columna de éxito.
Tabla de frecuencia de los datos
tabla <- table(datos$Exito)
tabla##
## 0 1
## 9 3
Según la tabla anterior se presentarón 7 fallas de 23 ensayos o inspeciones y 16 de 23 no presentaron fallas en las turbinas
Viendo esta relación gráficamente
colores <- NULL
colores[datos$Exito == 0] <- "green"
colores[datos$Exito == 1] <- "red"
plot(datos$ï..AlimentoDiario, datos$Exito, pch=21, bg= colores,
xlab = "Cantidad de alimento diario", ylab = "Camarones Aprobados")
legend("topleft", c("Fracaso", "Aprobado"), pch=21, col = c("green", "red") )regresion <- glm(Exito ~ ï..AlimentoDiario , data=datos, family= binomial)
summary(regresion)##
## Call:
## glm(formula = Exito ~ ï..AlimentoDiario, family = binomial, data = datos)
##
## 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
en el modelo de regresión 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 Yi−p̂i . En la salida anterior estas cantidades se denominan deviance residuals.
Análisis del modelo
Formulación matemática del modelo de regresión logísitca
\[ P(Y=1|X)=\dfrac{e^{15.0429-0.2322x}}{1+e^{15.0429-0.2322x}} \] Existen varias funciones que cumplen esta descripción, una de las más utilizadas es la función logística (también conocida como función sigmoide): \[ \text{función sigmoide} = \sigma(x) = \dfrac{1}{1 + e^{-x}} \tag{1} \]
Predicción para valores nuevos con el modelo ajustado
Para representar gráficamente la función logística estimada, calculamos las probabilidades de éxtio estimadas para obtener un rango en cuanto a las cantidades de comida, aproximadamente de 250 a 280kg.
datos_nuevos <- data.frame(ï..AlimentoDiario = seq (250, 280, 0.1))
datos_nuevos ## ï..AlimentoDiario
## 1 250.0
## 2 250.1
## 3 250.2
## 4 250.3
## 5 250.4
## 6 250.5
## 7 250.6
## 8 250.7
## 9 250.8
## 10 250.9
## 11 251.0
## 12 251.1
## 13 251.2
## 14 251.3
## 15 251.4
## 16 251.5
## 17 251.6
## 18 251.7
## 19 251.8
## 20 251.9
## 21 252.0
## 22 252.1
## 23 252.2
## 24 252.3
## 25 252.4
## 26 252.5
## 27 252.6
## 28 252.7
## 29 252.8
## 30 252.9
## 31 253.0
## 32 253.1
## 33 253.2
## 34 253.3
## 35 253.4
## 36 253.5
## 37 253.6
## 38 253.7
## 39 253.8
## 40 253.9
## 41 254.0
## 42 254.1
## 43 254.2
## 44 254.3
## 45 254.4
## 46 254.5
## 47 254.6
## 48 254.7
## 49 254.8
## 50 254.9
## 51 255.0
## 52 255.1
## 53 255.2
## 54 255.3
## 55 255.4
## 56 255.5
## 57 255.6
## 58 255.7
## 59 255.8
## 60 255.9
## 61 256.0
## 62 256.1
## 63 256.2
## 64 256.3
## 65 256.4
## 66 256.5
## 67 256.6
## 68 256.7
## 69 256.8
## 70 256.9
## 71 257.0
## 72 257.1
## 73 257.2
## 74 257.3
## 75 257.4
## 76 257.5
## 77 257.6
## 78 257.7
## 79 257.8
## 80 257.9
## 81 258.0
## 82 258.1
## 83 258.2
## 84 258.3
## 85 258.4
## 86 258.5
## 87 258.6
## 88 258.7
## 89 258.8
## 90 258.9
## 91 259.0
## 92 259.1
## 93 259.2
## 94 259.3
## 95 259.4
## 96 259.5
## 97 259.6
## 98 259.7
## 99 259.8
## 100 259.9
## 101 260.0
## 102 260.1
## 103 260.2
## 104 260.3
## 105 260.4
## 106 260.5
## 107 260.6
## 108 260.7
## 109 260.8
## 110 260.9
## 111 261.0
## 112 261.1
## 113 261.2
## 114 261.3
## 115 261.4
## 116 261.5
## 117 261.6
## 118 261.7
## 119 261.8
## 120 261.9
## 121 262.0
## 122 262.1
## 123 262.2
## 124 262.3
## 125 262.4
## 126 262.5
## 127 262.6
## 128 262.7
## 129 262.8
## 130 262.9
## 131 263.0
## 132 263.1
## 133 263.2
## 134 263.3
## 135 263.4
## 136 263.5
## 137 263.6
## 138 263.7
## 139 263.8
## 140 263.9
## 141 264.0
## 142 264.1
## 143 264.2
## 144 264.3
## 145 264.4
## 146 264.5
## 147 264.6
## 148 264.7
## 149 264.8
## 150 264.9
## 151 265.0
## 152 265.1
## 153 265.2
## 154 265.3
## 155 265.4
## 156 265.5
## 157 265.6
## 158 265.7
## 159 265.8
## 160 265.9
## 161 266.0
## 162 266.1
## 163 266.2
## 164 266.3
## 165 266.4
## 166 266.5
## 167 266.6
## 168 266.7
## 169 266.8
## 170 266.9
## 171 267.0
## 172 267.1
## 173 267.2
## 174 267.3
## 175 267.4
## 176 267.5
## 177 267.6
## 178 267.7
## 179 267.8
## 180 267.9
## 181 268.0
## 182 268.1
## 183 268.2
## 184 268.3
## 185 268.4
## 186 268.5
## 187 268.6
## 188 268.7
## 189 268.8
## 190 268.9
## 191 269.0
## 192 269.1
## 193 269.2
## 194 269.3
## 195 269.4
## 196 269.5
## 197 269.6
## 198 269.7
## 199 269.8
## 200 269.9
## 201 270.0
## 202 270.1
## 203 270.2
## 204 270.3
## 205 270.4
## 206 270.5
## 207 270.6
## 208 270.7
## 209 270.8
## 210 270.9
## 211 271.0
## 212 271.1
## 213 271.2
## 214 271.3
## 215 271.4
## 216 271.5
## 217 271.6
## 218 271.7
## 219 271.8
## 220 271.9
## 221 272.0
## 222 272.1
## 223 272.2
## 224 272.3
## 225 272.4
## 226 272.5
## 227 272.6
## 228 272.7
## 229 272.8
## 230 272.9
## 231 273.0
## 232 273.1
## 233 273.2
## 234 273.3
## 235 273.4
## 236 273.5
## 237 273.6
## 238 273.7
## 239 273.8
## 240 273.9
## 241 274.0
## 242 274.1
## 243 274.2
## 244 274.3
## 245 274.4
## 246 274.5
## 247 274.6
## 248 274.7
## 249 274.8
## 250 274.9
## 251 275.0
## 252 275.1
## 253 275.2
## 254 275.3
## 255 275.4
## 256 275.5
## 257 275.6
## 258 275.7
## 259 275.8
## 260 275.9
## 261 276.0
## 262 276.1
## 263 276.2
## 264 276.3
## 265 276.4
## 266 276.5
## 267 276.6
## 268 276.7
## 269 276.8
## 270 276.9
## 271 277.0
## 272 277.1
## 273 277.2
## 274 277.3
## 275 277.4
## 276 277.5
## 277 277.6
## 278 277.7
## 279 277.8
## 280 277.9
## 281 278.0
## 282 278.1
## 283 278.2
## 284 278.3
## 285 278.4
## 286 278.5
## 287 278.6
## 288 278.7
## 289 278.8
## 290 278.9
## 291 279.0
## 292 279.1
## 293 279.2
## 294 279.3
## 295 279.4
## 296 279.5
## 297 279.6
## 298 279.7
## 299 279.8
## 300 279.9
## 301 280.0
Cálculo de nuevas probabilidades
probabilidades <- predict(regresion, datos_nuevos, type="response")
probabilidades## 1 2 3 4 5 6
## 0.005057235 0.005117652 0.005178786 0.005240646 0.005303242 0.005366581
## 7 8 9 10 11 12
## 0.005430673 0.005495526 0.005561148 0.005627551 0.005694741 0.005762729
## 13 14 15 16 17 18
## 0.005831524 0.005901135 0.005971573 0.006042846 0.006114964 0.006187938
## 19 20 21 22 23 24
## 0.006261777 0.006336491 0.006412092 0.006488588 0.006565991 0.006644311
## 25 26 27 28 29 30
## 0.006723559 0.006803746 0.006884882 0.006966979 0.007050049 0.007134101
## 31 32 33 34 35 36
## 0.007219149 0.007305202 0.007392274 0.007480376 0.007569520 0.007659718
## 37 38 39 40 41 42
## 0.007750983 0.007843326 0.007936761 0.008031299 0.008126955 0.008223740
## 43 44 45 46 47 48
## 0.008321669 0.008420753 0.008521008 0.008622445 0.008725080 0.008828925
## 49 50 51 52 53 54
## 0.008933995 0.009040304 0.009147866 0.009256696 0.009366809 0.009478219
## 55 56 57 58 59 60
## 0.009590941 0.009704991 0.009820383 0.009937134 0.010055258 0.010174773
## 61 62 63 64 65 66
## 0.010295693 0.010418035 0.010541815 0.010667050 0.010793757 0.010921952
## 67 68 69 70 71 72
## 0.011051653 0.011182876 0.011315640 0.011449962 0.011585859 0.011723350
## 73 74 75 76 77 78
## 0.011862454 0.012003188 0.012145570 0.012289621 0.012435359 0.012582803
## 79 80 81 82 83 84
## 0.012731973 0.012882888 0.013035568 0.013190034 0.013346305 0.013504403
## 85 86 87 88 89 90
## 0.013664347 0.013826159 0.013989860 0.014155471 0.014323014 0.014492512
## 91 92 93 94 95 96
## 0.014663985 0.014837456 0.015012949 0.015190485 0.015370088 0.015551780
## 97 98 99 100 101 102
## 0.015735586 0.015921530 0.016109635 0.016299925 0.016492425 0.016687161
## 103 104 105 106 107 108
## 0.016884156 0.017083436 0.017285027 0.017488955 0.017695245 0.017903924
## 109 110 111 112 113 114
## 0.018115018 0.018328555 0.018544562 0.018763066 0.018984094 0.019207675
## 115 116 117 118 119 120
## 0.019433838 0.019662609 0.019894020 0.020128097 0.020364872 0.020604374
## 121 122 123 124 125 126
## 0.020846632 0.021091677 0.021339540 0.021590251 0.021843842 0.022100345
## 127 128 129 130 131 132
## 0.022359791 0.022622212 0.022887640 0.023156109 0.023427652 0.023702302
## 133 134 135 136 137 138
## 0.023980093 0.024261058 0.024545233 0.024832651 0.025123348 0.025417360
## 139 140 141 142 143 144
## 0.025714722 0.026015469 0.026319639 0.026627268 0.026938394 0.027253052
## 145 146 147 148 149 150
## 0.027571282 0.027893122 0.028218609 0.028547783 0.028880683 0.029217348
## 151 152 153 154 155 156
## 0.029557817 0.029902133 0.030250334 0.030602462 0.030958558 0.031318663
## 157 158 159 160 161 162
## 0.031682821 0.032051072 0.032423461 0.032800029 0.033180821 0.033565881
## 163 164 165 166 167 168
## 0.033955252 0.034348979 0.034747107 0.035149682 0.035556750 0.035968356
## 169 170 171 172 173 174
## 0.036384547 0.036805370 0.037230872 0.037661101 0.038096104 0.038535931
## 175 176 177 178 179 180
## 0.038980630 0.039430251 0.039884842 0.040344455 0.040809138 0.041278943
## 181 182 183 184 185 186
## 0.041753922 0.042234125 0.042719604 0.043210413 0.043706603 0.044208227
## 187 188 189 190 191 192
## 0.044715340 0.045227994 0.045746244 0.046270145 0.046799751 0.047335119
## 193 194 195 196 197 198
## 0.047876303 0.048423360 0.048976347 0.049535320 0.050100336 0.050671453
## 199 200 201 202 203 204
## 0.051248730 0.051832224 0.052421995 0.053018101 0.053620602 0.054229558
## 205 206 207 208 209 210
## 0.054845029 0.055467075 0.056095758 0.056731138 0.057373278 0.058022239
## 211 212 213 214 215 216
## 0.058678083 0.059340874 0.060010674 0.060687547 0.061371555 0.062062764
## 217 218 219 220 221 222
## 0.062761237 0.063467039 0.064180235 0.064900890 0.065629069 0.066364838
## 223 224 225 226 227 228
## 0.067108264 0.067859413 0.068618350 0.069385144 0.070159861 0.070942569
## 229 230 231 232 233 234
## 0.071733335 0.072532227 0.073339313 0.074154662 0.074978343 0.075810423
## 235 236 237 238 239 240
## 0.076650972 0.077500059 0.078357754 0.079224126 0.080099244 0.080983178
## 241 242 243 244 245 246
## 0.081876000 0.082777777 0.083688582 0.084608483 0.085537553 0.086475860
## 247 248 249 250 251 252
## 0.087423476 0.088380472 0.089346918 0.090322885 0.091308444 0.092303666
## 253 254 255 256 257 258
## 0.093308621 0.094323381 0.095348016 0.096382597 0.097427195 0.098481881
## 259 260 261 262 263 264
## 0.099546724 0.100621796 0.101707168 0.102802908 0.103909088 0.105025777
## 265 266 267 268 269 270
## 0.106153046 0.107290963 0.108439599 0.109599021 0.110769301 0.111950505
## 271 272 273 274 275 276
## 0.113142702 0.114345960 0.115560348 0.116785932 0.118022780 0.119270958
## 277 278 279 280 281 282
## 0.120530532 0.121801568 0.123084133 0.124378289 0.125684103 0.127001638
## 283 284 285 286 287 288
## 0.128330957 0.129672123 0.131025199 0.132390246 0.133767325 0.135156497
## 289 290 291 292 293 294
## 0.136557820 0.137971356 0.139397160 0.140835292 0.142285807 0.143748763
## 295 296 297 298 299 300
## 0.145224213 0.146712213 0.148212816 0.149726074 0.151252038 0.152790760
## 301
## 0.154342290
colores <- NULL
colores[datos$Exito == 0] <- "green"
colores[datos$Exito == 1] <- "red"
plot(datos$ï..AlimentoDiario, datos$Exito, pch=21, bg= colores,
xlab = "Cantidad de alimento diario", ylab = "Camarones Aprobados")
legend("topleft", c("Fracaso", "Aprobado"), pch=21, col = c("green", "red") )
lines(datos_nuevos$ï..AlimentoDiario, probabilidades, col="purple", lwd=2)