Analítica de datos aplicada a estudios sobre desarrollo

Análisis de Regresión

REGRESIÓN Y CORRELACIÓN

Proceso de analítica

Wickham, H. y otros (2023)

Modelos de analítica

Modelando datos

Diagramas de dispersión

Permite identificar de una forma visual la existencia de asociación lineal o no lineal.

Diagramas de dispersión

Descripción del conjunto de datos precios_viviendas.rds1

Variable Descripción
precio Precio de venta, en millones de COP
habitaciones Número de habitaciones
banos Número de baños
área Área de la casa en metros cuadrados
anio_construccion Año en que se construyó la casa
lote Área de toda la propiedad
zona Ubicación de la vivienda

Diagramas de dispersión

# Cargar librerías necesarias
library(pacman)
p_load(tidyverse, psych)

# Cargar datos desde GitHub
datos <- "https://github.com/jgbabativam/AnaDatos/raw/main/datos/precios_viviendas.rds"
load(url(datos))

precios_viviendas |> 
  select(-zona) |> 
  pivot_longer(cols = -precio, names_to = "variable", values_to = "valores") |> 
  ggplot(aes(y = precio, x = valores)) +
  geom_point(color = "steelblue") +
  facet_wrap(~variable, scales = "free") +
  theme_bw()

Diagramas de dispersión

pairs.panels(precios_viviendas |> select(-zona))

Matriz de correlación

(r <- cor(precios_viviendas |> select(-zona), use="complete.obs"))
                     precio habitaciones     banos   area_m2 anio_construccion
precio            1.0000000    0.4116450 0.5938422 0.6713036         0.2176978
habitaciones      0.4116450    1.0000000 0.6203248 0.5671557         0.1817579
banos             0.5938422    0.6203248 1.0000000 0.6753512         0.3230760
area_m2           0.6713036    0.5671557 0.6753512 1.0000000         0.1706275
anio_construccion 0.2176978    0.1817579 0.3230760 0.1706275         1.0000000
tamano_lote       0.5340996    0.2218340 0.2927030 0.4020429        -0.1221422
                  tamano_lote
precio              0.5340996
habitaciones        0.2218340
banos               0.2927030
area_m2             0.4020429
anio_construccion  -0.1221422
tamano_lote         1.0000000
r <- cor(precios_viviendas |> select(-zona), use="complete.obs")
corrplot(r, method = "color", type = "upper", addCoef.col = "black", tl.col = "black", tl.srt = 45,
         col = colorRampPalette(c("darkred", "white", "steelblue"))(200))

ESTUDIO DE CASO

DASS 21

  • El instrumento del DASS 21 permite construir una escala de Depresión, Ansiedad y Estrés (DASS-21). Investigue más sobre su contrucción y propiedades psicométricas. Una versión del instrumento puede ser consultada aquí

  • Explore el conjunto del datos DASS21.sav el cual contiene los resultados para una muestra de 800 personas de Colombia realizada en el año 2022.

library(haven)

url <- "https://github.com/jgbabativam/AnaDatos/raw/main/datos/DASS21.sav"
dass <- read_sav(url)

Puede usar lapply(dass, function(x) attributes(x)$label) para ver las etiquetas de las preguntas.

Punto 1 - Taller

  1. Haga el diagrama de dispersión y calcule la correlación entre las variables cuantitativas de nivel de depresión, estrés y ansiedad.

  2. ¿Considera que el grado de asociación se diferencia entre hombres y mujeres?, haga los gráficos de dispersión segmentados por sexo

  3. Realice los análisis que le permitan concluir sobre la asociación entre la depresión y la satisfacción con la vivienda, trabajo, amigos, vecinos y el barrio.

  4. Teniendo en cuenta que las variables sobre la participación en actividades no son cuantitativas, investigue y discuta sobre la forma en que podría identificarse alguna asociación con la depresión.

ANÁLISIS DE REGRESIÓN

Conceptos básicos

Denominar a \(\mathbf{X}=(\mathbf{x}_1,\ldots, \mathbf{x}_p)\) como variables explicativas o predictoras se debe a uno de dos propósitos:

  • Modelo para explicar las relaciones: busca describir y cuantificar explícitamente la relación entre la variable de resultado \(y\) y un conjunto de variables explicativas \(\mathbf{X}\), así como determinar la importancia de cualquier relación.

  • Modelo para la predicción: busca predecir una variable de resultado \(y\) basado en la información contenida en un conjunto de variables predictivas \(\mathbf{X}\). Acá no necesariamente importa comprender cómo se relacionan e interactúan todas las variables entre sí, sólo lograr buenas predicciones sobre \(y\) utilizando la información en \(\mathbf{X}\).

¿Qué objetivo se persigue?

En cada caso indique si el objetivo del modelo debe ser explicativo o predictivo. Suponga que tenemos interés en identificar:

  1. ¿Cuáles son los factores de riesgo (como el hábito de fumar, la edad, etc) que se asocian con el cáncer de pulmón?.
  1. ¿Qué contenido que le gustaría ver con mayor probabilidad a un usuario de una plataforma digital?.
  1. ¿Cuál es el efecto de la edad de la mujer, nivel educativo de la pareja, hábitos de la pareja (fuma, toma, etc), composición del hogar (con hijos/sin hijos) sobre la violencia doméstica?.

Especificación del modelo

. . .

¿Es un buen modelo?

. . .

Análisis de los residuales

. . .

Pronósticos

. . .

Algoritmos

Algunos modelos son:

  • Lineales: lm().

  • Generalizados: glm().

  • Bayesianos: stan_glm()

  • Penalizados: glmnet()

  • ML: tidymodels

Formulación del modelo

Sea \(\mathcal{D}=\{(y_i, \mathbf{x}_i): i=1,\ldots,n\}\), con \(y_i\) la \(i\)-ésima respuesta medida en una escala continua; \(\mathbf{x}_i=(x_{i1},\ldots,x_{ip})^t \in \mathbb{R}^p\) es el vector de variables predictoras; y \(n\) \((\gg p)\) es el tamaño de la muestra. El modelo lineal se especifica así:


\[y_i = \beta_0 + \beta_1 x_{i1} + \cdots + \beta_p x_{ip} + \varepsilon_i \hspace{0.25cm} \text{con } \varepsilon_i \overset{\text{iid}}{\sim} \mathcal{N}(\mu, \sigma^2)\]

Supuestos del modelo



  • Linealidad: \(\mu = E(y|\mathbf{x}_i) = \mathbf{X}\mathbf{\beta}\)
  • Independencia: \(Cov(\varepsilon_i, \varepsilon_j) = 0\) para \(i\ne j\)
  • Homocedasticidad: \(V(\varepsilon_i|\mathbf{X}) = \sigma^2\)
  • Normalidad: \(\varepsilon_i \overset{\text{iid}}{\sim} \mathcal{N}(\mu, \sigma^2)\)

Estimación de los parámetros

  • Puntos: Valores observados
  • Recta: Valores ajustados

Estimación de los parámetros

\[Y_i = \hat{Y}_i + (Y_i - \hat{Y}_i) = \hat{Y}_i + e_i \]

El objetivo entonces es minimizar

\[\sum_{i=1}^n e_i^2 = \sum_{i=1}^n (Y_i - \hat{Y}_i)^2 = \sum_{i=1}^n (Y_i - [\beta_0 + \beta_1 X_i])^2\]

El procedimiento se conoce como Mínimos Cuadrados Ordinarios (MCO).

EJEMPLO

Variable explicativa cuantitativa

Considere los datos sobre los precios_vivienda

library(pacman)
p_load(tidyverse, psych, corrplot, broom, performance)

datos <- "https://github.com/jgbabativam/AnaDatos/raw/main/datos/precios_viviendas.rds"
load(url(datos))

modelo1 <- lm(precio ~ area_m2, data = precios_viviendas)
  • Escriba la ecuación del modelo
  • ¿Cuál debe ser la interpretación de los coeficientes de regresión?
  • ¿Cómo debe interpretarse el p-value?

Variable explicativa cuantitativa

library(gtsummary)
tbl_regression(modelo1, intercept = TRUE)
Characteristic Beta 95% CI p-value
(Intercept) 156 14, 297 0.031
area_m2 2.1 1.6, 2.6 <0.001
Abbreviation: CI = Confidence Interval

Visualización

ggplot(precios_viviendas, aes(x = area_m2, y = precio)) +
  geom_point() +
  labs(x = "Área en m2", y = "Precio (en millones)",
       title = "Relación entre el precio y el área") +  
  geom_smooth(method = "lm", formula = 'y ~ x') + 
  theme_classic()

¿Es bueno el modelo?

Descomposición de la varianza

\[\sum_{i=1}^n (Y_i - \bar{Y}_i)^2 = \sum_{i=1}^n (Y_i - \hat{Y})^2+ \sum_{i=1}^n e^2\]

Bondad del ajuste

El coeficiente de determinación es un indicador entre 0 y 1:

\[\sum_{i=1}^n (Y_i - \bar{Y}_i)^2 = \sum_{i=1}^n (Y_i - \hat{Y})^2+ \sum_{i=1}^n e^2\]

\[SCT = SCR + SCE\]

Se deduce que:

\[R^2 = \frac{SCR}{SCT} = 1 - \frac{SCE}{SCT}\]

El valor de \(R^2\) está entre 0 y 1.

Bondad del ajuste

tab <- modelo1 |> glance() 
gt::gt(tab)
r.squared adj.r.squared sigma statistic p.value df logLik AIC BIC deviance df.residual nobs
0.4451945 0.4394153 225.0639 77.03361 6.292077e-14 1 -668.8513 1343.703 1351.458 4862761 96 98
  • ¿cuál debe ser la interpretación del \(R^2\)?
  • ¿cómo debe interpretarse el p-value?
tbl_regression(modelo1, intercept = TRUE) |> 
  add_glance_table(include = c(r.squared, p.value))
Characteristic Beta 95% CI p-value
(Intercept) 156 14, 297 0.031
area_m2 2.1 1.6, 2.6 <0.001
0.445

p-value <0.001

Abbreviation: CI = Confidence Interval

Variable explicativa categórica

Podemos preguntarnos si el precio depende de la zona donde se ubica la vivenda. Para ello es clave que la variable categórica sea de clase factor.

Revisión de la clase de la variable zona

class(precios_viviendas$zona)
[1] "character"
precios_viviendas$zona <- as.factor(precios_viviendas$zona)

Análisis exploratorio de la asociación

precios_viviendas |> 
  ggplot(aes(x = zona, y = precio, fill = zona)) +
  geom_boxplot()+
  theme_bw()

Modelo de regresión


El ajuste del modelo no cambia:


levels(precios_viviendas$zona)
[1] "rural"  "urbana"
precios_viviendas$zona <-relevel(precios_viviendas$zona, ref = "rural") 

modelo2 <- lm(precio ~ zona, data = precios_viviendas)

Resultados

tbl_regression(modelo2, intercept = TRUE) |> 
  add_glance_table(include = c(r.squared, p.value))
Characteristic Beta 95% CI p-value
(Intercept) 689 608, 770 <0.001
zona


    rural
    urbana 125 6.3, 244 0.039
0.044

p-value 0.039

Abbreviation: CI = Confidence Interval
  • ¿Por qué no hay coeficiente para la zona rural?. Escriba la ecuación del modelo.
  • ¿Cómo se interpreta el coeficiente de la zona urbana y su valor p?
  • ¿Qué se puede decir acerca de la bondad del ajuste?

Punto 2 - Taller: DASS21

  1. Ajuste los modelos de regresión simple para interpretar las relaciones entre:
  • Depresión y Satisfacción con su trabajo
  • Ansiedad y Satisfacción con su trabajo
  • Estrés y Satisfacción con su trabajo
  1. Interprete los coeficientes de los modelos ajustados y discuta la bondad del ajuste.

Punto 3 - Taller: DASS21

Ajuste un modelo de regresión simple que le permita identificar si el puntaje de depresión se relaciona con el sexo

  1. Convierta la variable sexo en factor así: dass$sexo <- as_factor(dass$sexo)

  2. Ajuste el modelo de regresión y presente los resultados.

  3. Interprete los coeficientes y el valor p.

Variables explicativas mixtas

Ajuste un modelo para el ingreso en función de las variables de área m2, zona, número de habitaciones, cantidad de baños y año de construcción.

\[Precio_i = \beta_0 + \beta_1 \cdot area_i + \beta_2 \cdot zona_{urb} + \beta_3 \cdot habitaciones_i \\ + \beta_4 \cdot baños_i + \beta_5 \cdot añoc_i + \varepsilon_i\]

Observe que las variables explicativas son cuantitativas y cualitativas. Verifique que la clase esté bien definida.

modelo3 <-lm(precio ~ area_m2 + zona + habitaciones + banos + anio_construccion, 
             data = precios_viviendas)

Escriba la ecuación del modelo e interprete los resultados.

Resultados

tbl_regression(modelo3, intercept = TRUE) |> 
  add_glance_table(include = c(r.squared, p.value))
Characteristic Beta 95% CI p-value
(Intercept) -434 -5,329, 4,460 0.9
area_m2 1.7 1.1, 2.4 <0.001
zona


    rural
    urbana 113 25, 202 0.013
habitaciones -21 -97, 55 0.6
banos 70 -1.9, 142 0.056
anio_construccion 0.26 -2.3, 2.8 0.8
0.516

p-value <0.001

Abbreviation: CI = Confidence Interval

Análisis de los supuestos

Que no se cumplan los supuestos puede afectar varios aspectos: sesgos, problemas de pronóstico, error de contraste.

Análisis de los supuestos

par(mfrow = c(2, 2))
plot(modelo3)

Análisis de los supuestos

plot(modelo3, 1)

No se debe presentar un patrón, así que la línea roja debe estar aproximadamente de forma horizontal en cero.

plot(modelo3, 3)

Se espera que que los puntos queden igualmente distribuidos dentro de una banda estable.

plot(modelo3, 2)

Un ajuste cercano a la línea de 45 grados es indice de que el supuesto de normalidad se satisface.

plot(modelo3, 4)

Análisis de los supuestos

ajuste <- augment(modelo3)

ajuste |> 
  slice_max(.cooksd, n=3)
# A tibble: 3 × 12
  precio area_m2 zona   habitaciones banos anio_construccion .fitted .resid
   <dbl>   <dbl> <fct>         <dbl> <dbl>             <dbl>   <dbl>  <dbl>
1  2027.    604  urbana            3   4                1972   1458.   569.
2   127.    489. rural             4   4.5              1970   1158. -1031.
3   533.    477. rural             4   3                1962   1028.  -495.
# ℹ 4 more variables: .hat <dbl>, .sigma <dbl>, .cooksd <dbl>, .std.resid <dbl>

Revisar que no existan distancias mayores que 1.

Paquete performance

Revise los supuestos del modelo, verifique que VIF es menor que 10 (no hay multicolinealidad) y haga las pruebas de los supuestos


check_model(modelo3)
check_collinearity(modelo3)
check_outliers(modelo3)
check_autocorrelation(modelo3)
check_heteroscedasticity(modelo3)
check_normality(modelo3)

Punto 4 - Taller: DASS 21

  1. Ajuste un modelo de regresión lineal múltiple con al menos 3 variables explicativas que resulten significativas para modelar el puntaje de depresión. Escriba la ecuación, interprete los coeficientes, revise los supuestos y concluya.

  2. Ajuste un modelo de regresión lineal múltiple con al menos 3 variables explicativas que resulten significativas para modelar el puntaje de estrés. Escriba la ecuación, interprete los coeficientes, revise los supuestos y concluya.

TALLER

Reglas

  • Se deben realizar los 4 ejercicios propuestos a lo largo de la sesión.
  • El trabajo se debe realizar con los grupos que han conformado.
  • Si existen dos trabajos exactamente iguales, la nota será dividida entre el número de personas involucradas.

GRACIAS!

Referencias

  • Çetinkaya-Rundel, M. and Hardin, J. (2021) Introduction to modern statistics. Sections of Regression modeling: 7, 8, 9 y 10. Disponible aquí: https://openintro-ims.netlify.app/

  • Ismay, C., & Kim, A.Y. (2019). Statistical Inference via Data Science: A ModernDive into R and the Tidyverse (1st ed.). Chapman and Hall/CRC. https://doi.org/10.1201/9780367409913

  • Thompson, J. (2019). Tidy Data Science with the tidyverse and tidymodels. https://tidyds-2021.wjakethompson.com