library(haven)
data <- read_dta("Data1_R.dta")
View(data)Clase 25 de mayo
Cambiar directorio de trabajo. setwd(“C:/Users/JEFERSON CHAMBA/Desktop/Clase 25 de mayo”)
Instalar paquetes si es necesario
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ón
Cargar librerías library(MASS) library(car) library(margins) library(pROC) library(ggplot2) library(ggplot) library(caret) library(haven)
Leer los datos (ajustar la ruta del archivo)
Ver las primeras filas de la base de datos
head(data)# A tibble: 6 × 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]
# ℹ 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: 4Análisis
-El ingreso y los años de escolaridad reflejan no tener relevancia estadística para este modelo es decir estas variables o datos no ayudan a entender o comprender la probabilidad de que una mujer ecuatoriana sufra depresión post parto o ayude a contrarestar este efecto.
La edad es estadisticamente significtiva, lo que se traduce en que, las mujeres con mayor edad tienen probabilidad de sufrir de depresión post parto.
El número de hijos también es estadisticamente significativa por lo que se deduce que al aumentar el número de hijos las mujeres tienden a sufrir depresión post parto.
La etnia es estadisticamente significativa, en el modelo tenemo primero la etnia indigena entonces a partir de ella se muestra que las mujeres indigenas tienen mayor probabilidad de sufrir depresión post parto.
En nuestro modelo por último tenemos la varible área rural la cual presenta evidencia estadistica significativa mostrando que las mujeres del area rural tienen mas probabilidad de presentar depresión post parto a comparación de las mujeres del area 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)
aic_logit <- AIC(modelo_logit)
aic_probit <- AIC(modelo_probit)El criterio de Información Bayesiano (BIC)
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 Conclusión
El modelo con menor AIC/BIC es el preferido, ya que tiene mejor ajuste, es decir 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.0129Análisis
Los años de escolaridad y los ingresos para este modelo arrojan no tener evidencia estadistica significativa. Por lo cual no ayudan a inferir si ayuda o empeora el tema de la depresión post parto de las mujeres ecuatorianas.
El area si es estadisticamente significativa, nos presenta que las mujeres del area rural tienen 1.84% más probabilidad de sufrir depresión post parto.
Edad es una variable igualmente con evidencia significativa en terminos estadisticos muestra que las mujeres aumentan la probabilidad de sufrir depresión post parto en un 0.56%.
Las mujeres indigenas tambien con la variable etnia muestran esta preocupación de sufrir depresión post parto con una probabilidad del 5.92%
-Por último y no menos importante tenemos que la variable número de hijos la cual tambien es estadisticamente significativa aumenta la probabilidad de sufrir depresión post parto de las mujeres en un 0.66%.
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.0132Análisis
Los años de escolaridad y los ingresos para este modelo arrojan no tener evidencia estadistica significativa. Por lo cual no ayudan a inferir si ayuda o empeora el tema de la depresión post parto de las mujeres ecuatorianas.
El area si es estadisticamente significativa, nos presenta que las mujeres del area rural tienen 1.88% más probabilidad de sufrir depresión post parto.
Edad es una variable igualmente con evidencia significativa en terminos estadisticos muestra que las mujeres aumentan la probabilidad de sufrir depresión post parto en un 0.57%.
Las mujeres indigenas tambien con la variable etnia muestran esta preocupación de sufrir depresión post parto con una probabilidad del 6.06%
-Por último y no menos importante tenemos que la variable número de hijos la cual tambien es estadisticamente significativa aumenta la probabilidad de sufrir depresión post parto de las mujeres en un 0.68%.
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 con depresión post parto ocurra de acuerdo en el modelo logit es del 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 Deber aumentar área bajo la curva (AUC) - Modelo logit 0.5836 Estado actual
modelo_logit2 <- glm(depresion_pp ~ lingrl + anios_esc + edad + etnia + area + planf + form_parto + c_posparto + regionx1 + regionx3,
data = data, family = binomial(link = "logit"))Resumen modelo logit
summary(modelo_logit2)
Call:
glm(formula = depresion_pp ~ lingrl + anios_esc + edad + etnia +
area + planf + form_parto + c_posparto + regionx1 + regionx3,
family = binomial(link = "logit"), data = data)
Coefficients:
Estimate Std. Error z value Pr(>|z|)
(Intercept) -2.016738 0.113114 -17.829 < 2e-16 ***
lingrl -0.010240 0.007292 -1.404 0.1602
anios_esc -0.010391 0.004840 -2.147 0.0318 *
edad 0.038099 0.002797 13.622 < 2e-16 ***
etnia 0.324727 0.064052 5.070 3.98e-07 ***
area 0.096475 0.043506 2.217 0.0266 *
planf -0.617628 0.039673 -15.568 < 2e-16 ***
form_parto -0.264505 0.040322 -6.560 5.39e-11 ***
c_posparto 0.049042 0.040276 1.218 0.2234
regionx1 0.346036 0.044435 7.788 6.83e-15 ***
regionx3 0.279774 0.057210 4.890 1.01e-06 ***
---
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: 16770 on 16440 degrees of freedom
AIC: 16792
Number of Fisher Scoring iterations: 4Calcular efectos marginales para Logit
library(margins)
marg_logit2 <- margins(modelo_logit2)
summary(marg_logit2) factor AME SE z p lower upper
anios_esc -0.0017 0.0008 -2.1475 0.0318 -0.0033 -0.0002
area 0.0160 0.0072 2.2183 0.0265 0.0019 0.0301
c_posparto 0.0081 0.0067 1.2178 0.2233 -0.0049 0.0212
edad 0.0063 0.0005 13.8204 0.0000 0.0054 0.0072
etnia 0.0538 0.0106 5.0804 0.0000 0.0330 0.0745
form_parto -0.0438 0.0067 -6.5802 0.0000 -0.0568 -0.0307
lingrl -0.0017 0.0012 -1.4046 0.1602 -0.0041 0.0007
planf -0.1023 0.0064 -15.8701 0.0000 -0.1149 -0.0896
regionx1 0.0573 0.0073 7.8210 0.0000 0.0429 0.0716
regionx3 0.0463 0.0095 4.8981 0.0000 0.0278 0.0649Obtener predicciones del modelo logit
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 de confusión
print(conf_matrix2) Real
Predicho 0 1
0 12803 3601
1 25 22Calcular exactitud
exactitud2 <- sum(diag(conf_matrix2)) / sum(conf_matrix2)
cat("Exactitud del modelo logit:", exactitud2, "\n")Exactitud del modelo logit: 0.7795879 Las variables consideradas para explicar la probabilidad de que las mujeres con depresión post parto ocurra de acuerdo en el modelo logit es del 77.96%
Curva
library(pROC)
roc_logit2 <- roc(data$depresion_pp, predict(modelo_logit2, type = "response"))Setting levels: control = 0, case = 1Setting direction: controls < caseslibrary(ggplot2)
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",
x = "1 - Especificidad",
y = "Sensibilidad") +
theme_minimal()
Mostrar el area bajo la curva
auc_logit2 <- auc(roc_logit2)
cat("Área bajo la curva 2 (AUC) - Modelo Logit 2:", auc_logit2, "\n")Área bajo la curva 2 (AUC) - Modelo Logit 2: 0.6265065 Se añadió las siguientes variables:
- planf que se refiere a que si planificaron familiarmente su cuidado para tener hijos esta variable es estadisticamente significativa y tiene la probabilidad de restar la depresión post parto.
- c_posparto que se refiere a si tuvo algún control despues del parto esta variable no es estadisticamente significativa por tanto no ayuda a inferir si suma o resta.
- region x1 que se refiere a las mujeres de la sierra se tomo esta variable ya que en la sierra se encuentra un número considerables de mujeres indigenas esta variables si es estadisticamente significativa y tiene la probabilidad de aumentar la depresión post parto.
- region x3 que se refiere a las mujeres de la amazonia se tomo esta variable ya que en la sierra se encuentra un número considerables de mujeres indigenas esta variables si es estadisticamente significativa y tiene la probabilidad de aumentar la depresión post parto.
En conclusión se pudo aumentar el area bajo la curva de 0.5836 a 0.6265