# Estos paquetes solo se instalan una vez
install.packages("MASS") # Para modelo probit
install.packages("car") # Para pruebas estadísticas
install.packages("margins") # Para calcular efectos marginales
install.packages("pROC") # Para curva ROC
install.packages("ggplot2") # Para visualización
install.packages("caret") # Para matriz de confusión
install.packages("ggplot") # Para visualizaciónCLASE 25 DE MAYO YURY
Quarto
Quarto enables you to weave together content and executable code into a finished document. To learn more about Quarto see https://quarto.org.
Instalar paquetes si es necesario
# Cargar librerías
library(MASS)
library(car)
library(margins)
library(pROC)
library(ggplot2)
library(ggplot)
library(caret)
library(haven)#Leer la base de datos
library(haven)
data <- read_dta("Data1_R.dta")
View(data)Ver las primeras filas de la base de datos
head(data)# A tibble: 6 x 50
area empleo region edad t_hijos nac_vivo_murieron mortinato_2
<dbl+lbl> <dbl+lbl> <dbl+l> <dbl> <dbl> <dbl+lbl> <dbl+lbl>
1 1 [Urbano] 1 [Trabajó al ~ 1 [Sie~ 19 1 0 [No] 0 [No]
2 1 [Urbano] 0 [No trabajó] 1 [Sie~ 23 1 0 [No] 0 [No]
3 1 [Urbano] 1 [Trabajó al ~ 1 [Sie~ 38 5 0 [No] 0 [No]
4 1 [Urbano] 0 [No trabajó] 1 [Sie~ 18 1 0 [No] 0 [No]
5 1 [Urbano] 0 [No trabajó] 1 [Sie~ 21 1 0 [No] 0 [No]
6 1 [Urbano] 1 [Trabajó al ~ 1 [Sie~ 22 1 0 [No] 0 [No]
# i 43 more variables: depresion_pp <dbl+lbl>, intensidad_dpp <dbl+lbl>,
# etnia <dbl+lbl>, f2_s2_216_1 <dbl+lbl>, f2_s2_216_2 <dbl>,
# f2_s2_218_1_a <dbl+lbl>, tiempo_dpp <dbl+lbl>, f2_s5_504a_1 <dbl+lbl>,
# f2_s5_504b_1 <dbl+lbl>, f2_s5_504c_1 <dbl+lbl>, f2_s5_504d_1 <dbl+lbl>,
# f2_s5_504e_1 <dbl+lbl>, f2_s5_504f_1 <dbl+lbl>, f2_s5_504g_1 <dbl+lbl>,
# f2_s5_504h_1 <dbl+lbl>, f2_s5_504i_1 <dbl+lbl>, f2_s5_504j_1 <dbl+lbl>,
# f2_s5_504k_1 <dbl+lbl>, est_civil <dbl+lbl>, q_usted <dbl+lbl>, ...Revisar estructura de los datos
str(data)##Ejemplo 1: Modelos con variable dependiente dicotómica
##Modelos logit y probit
Ajustar el modelo logit
modelo_logit <- glm(depresion_pp ~ lingrl + anios_esc + edad + t_hijos + etnia + area,
data = data, family = binomial(link = "logit"))Ajustar el modelo probit
modelo_probit <- glm(depresion_pp ~ lingrl + anios_esc + edad + t_hijos + etnia + area,
data = data, family = binomial(link = "probit"))Resumen modelo Logit
summary(modelo_logit)
Call:
glm(formula = depresion_pp ~ lingrl + anios_esc + edad + t_hijos +
etnia + area, family = binomial(link = "logit"), data = data)
Coefficients:
Estimate Std. Error z value Pr(>|z|)
(Intercept) -2.3377859 0.1015521 -23.021 < 2e-16 ***
lingrl 0.0006157 0.0071763 0.086 0.9316
anios_esc -0.0078052 0.0049109 -1.589 0.1120
edad 0.0333503 0.0032243 10.344 < 2e-16 ***
t_hijos 0.0391392 0.0189765 2.063 0.0392 *
etnia 0.3502255 0.0605997 5.779 7.5e-09 ***
area 0.1089295 0.0425378 2.561 0.0104 *
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
(Dispersion parameter for binomial family taken to be 1)
Null deviance: 17346 on 16450 degrees of freedom
Residual deviance: 17105 on 16444 degrees of freedom
AIC: 17119
Number of Fisher Scoring iterations: 4#Análisis
El ingreso del hogar y los años de escolaridad no resultan estadísticamente significativos, lo que significa que no contribuyen de manera relevante a explicar la probabilidad de que las mujeres ecuatorianas presenten depresión posparto, según los datos analizados en este modelo. La edad es estadisticamente significativa, lo que se traduce en que las mujeres de mayor edad, tienen probabilidad de sufrir depresión posparto. El número de hijos presenta significancia estadistica, es decir que a medida que aunmenta el número de hijos, aunmenta la probabilidad de sufran depresión psparto La etnia es una variable estadisticamente significativa en el modelo, lo que significa que, las mujeres indigenas tienen mayor probabilidad de tener depresión posparto. Las mujeres del área rural tienen mayor probabilidad de sufrir depresipon posparto, en comparación a las mujeres del área urbana.
Resumen modelo probit
summary(modelo_probit)
Call:
glm(formula = depresion_pp ~ lingrl + anios_esc + edad + t_hijos +
etnia + area, family = binomial(link = "probit"), data = data)
Coefficients:
Estimate Std. Error z value Pr(>|z|)
(Intercept) -1.401e+00 5.852e-02 -23.942 < 2e-16 ***
lingrl 3.942e-05 4.170e-03 0.009 0.99246
anios_esc -4.481e-03 2.861e-03 -1.566 0.11733
edad 1.958e-02 1.890e-03 10.363 < 2e-16 ***
t_hijos 2.334e-02 1.123e-02 2.078 0.03774 *
etnia 2.078e-01 3.585e-02 5.796 6.8e-09 ***
area 6.431e-02 2.452e-02 2.623 0.00872 **
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
(Dispersion parameter for binomial family taken to be 1)
Null deviance: 17346 on 16450 degrees of freedom
Residual deviance: 17103 on 16444 degrees of freedom
AIC: 17117
Number of Fisher Scoring iterations: 4Comparar AIC y BIC de ambos modelos
El criterio de información de Akaike (AIC) es un método matemático para evaluar el ajuste de un modelo a los datos a partir de los cuales se generó.
El criterio de Información Bayesiano (BIC) es una medida de bondad de ajuste de un modelo estadístico, y es a muenudo utilizado como un criterio para para la selección de modelos entre un conjunto finito de modelos.
##Comparar AIC de los modelos
aic_logit <- AIC(modelo_logit)
aic_probit <- AIC(modelo_probit)##Comparar BIC de los modelos
bic_logit <- BIC(modelo_logit)
bic_probit <- BIC(modelo_probit)Mostrar resultados
cat("AIC Logit:", aic_logit, " | AIC Probit:", aic_probit, "\n")AIC Logit: 17119.38 | AIC Probit: 17117.13 cat("BIC Logit:", bic_logit, " | BIC Probit:", bic_probit, "\n")BIC Logit: 17173.34 | BIC Probit: 17171.09 Cinclusión
El modelo con menor criterio de información, tanto AIC como BIC, es el modelo probit.
Calcular efectos marginales para Logit
library(margins)
marg_logit <- margins(modelo_logit)
summary(marg_logit) factor AME SE z p lower upper
anios_esc -0.0013 0.0008 -1.5897 0.1119 -0.0029 0.0003
area 0.0184 0.0072 2.5619 0.0104 0.0043 0.0325
edad 0.0056 0.0005 10.4239 0.0000 0.0046 0.0067
etnia 0.0592 0.0102 5.7944 0.0000 0.0392 0.0793
lingrl 0.0001 0.0012 0.0858 0.9316 -0.0023 0.0025
t_hijos 0.0066 0.0032 2.0632 0.0391 0.0003 0.0129Analisis
Los años de escolaridad y el ingreso no son estadisticamente significativos, es decir no ayudan a explicar la probabilidad de que las mujeres ecuatorianas padezcan de depresión posparto.
El área al ser una variable estadísticamente significativa 0.0184*100, nos explica que las mujeres del área rural tienen el 1.84% tienen más probabilidades de sufrir depresión posparto.
Un año adicional en las mujeres ecuatorias aunmentan la probabilidad en 0,56% de sufrir de depresipon posparto.
Una mujer indígena tiene el 5.92% de probabilidad de sufrir de depresión posparto.
A medida que aunmenta el número de hijos aunmenta la probabilidad en 0.66% que padescan de depresipon posparto
Calcular efectos marginales para Probit
marg_probit <- margins(modelo_probit)
summary(marg_probit) factor AME SE z p lower upper
anios_esc -0.0013 0.0008 -1.5664 0.1173 -0.0029 0.0003
area 0.0188 0.0072 2.6238 0.0087 0.0047 0.0328
edad 0.0057 0.0005 10.4392 0.0000 0.0046 0.0068
etnia 0.0606 0.0104 5.8096 0.0000 0.0402 0.0811
lingrl 0.0000 0.0012 0.0095 0.9925 -0.0024 0.0024
t_hijos 0.0068 0.0033 2.0783 0.0377 0.0004 0.0132Analisis
Los años de escolaridad y el ingreso no son estadisticamente significativos, es decir no ayudan a explicar la probabilidad de que las mujeres ecuatorianas padezcan de depresión posparto.
El área al ser una variable aunmenta en 1.88% tienen más probabilidades de sufrir depresión posparto.
Un año adicional en las mujeres ecuatorias aunmentan la probabilidad en 0,57% de sufrir de depresipon posparto.
Una mujer indígena tiene el 6.06% de probabilidad de sufrir de depresión posparto.
A medida que aunmenta el número de hijos aunmenta la probabilidad en 0.68% que padescan de depresipon posparto
###Para el modelo logit podemos calcular la matriz de confusión
Obtener predicciones del modelo logit
pred_logit <- ifelse(predict(modelo_logit, type = "response") > 0.5, 1, 0)Crear matriz de confusión
conf_matrix <- table(Predicho = pred_logit, Real = data$depresion_pp)Mostrar matriz de confusión
print(conf_matrix) Real
Predicho 0 1
0 12828 3623Calcular exactitud
exactitud <- sum(diag(conf_matrix)) / sum(conf_matrix)
cat("Exactitud del modelo logit:", exactitud, "\n")Exactitud del modelo logit: 0.7797702 Las variables consideradas para explicar la probabilidad de que las mujeres de depresión posparto ocurra de acuerdo al modelo ligin es de 77,98%
Calcular la curva ROC para el modelo logit (ayuda a ver con cuanta exactitud el modelo está prediciendo los resultados)
library(pROC)Type 'citation("pROC")' for a citation.
Adjuntando el paquete: 'pROC'The following objects are masked from 'package:stats':
cov, smooth, varroc_logit <- roc(data$depresion_pp, predict(modelo_logit, type = "response"))Setting levels: control = 0, case = 1Setting direction: controls < casesGraficar la curva ROC
library(ggplot2)
ggplot() +
geom_line(aes(x = roc_logit$specificities, y = roc_logit$sensitivities), color = "blue") +
geom_abline(linetype = "dashed", color = "red") +
labs(title = "Curva ROC - Modelo Logit",
x = "1 - Especificidad",
y = "Sensibilidad") +
theme_minimal()
Mostrar el área bajo la curva (AUC)
auc_logit <- auc(roc_logit)
cat("Área bajo la curva (AUC) - Modelo Logit:", auc_logit, "\n")Área bajo la curva (AUC) - Modelo Logit: 0.5836005 #Curva ROC: Muestra el rendimiento del modelo en diferentes umbrales de clasificación. #AUC > 0.7 indica un buen modelo predictivo.
Convertir variables haven_labelled a numeric usando haven::as_factor
data$empleo <- as.numeric(haven::as_factor(data$empleo))
data$nac_vivo_murieron <- as.numeric(haven::as_factor(data$nac_vivo_murieron))
data$mortinato_2 <- as.numeric(haven::as_factor(data$mortinato_2))
data$est_civil <- as.numeric(haven::as_factor(data$est_civil))Ajustar un nuevo modelo logit con más variables
modelo_logit2 <- glm(depresion_pp ~ lingrl + anios_esc + edad + t_hijos + etnia + area +
empleo + nac_vivo_murieron + mortinato_2 + est_civil,
data = data, family = binomial(link = "logit"))Resumen del nuevo modelo logit con variables adicionales
summary(modelo_logit2)
Call:
glm(formula = depresion_pp ~ lingrl + anios_esc + edad + t_hijos +
etnia + area + empleo + nac_vivo_murieron + mortinato_2 +
est_civil, family = binomial(link = "logit"), data = data)
Coefficients:
Estimate Std. Error z value Pr(>|z|)
(Intercept) -4.915478 0.192401 -25.548 < 2e-16 ***
lingrl -0.029765 0.012826 -2.321 0.02030 *
anios_esc 0.002179 0.005060 0.431 0.66669
edad 0.027485 0.003387 8.115 4.85e-16 ***
t_hijos 0.037720 0.019392 1.945 0.05176 .
etnia 0.295325 0.063081 4.682 2.85e-06 ***
area 0.136516 0.043461 3.141 0.00168 **
empleo 0.223930 0.071417 3.136 0.00172 **
nac_vivo_murieron 1.502174 0.091312 16.451 < 2e-16 ***
mortinato_2 0.844528 0.105359 8.016 1.09e-15 ***
est_civil -0.051580 0.043270 -1.192 0.23325
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
(Dispersion parameter for binomial family taken to be 1)
Null deviance: 17346 on 16450 degrees of freedom
Residual deviance: 16744 on 16440 degrees of freedom
AIC: 16766
Number of Fisher Scoring iterations: 4Cinclusión
En el nuevo modelo logit con variables ampliadas, se observa una mejora en el ajuste general (AIC bajó a 16766). Además, varias variables resultaron estadísticamente significativas:
La edad, el número de hijos, el área de residencia, etnia, empleo, y la pérdida de hijos vivos o mortinatos se asocian significativamente con una mayor probabilidad de depresión posparto.
En cambio, los años de escolaridad y el estado civil no resultaron significativos en este modelo.
De forma interesante, el ingreso del hogar (lingrl) sí fue significativo y con efecto negativo, lo que indica que mayores ingresos estarían asociados con una menor probabilidad de sufrir depresión posparto.
Calcular efectos marginales para el modelo logit extendido
marg_logit2 <- margins(modelo_logit2)Ver resumen de los efectos marginales
summary(marg_logit2) factor AME SE z p lower upper
anios_esc 0.0004 0.0008 0.4307 0.6667 -0.0013 0.0020
area 0.0224 0.0071 3.1431 0.0017 0.0084 0.0364
edad 0.0045 0.0006 8.1515 0.0000 0.0034 0.0056
empleo 0.0368 0.0117 3.1377 0.0017 0.0138 0.0598
est_civil -0.0085 0.0071 -1.1922 0.2332 -0.0224 0.0055
etnia 0.0486 0.0104 4.6892 0.0000 0.0283 0.0689
lingrl -0.0049 0.0021 -2.3216 0.0203 -0.0090 -0.0008
mortinato_2 0.1389 0.0172 8.0642 0.0000 0.1051 0.1726
nac_vivo_murieron 0.2470 0.0146 16.9456 0.0000 0.2184 0.2756
t_hijos 0.0062 0.0032 1.9457 0.0517 -0.0000 0.0125Analisis efectos (multiplicarlos por 100 y analizarlos)
Los efectos marginales del modelo logit ampliado permiten observar el impacto porcentual de cada variable sobre la probabilidad de que una mujer ecuatoriana presente depresión posparto.
Las pérdidas de hijos vivos y mortinatos son los factores más influyentes, aumentando la probabilidad en 24.7% y 13.9% respectivamente.
Otras variables significativas incluyen: etnia indígena (↑ 4.86%), estar empleada (↑ 3.68%), vivir en zona rural (↑ 2.24%), y la edad (↑ 0.45% por año adicional).
El ingreso del hogar tiene un efecto negativo leve, reduciendo la probabilidad en 0.49%, mientras que el número de hijos tiene un efecto marginal con un valor p ligeramente por encima de 0.05.
En contraste, los años de escolaridad y el estado civil no mostraron efectos significativos.
Obtener predicciones del modelo logit extendido (umbral de 0.5)
pred_logit2 <- ifelse(predict(modelo_logit2, type = "response") > 0.5, 1, 0)Crear matriz de confusión
conf_matrix2 <- table(Predicho = pred_logit2, Real = data$depresion_pp)Mostrar matriz
print(conf_matrix2) Real
Predicho 0 1
0 12626 3344
1 202 279Calcular exactitud
exactitud2 <- sum(diag(conf_matrix2)) / sum(conf_matrix2)
cat("Exactitud del modelo logit 2:", exactitud2, "\n")Exactitud del modelo logit 2: 0.7844508 Analisis
La exactitud del modelo logit ampliado es del 78.45%, lo que indica que el modelo logra clasificar correctamente a aproximadamente 8 de cada 10 mujeres en cuanto a si presentan o no depresión posparto. Esta exactitud representa una ligera mejora respecto al modelo inicial, lo que sugiere que las variables adicionales incorporadas aportan valor predictivo.
Calcular la curva ROC para modelo_logit2
roc_logit2 <- roc(data$depresion_pp, predict(modelo_logit2, type = "response"))Setting levels: control = 0, case = 1Setting direction: controls < casesGraficar la curva ROC
ggplot() +
geom_line(aes(x = roc_logit2$specificities, y = roc_logit2$sensitivities), color = "blue") +
geom_abline(linetype = "dashed", color = "red") +
labs(title = "Curva ROC - Modelo Logit Ampliado",
x = "1 - Especificidad",
y = "Sensibilidad") +
theme_minimal()
Calcular AUC
auc_logit2 <- auc(roc_logit2)
cat("Área bajo la curva (AUC) - Modelo Logit 2:", auc_logit2, "\n")Área bajo la curva (AUC) - Modelo Logit 2: 0.6152535 Análisis
El área bajo la curva (AUC) del modelo logit ampliado es de 0.6152, lo que indica que el modelo tiene una capacidad moderada para discriminar entre mujeres con y sin depresión posparto. Si bien este valor representa una ligera mejora frente al modelo inicial (AUC ≈ 0.58), aún se encuentra por debajo del umbral de 0.7, considerado como un nivel aceptable de discriminación. Esto sugiere que, aunque las variables adicionales aportan información relevante.
Selección automática de variables para mejorar el modelo
modelo_logit_step <- step(modelo_logit2)Start: AIC=16765.9
depresion_pp ~ lingrl + anios_esc + edad + t_hijos + etnia +
area + empleo + nac_vivo_murieron + mortinato_2 + est_civil
Df Deviance AIC
- anios_esc 1 16744 16764
- est_civil 1 16745 16765
<none> 16744 16766
- t_hijos 1 16748 16768
- lingrl 1 16749 16769
- empleo 1 16754 16774
- area 1 16754 16774
- etnia 1 16765 16785
- mortinato_2 1 16805 16825
- edad 1 16809 16829
- nac_vivo_murieron 1 17013 17033
Step: AIC=16764.08
depresion_pp ~ lingrl + edad + t_hijos + etnia + area + empleo +
nac_vivo_murieron + mortinato_2 + est_civil
Df Deviance AIC
- est_civil 1 16745 16763
<none> 16744 16764
- t_hijos 1 16748 16766
- lingrl 1 16749 16767
- empleo 1 16754 16772
- area 1 16755 16773
- etnia 1 16765 16783
- mortinato_2 1 16805 16823
- edad 1 16810 16828
- nac_vivo_murieron 1 17014 17032
Step: AIC=16763.4
depresion_pp ~ lingrl + edad + t_hijos + etnia + area + empleo +
nac_vivo_murieron + mortinato_2
Df Deviance AIC
<none> 16745 16763
- t_hijos 1 16749 16765
- lingrl 1 16751 16767
- empleo 1 16755 16771
- area 1 16756 16772
- etnia 1 16766 16782
- mortinato_2 1 16806 16822
- edad 1 16810 16826
- nac_vivo_murieron 1 17016 17032Ver resumen del modelo optimizado
summary(modelo_logit_step)
Call:
glm(formula = depresion_pp ~ lingrl + edad + t_hijos + etnia +
area + empleo + nac_vivo_murieron + mortinato_2, family = binomial(link = "logit"),
data = data)
Coefficients:
Estimate Std. Error z value Pr(>|z|)
(Intercept) -4.923602 0.176905 -27.832 < 2e-16 ***
lingrl -0.028929 0.012624 -2.292 0.02193 *
edad 0.026656 0.003295 8.090 5.96e-16 ***
t_hijos 0.035072 0.018609 1.885 0.05947 .
etnia 0.288404 0.062807 4.592 4.39e-06 ***
area 0.139186 0.042836 3.249 0.00116 **
empleo 0.221076 0.071259 3.102 0.00192 **
nac_vivo_murieron 1.500702 0.091009 16.490 < 2e-16 ***
mortinato_2 0.842235 0.105228 8.004 1.21e-15 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
(Dispersion parameter for binomial family taken to be 1)
Null deviance: 17346 on 16450 degrees of freedom
Residual deviance: 16745 on 16442 degrees of freedom
AIC: 16763
Number of Fisher Scoring iterations: 4Calcular efectos marginales para modelo optimizado
marg_step <- margins(modelo_logit_step)
summary(marg_step) factor AME SE z p lower upper
area 0.0229 0.0070 3.2515 0.0011 0.0091 0.0367
edad 0.0044 0.0005 8.1258 0.0000 0.0033 0.0054
empleo 0.0364 0.0117 3.1046 0.0019 0.0134 0.0593
etnia 0.0474 0.0103 4.5989 0.0000 0.0272 0.0676
lingrl -0.0048 0.0021 -2.2925 0.0219 -0.0088 -0.0007
mortinato_2 0.1385 0.0172 8.0522 0.0000 0.1048 0.1722
nac_vivo_murieron 0.2468 0.0145 16.9883 0.0000 0.2183 0.2753
t_hijos 0.0058 0.0031 1.8852 0.0594 -0.0002 0.0118Obtener predicciones para el modelo optimizado
pred_step <- ifelse(predict(modelo_logit_step, type = "response") > 0.5, 1, 0)Crear matriz de confusión
conf_matrix_step <- table(Predicho = pred_step, Real = data$depresion_pp)Mostrar matriz
print(conf_matrix_step) Real
Predicho 0 1
0 12628 3343
1 200 280Calcular exactitud
exactitud_step <- sum(diag(conf_matrix_step)) / sum(conf_matrix_step)
cat("Exactitud del modelo optimizado (step):", exactitud_step, "\n")Exactitud del modelo optimizado (step): 0.7846332 Calcular la curva ROC para el modelo optimizado
roc_step <- roc(data$depresion_pp, predict(modelo_logit_step, type = "response"))Setting levels: control = 0, case = 1Setting direction: controls < casesGraficar la curva ROC
ggplot() +
geom_line(aes(x = roc_step$specificities, y = roc_step$sensitivities), color = "blue") +
geom_abline(linetype = "dashed", color = "red") +
labs(title = "Curva ROC - Modelo Logit Optimizado",
x = "1 - Especificidad",
y = "Sensibilidad") +
theme_minimal()
Calcular AUC
auc_step <- auc(roc_step)
cat("Área bajo la curva (AUC) - Modelo Logit Optimizado:", auc_step, "\n")Área bajo la curva (AUC) - Modelo Logit Optimizado: 0.6146826 Analisis
El modelo logit optimizado mediante la función step() logró seleccionar las variables con mayor capacidad explicativa, logrando un ajuste eficiente con el menor AIC registrado (16763). Si bien la exactitud del modelo se mantiene alrededor del 78.4%, la área bajo la curva ROC (AUC = 0.6147) muestra una mejora respecto al modelo inicial (AUC = 0.58), aunque todavía se encuentra por debajo del umbral de 0.7 deseado.