alt text

Métodos estadísticos I

Tarea 11 - Regresión Logística

En este ejercicio aplicaremos una regresión logística a los siguientes datos:


datos <- data.frame(f_cond = c(1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0, 0, 0, 0, 0, 
    0, 0, 0, 0, 0), loans_assets = c(0.64, 1.04, 0.66, 0.8, 0.69, 0.74, 0.63, 
    0.75, 0.56, 0.65, 0.55, 0.46, 0.72, 0.43, 0.52, 0.54, 0.3, 0.67, 0.51, 0.79), 
    expenses_assets = c(0.13, 0.1, 0.11, 0.09, 0.11, 0.14, 0.12, 0.12, 0.16, 
        0.12, 0.1, 0.08, 0.08, 0.08, 0.07, 0.08, 0.09, 0.07, 0.09, 0.13))

datos
##    f_cond loans_assets expenses_assets
## 1       1         0.64            0.13
## 2       1         1.04            0.10
## 3       1         0.66            0.11
## 4       1         0.80            0.09
## 5       1         0.69            0.11
## 6       1         0.74            0.14
## 7       1         0.63            0.12
## 8       1         0.75            0.12
## 9       1         0.56            0.16
## 10      1         0.65            0.12
## 11      0         0.55            0.10
## 12      0         0.46            0.08
## 13      0         0.72            0.08
## 14      0         0.43            0.08
## 15      0         0.52            0.07
## 16      0         0.54            0.08
## 17      0         0.30            0.09
## 18      0         0.67            0.07
## 19      0         0.51            0.09
## 20      0         0.79            0.13

La idea es modelar la condición financiera (f_cond) que es una variable dicotómica, a partir de las variables:

*Razón de prestámos totales entre el total de activos (loans_assets)

*Razón de gastos totales entre el total de activos (expensses_assets )

Antes de ajustar el modelo logístico veamos gráficamente los datos:


par(mfrow = c(1, 2))
plot(datos$loans_assets, datos$f_cond, main = "Prestamos y Rentas", xlab = expression(frac("Prestamos+Rentas", 
    "Activos")), ylab = "Condición financiera", cex.lab = 0.8)

plot(datos$expenses_assets, datos$f_cond, main = "Gastos", xlab = expression(frac("Gastos", 
    "Activos")), ylab = "Condición financiera", cex.lab = 0.8)

plot of chunk plot1

Es difícil ver algo así como que sí, quizá, quien sabe; veamos otra representación gráfica donde quizá veamos algo más:


par(mfrow = c(1, 1), mai = c(1.5, 1.5, 1, 1))
plot(datos$loans_assets, datos$expenses_assets, main = "Prestamos y Rentas", 
    xlab = expression(frac("Prestamos+Rentas", "Activos")), ylab = expression(frac("Gastos", 
        "Activos")), cex.lab = 0.8, pch = 16, col = c(2, 4)[as.factor(datos$f_cond)])
legend("topleft", pch = 16, col = c(2, 4), legend = c("1", "0"), bty = "n")

plot of chunk plot2

Aparentemente si son relevantes, ¿pero lo són? y ¿qué tanto?. Ahora utilizaremos un modelo de regresión logística para responder estas preguntas:

mod_log <- glm(f_cond ~ loans_assets + expenses_assets, data = datos, family = binomial(logit))
mod_log
## 
## Call:  glm(formula = f_cond ~ loans_assets + expenses_assets, family = binomial(logit), 
##     data = datos)
## 
## Coefficients:
##     (Intercept)     loans_assets  expenses_assets  
##          -14.19             9.17            79.96  
## 
## Degrees of Freedom: 19 Total (i.e. Null);  17 Residual
## Null Deviance:       27.7 
## Residual Deviance: 12.8  AIC: 18.8

De tal manera que nuestro modelo es:

\[ \frac{P(Condición)}{1-P(Condición)}=\frac{e^{\beta _0+ \beta _1 x_1 + \beta _2 x_2 }}{1+e^{\beta_0+ \beta_1 x_1 + \beta_2 x_2 }} \]

Donde \( x_1=\frac{Prestamos+Rentas}{Activos} \) y \( x_2=\frac{Gastos}{Activos} \)

Ya con lo estimadores de los coeficientes obtenemos el modelo:

\[ \frac{P(Condición)}{1-P(Condición)}=\frac{e^{-14.1876+ 9.1732x_1 + 79.9639 x_2 }}{1+e^{-14.1876+ 9.1732x_1 + 79.9639 x_2 }} \]

Antes de pasar a algo más conviene verificar que el modelo sea significativo, y que cada una de las variables incluidas también lo sean, si bien podemos obtener los estasdísticoas t de Wald, para este caso siendo un modelo múltiple preferiremos una \( X^2 \)

anova(mod_log, test = "Chisq")
## Analysis of Deviance Table
## 
## Model: binomial, link: logit
## 
## Response: f_cond
## 
## Terms added sequentially (first to last)
## 
## 
##                 Df Deviance Resid. Df Resid. Dev Pr(>Chi)   
## NULL                               19       27.7            
## loans_assets     1     7.17        18       20.6   0.0074 **
## expenses_assets  1     7.73        17       12.8   0.0054 **
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Esta es una prueba por orden de inclusión, de tal manera que para este caso prueba con respecto a expenses_assets dado que loans_assets ya está incluido en el modelo, y en este caso ambas variables resultan significativas.

Pasemos a intepretar los coeficientes del modelo, para esto es necesario tomar la exponencial de los coeficientes en su forma linealizada donde tenemos:

\( \beta_1= \) 9635.5494

\( \beta_2= \) 5.3444 × 1034

Es decir que cuando la razón de prestamos aumenta un punto la probabilidad de tener una condición 1 aumenta exp(mod_log$coefficients)[2] veces. Y la interpretación para la razón de gastos es la misma. Hay otra manera de ver estos resultados, utilizando prediciónes




pred_mod <- function(modelo, x1, x2) {
    pred_num <- exp(coef(modelo)[[1]] + (coef(modelo)[[2]] * x1) + (coef(modelo)[[3]] * 
        x2))
    pred <- pred_num/(1 + pred_num)
    p <- pred/(1 + pred)
    p
}

predicciones <- matrix(nrow = 40, ncol = 40, dimnames = list(seq(from = min(datos$loans_assets), 
    to = max(datos$loans_assets), length = 40), seq(from = min(datos$expenses_assets), 
    to = max(datos$expenses_assets), length = 40)))

for (i in 1:40) {
    for (j in 1:40) {
        predicciones[i, j] <- pred_mod(mod_log, as.numeric(dimnames(predicciones)[[1]][i]), 
            as.numeric(dimnames(predicciones)[[2]][j]))
    }
}

heatmap(t(predicciones), Rowv = NA, Colv = NA)

plot of chunk plot3

Básicamente lo que estamos haciendo es predicciones dentro de los rangos de las variables en este caso \( 40\times 40=1600 \) predicciones, sin embargo no son predicciones sobre la razón de momios, sino utilizando el modelo se buscó realizar predicciones sobre la probabilidad de éxito, i.e. f_cond==1. Para esto se utilizó \( p=\frac{modelo}{1+modelo} \). La gráfica que se presenta es un heat map (¿mapa de calor?), básicamente los colores más obscuros representan una mayor probabilidad.

Cabe señalar que además de la regresión logística se pueden utilizar técnicas más nuevas o más viejas, que en términos de predicción pueden ser más flexibles, por ejemplo: Ánálisis Discriminante, Árboles de regresión y bosques.

nota: Las probabilidades obtenidas en este caso quedaron todas debajo de 0.5, nos parece que hubo un error y quedaron a la mitad las probabilidades, sin embarga si este es el caso no afecta el heat map.

mod_log