Inferencia causal aplicada

1 Introduccion

La inferencia causal es una metodología estadística y analítica que busca identificar relaciones de causa y efecto entre variables, más allá de simples correlaciones. Su objetivo es responder preguntas del tipo “¿qué pasaría si…?”, permitiendo estimar el impacto de una intervención o cambio en una variable sobre otra. Esta técnica es ampliamente utilizada en campos como la economía, la medicina, las ciencias sociales, el marketing y el aprendizaje automático, donde es crucial entender cómo ciertas acciones o políticas influyen en los resultados observados. Por ejemplo, puede emplearse para evaluar el efecto de una política pública, el impacto de un tratamiento médico o el rendimiento de una estrategia comercial.

En este caso, asumiremos que a unos investigadores les interesa saber si el uso de mosquiteros disminuye el riesgo de que una persona contraiga malaria. Para ello han recopilado datos de 1.752 hogares en un país sin nombre y tienen variables relacionadas con los factores ambientales, salud individual y las características de los hogares. Además, este país tiene un programa especial que proporciona mosquiteros gratuitos a los hogares que cumplen con requisitos específicos: para calificar para el programa, debe haber más de 4 miembros del hogar, y el ingreso mensual del hogar debe ser inferior a 700 dólares al mes. Los hogares no se inscriben automáticamente en el programa, y muchos no lo usan. Los datos no son experimentales - los investigadores no tienen control sobre quién usa los mosquiteros, y los hogares individuales toman sus propias decisiones sobre si solicitan mosquiteros gratuitos o compran sus propios mosquiteros, así como si usan los mosquiteros si los tienen.

2 Importamos librerias

library(dagitty)   # Para crear, analizar y evaluar diagramas causales (DAGs)
library(ggdag)     # Extensión de ggplot2 para visualizar y personalizar DAGs
library(tidyverse) # Conjunto de herramientas para manipulación, transformación y análisis de datos
library(ggplot2)   # Para crear gráficos elegantes y personalizables en R
library(MatchIt)   # Para realizar emparejamiento de muestras (matching) en análisis causal
library(XReg) # Salida de modelos de regresion
library(texreg) # Resumen y comparaciones de modelos

3 Modelamos DAG (Diagram Aclyclic Graph)

mosquito_dag <- dagify(
  riesgo_malaria ~ red + ingreso + salud + temperatura + resistencia,
  red  ~ ingreso + salud + temperatura + elegible + hogar,
  elegible ~ ingreso + hogar,
  salud ~ ingreso,
  exposure = "red",
  outcome = "riesgo_malaria",
  coords = list(x = c(riesgo_malaria = 7, red  = 3, ingreso = 4, salud = 5,
                      temperatura = 6, resistencia = 8.5, elegible = 2, hogar = 1),
                y = c(riesgo_malaria = 2, red  = 2, ingreso = 3, salud = 1,
                      temperatura = 3, resistencia = 2, elegible = 3, hogar = 2)),
  labels = c(riesgo_malaria = "Riesgo de malaria", red  = "Mosquitero", ingreso = "Ingreso",
             salud = "Salud", temperatura = "temperaturas nocturnas", 
             resistencia = "Resistencia a insecticidas",
             elegible = "Eligible para el programa", hogar = "Miembros en el hogar")
)

ggdag_status(mosquito_dag, use_labels = "label", text = FALSE) + 
  theme_dag()

Este código en R construye y visualiza un grafo acíclico dirigido (DAG, por sus siglas en inglés) que representa relaciones causales hipotéticas entre variables relacionadas con el riesgo de malaria. Utiliza la función dagify() para definir cómo se cree que las variables están causalmente conectadas: por ejemplo, el uso de mosquiteros (red) depende de factores como ingreso, salud, temperatura, elegibilidad y tamaño del hogar, y a su vez influye en el riesgo de malaria (riesgo_malaria). Se especifican las variables de exposición (red) y resultado (riesgo_malaria), así como las coordenadas para su visualización y etiquetas descriptivas. Finalmente, ggdag_status() genera un gráfico que muestra el estado de las variables en el DAG, útil para analizar supuestos de identificación causal y planear estrategias de estimación.

Cada nodo en el DAG es una columna en la base de datos recogida por los investigadores, e incluye lo siguiente:

Riesgo de malaria (malaria_risk): La probabilidad de que alguien en el hogar se infecte con malaria. Se mide en una escala de 0-100, con valores más altos que indican un mayor riesgo.

Mosquitero ((net and net_num): Una variable binaria que indica si el hogar utiliza mosquiteros.

Elegible para el programa (eligible): Una variable binaria que indica si el hogar es elegible para el programa de mosquiteros gratuitos.

Ingresos (income): Los ingresos mensuales del hogar, en dólares americanos.

Temperaturas nocturnas (temperature): La temperatura media nocturna, en grados centígrados.

Salud (health): La salud de la familia se declarada por sus propios miembros. Medido en una escala de 0-100, con valores más altos que indican mejor salud.

Número en el hogar (household): Número de personas que viven en el hogar.

Resistencia a los insecticidas (resistance): Algunas cepas de mosquitos son más resistentes a los insecticidas y, por lo tanto, suponen un mayor riesgo de infectar a las personas con malaria. Esto se mide en una escala de 0-100, con valores más altos que indican una mayor resistencia.

Según el DAG, el riesgo de contraer malaria es causado por los ingresos, las temperaturas, la salud, la resistencia a los insecticidas y el uso de mosquiteros. Las personas que viven en regiones más cálidas, tienen ingresos más bajos, tienen peor salud, están rodeadas de mosquitos con alta resistencia a los insecticidas y los que no utilizan mosquiteros corren un riesgo mayor de contraer malaria que los que los utilizan. El uso de los mosquiteros se explica por los ingresos, las temperaturas nocturnas, la salud, el número de personas que viven en la vivienda y la elegibilidad para el programa de mosquiteros gratuitos. Las personas que viven en zonas más frescas por la noche, tienen ingresos más altos, mejor salud, tienen más personas en el hogar, y son elegibles para los mosquiteros gratuitos del gobierno, tienen más probabilidades de usar regularmente los mosquiteros. El DAG también muestra que la elegibilidad para el programa de mosquiteros gratuitos se explica por los ingresos y el tamaño de los núcleos familiares, ya que los hogares deben cumplir con umbrales específicos para calificar.

4 Cargamos la base de datos

load("~/Library/Mobile Documents/com~apple~CloudDocs/Cursos/Causal Inference/Causal Inference in R/R for political data science/data/mosquito_nets.rda")

glimpse(mosquito_nets)

Rows: 1,752
Columns: 10
$ id           <dbl> 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22, 23, 24, 25, 26, 27, 28,…
$ net          <lgl> TRUE, FALSE, FALSE, TRUE, FALSE, FALSE, TRUE, FALSE, FALSE, FALSE, FALSE, TRUE, TRUE, FALSE, FALSE, TR…
$ net_num      <dbl> 1, 0, 0, 1, 0, 0, 1, 0, 0, 0, 0, 1, 1, 0, 0, 1, 1, 0, 0, 1, 0, 0, 0, 1, 0, 0, 1, 1, 0, 0, 1, 0, 1, 0, …
$ malaria_risk <dbl> 33, 42, 80, 34, 44, 25, 19, 35, 32, 40, 30, 14, 31, 34, 22, 54, 23, 41, 37, 65, 59, 59, 41, 21, 25, 27…
$ income       <dbl> 781, 974, 502, 671, 728, 1050, 1146, 1093, 1037, 828, 889, 1133, 807, 950, 1195, 594, 936, 816, 962, 4…
$ health       <dbl> 56, 57, 15, 20, 17, 48, 65, 75, 60, 36, 75, 62, 42, 53, 90, 29, 47, 30, 55, 6, 36, 30, 45, 64, 62, 75,…
$ household    <dbl> 2, 4, 3, 5, 5, 1, 3, 5, 3, 3, 6, 3, 4, 3, 1, 5, 2, 1, 2, 6, 5, 4, 3, 7, 3, 4, 3, 3, 2, 2, 6, 5, 4, 4, …
$ eligible     <lgl> FALSE, FALSE, FALSE, TRUE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE…
$ temperature  <dbl> 21.1, 26.5, 25.6, 21.3, 19.2, 25.3, 27.4, 29.8, 27.6, 21.3, 17.4, 18.9, 27.6, 24.0, 24.5, 25.5, 26.1, …
$ resistance   <dbl> 59, 73, 65, 46, 54, 34, 45, 65, 55, 54, 33, 39, 37, 53, 34, 55, 32, 73, 54, 38, 51, 43, 53, 51, 33, 37…

5 Analisamos las independencias condicionales

El comando impliedConditionalIndependencies() se utiliza en el paquete dagitty en R para analizar un DAG (Gráfico Acíclico Dirigido).Este comando identifica las independencias condicionales que deberían existir en el modelo, basándose en su estructura. Teniendo en cuenta que el signo se interpreta como independencia ⊥.

Cuando dos variables son independientes incondicionalmente, significa que no tienen relación directa o indirecta, y conocer el valor de una no proporciona información sobre la otra. Por ejemplo, si A y B son incondicionalmente independientes, el conocimiento de A no afecta nuestra predicción sobre B. Por otro lado, dos variables son independientes condicionalmente si dejan de estar relacionadas al conocer el valor de una tercera variable C. En este caso, A y B pueden estar correlacionadas inicialmente, pero esa correlación desaparece al controlar por C.

Primero, verificar independecias condicionales

impliedConditionalIndependencies(mosquito_dag)

elgb _||_ rsst
elgb _||_ rsg_ | ingr, red, sald, tmpr
elgb _||_ sald | ingr
elgb _||_ tmpr
hogr _||_ ingr
hogr _||_ rsst
hogr _||_ rsg_ | ingr, red, sald, tmpr
hogr _||_ sald
hogr _||_ tmpr
ingr _||_ rsst
ingr _||_ tmpr
red _||_ rsst
rsst _||_ sald
rsst _||_ tmpr
sald _||_ tmpr

las independencias incondicionales (como elgb || rsst) significan que, basándonos en la estructura del DAG, las variables “Eligible para el programa” (elegible) y “Resistencia a insecticidas” (resistencia) no tienen una relación directa o indirecta; conocer el estado de una no cambia nada sobre la otra. Por otro lado, las independencias condicionales (como elgb || sald | ingr) indican que, aunque podría haber una relación entre “Elegible” y “Salud”, esta relación desaparece si conocemos “Ingreso”, que actúa como variable de ajuste.

Segundo, validamos correlación entre las variables

# Elegible ⊥ resistencia
cor(mosquito_nets$eligible,mosquito_nets$resistance)

[1] 0.01234798

# Elegible ⊥ temperatura
cor(mosquito_nets$eligible,mosquito_nets$temperature)

[1] -0.02092867

# Miembros Hogar ⊥ Ingreso
cor(mosquito_nets$household,mosquito_nets$income)

[1] 0.008606828

# Miembros Hogar ⊥ Resistencia
cor(mosquito_nets$household,mosquito_nets$resistance)

[1] 0.02036479

# Miembros Hogar ⊥ Salud
cor(mosquito_nets$household,mosquito_nets$health)

[1] 0.00009785337

# Miembros Hogar ⊥ temperatura
cor(mosquito_nets$household,mosquito_nets$temperature)

[1] -0.02355615

#  Ingreso ⊥ resistencias
cor(mosquito_nets$income,mosquito_nets$resistance)

[1] 0.01371297

#  Ingreso ⊥ temperatura
cor(mosquito_nets$income,mosquito_nets$temperature)

[1] 0.1145676

# Red ⊥ resistencia
cor(mosquito_nets$net_num,mosquito_nets$resistance)

[1] 0.01677202

# resistencia ⊥ salud
cor(mosquito_nets$resistance,mosquito_nets$health)

[1] 0.007841029

# resistencia ⊥ temperatura
cor(mosquito_nets$resistance,mosquito_nets$temperature)

[1] -0.02198524

# salud ⊥ temperatura
cor(mosquito_nets$health,mosquito_nets$temperature)

[1] 0.1199855

Podemos que observar que la correlación entre las variables independientes incondicionadas es significativamente bajas. Sugiriendo la ausencia de dependencias entre ellas

6 Evaluamos la independencia condicional

El riesgo de malaria debe ser independiente del número de miembros de la familia, dado que los niveles de salud, ingresos, uso de mosquiteros y temperaturas nocturnas son similares.

No podemos utilizar cor() para probar esta implicación, ya que hay muchas variables en juego, pero podemos utilizar un modelo de regresión para comprobar si el número de miembros del hogar está significativamente relacionado con el riesgo de malaria

# Riesgo de malaria ⊥ Miembros en el hogar| Salud, Ingreso, Uso de la red, Temperatura:
lm(malaria_risk ~ household + health + income + net_num  + temperature, 
   data = mosquito_nets) %>% 
  summary()

Call:
lm(formula = malaria_risk ~ household + health + income + net_num + 
    temperature, data = mosquito_nets)

Residuals:
     Min       1Q   Median       3Q      Max 
-13.1458  -3.9117  -0.5715   3.3512  16.4438 

Coefficients:
              Estimate Std. Error t value            Pr(>|t|)    
(Intercept)  76.206730   0.965788  78.906 <0.0000000000000002 ***
household    -0.015453   0.089311  -0.173               0.863    
health        0.148337   0.010672  13.900 <0.0000000000000002 ***
income       -0.075144   0.001036 -72.563 <0.0000000000000002 ***
net_num     -10.437014   0.266500 -39.163 <0.0000000000000002 ***
temperature   1.005790   0.030964  32.483 <0.0000000000000002 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 5.245 on 1746 degrees of freedom
Multiple R-squared:  0.8852,    Adjusted R-squared:  0.8848 
F-statistic:  2692 on 5 and 1746 DF,  p-value: < 0.00000000000000022

No es significativo (t= −0.17, p = 0.863), lo que significa que los dos son independientes, como se esperaba. Por lo tanto, la variable riesgos de malaria es independiente a la variable miembros del hogar.

Si hay correlaciones sustanciales y significativas entre los nodos que deberían ser independientes, es probable que haya un problema con la especificación del DAG.

Para este ejemplo en concreto hemos validado que no hay dependencia significativas y que nuestro DAG calcula la realizad del modelo.

7 Econtramos el set de ajuste

Existe una vía directa entre el uso de mosquiteros y el riesgo de contraer malaria, pero el efecto no se identifica causalmente debido a otras varias vías abiertas. Podemos dejar que R encuentre los sets de ajuste apropiados automáticamente:

# Comando para evualuar los caminos del DAG
paths(mosquito_dag)

$paths
 [1] "red -> riesgo_malaria"                                         
 [2] "red <- elegible <- ingreso -> riesgo_malaria"                  
 [3] "red <- elegible <- ingreso -> salud -> riesgo_malaria"         
 [4] "red <- hogar -> elegible <- ingreso -> riesgo_malaria"         
 [5] "red <- hogar -> elegible <- ingreso -> salud -> riesgo_malaria"
 [6] "red <- ingreso -> riesgo_malaria"                              
 [7] "red <- ingreso -> salud -> riesgo_malaria"                     
 [8] "red <- salud -> riesgo_malaria"                                
 [9] "red <- salud <- ingreso -> riesgo_malaria"                     
[10] "red <- temperatura -> riesgo_malaria"                          

$open
 [1]  TRUE  TRUE  TRUE FALSE FALSE  TRUE  TRUE  TRUE  TRUE  TRUE

$paths

Muestra todas las rutas (causales, indirectas, espurias) entre variables en el DAG. Esto es útil para visualizar cómo las variables están conectadas directa e indirectamente.

$open

La lista $open indica si las rutas están abiertas o cerradas. Las rutas [4], [5] están cerradas (FALSE), lo que sugiere que están bloqueadas por el control de ciertas variables. Para entender, la ruta 4 está cerrada por la variable elegibles, ya que no permite el flujo de ingreso (variable confusa) a la variable tratamiento a pesar de que si tiene influencia en el outcome.

Una ruta se determina abierta o cerrada dependiendo si permite el flujo causal. Es decir, una ruta es abierta (TRUE) si hay relación entre la exposición (tratamiento) y el resultado (outcome), ya sea por medio de una variable confusor o de forma directa como la ruta 1.

El siguente comando identifica las variables que deben ajustarse para obtener una estimación no sesgada del efecto causal. Las variables de ajuste son aquellas que “cierran” las rutas de confusión, permitiendo analizar el efecto directo entre las variables de interés.

El efecto de confusión ocurre cuando la relación observada entre dos variables está influenciada por una tercera variable (llamada confusor) que afecta tanto a la variable independiente (exposición) como a la variable dependiente (resultado). Este efecto puede distorsionar o sesgar la estimación del efecto causal, llevando a conclusiones incorrectas.

# Variables necesarias para controlar efecto de confusion
adjustmentSets(mosquito_dag)

{ ingreso, salud, temperatura }

Basándose en las relaciones entre todos los nodos del DAG, el ajuste por salud, ingresos y temperatura es suficiente para cerrar todas las puertas traseras e identificar la relación entre el uso de la red y el riesgo de malaria. Es importante destacar que no tenemos que preocuparnos por ninguno de los nodos relacionados con el programa gubernamental de redes gratuitas, ya que esos nodos no están d conectados con el riesgo de malaria. Sólo tenemos que preocuparnos por las relaciones de confusión.

Esto significa que, para estimar correctamente el efecto causal de la exposición (por ejemplo, red) sobre el resultado (por ejemplo, riesgo_malaria), es necesario controlar (o ajustar) por las variables ingreso, salud, y temperatura, ya que estas variables influyen tanto en el tratamiento (red) como en el resultado (malaria), bajo el termino en ingles confounder.

Podemos confirmar esto gráficamente con ggdag_adjustment_set():

ggdag_adjustment_set(mosquito_dag, shadow = TRUE,
                     use_labels = "label", text = FALSE)

8 Estimación ingenua no ajustada

Primero, podemos ver cuál es la relación entre el uso de mosquiteros y el riesgo de malaria en ausencia de cualquier ajuste. Si creamos un cuadro de la distribución del riesgo de malaria entre las personas que usan y no usan mosquiteros, vemos que el riesgo medio es sustancialmente menor entre los que usan mosquiteros

ggplot(mosquito_nets, aes(x = net, y = malaria_risk)) +
  geom_boxplot()

Podemos ejecutar un simple modelo de regresión para medir la diferencia promedio exacta:

model_naive <- lm(malaria_risk ~ net, data = mosquito_nets) 
tbl_regression(model_naive)

G1;H1;Errorh in tbl_regression(model_naive) : 
  could not find function "tbl_regression"
Error during wrapup: not that many frames on the stack
Error: no more error handlers available (recursive errors?); invoking 'abort' restart
g

summary(model_naive)

Call:
lm(formula = malaria_risk ~ net, data = mosquito_nets)

Residuals:
    Min      1Q  Median      3Q     Max 
-26.937  -9.605  -1.937   7.063  55.395 

Coefficients:
            Estimate Std. Error t value            Pr(>|t|)    
(Intercept)  41.9365     0.4049  103.57 <0.0000000000000002 ***
netTRUE     -16.3315     0.6495  -25.15 <0.0000000000000002 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 13.25 on 1750 degrees of freedom
Multiple R-squared:  0.2654,    Adjusted R-squared:  0.265 
F-statistic: 632.3 on 1 and 1750 DF,  p-value: < 0.00000000000000022

Según este modelo, parece que el uso de un mosquitero está asociado con una disminución de 16 puntos en el riesgo de malaria. Sin embargo, este no es el efecto causal. Este es un caso en el que la correlación no es igual a la causalidad. Otros factores como los ingresos, la salud y las temperaturas confunden la relación entre el uso de los mosquiteros y el riesgo

9 Regresiones

Una forma rápida y fácil de tratar de ajustar las variables de confusión es incluirlas como covariables en una regresión lineal. A primera vista, esto tiene sentido intuitivo: el objetivo del ajuste es comparar los nodos de tratamiento y de resultado con los mismos valores de los diversos confusores, y el propósito de la regresión múltiple es explicar la variación del resultado manteniendo constantes las diferentes variables explicativas. Sin embargo, hacer ajustes de confusión con la regresión lineal dará como resultado relaciones causales correctamente identificadas sólo en circunstancias muy específicas. Para que el ajuste basado en la regresión funcione, las relaciones entre todos los nodos de tratamiento, resultado y confusión deben ser lineales, lo cual es difícil de probar y verificar con datos de observación reales. Casi siempre es mejor utilizar una de las otras técnicas de ajuste que se describen a continuación (matching o ponderación de probabilidad inversa), ya que esos métodos no se basan en el supuesto de la linealidad.

model_regression <- lm(malaria_risk ~ net + income + temperature + health,
data = mosquito_nets)
summary(model_regression)

Call:
lm(formula = malaria_risk ~ net + income + temperature + health, 
    data = mosquito_nets)

Residuals:
    Min      1Q  Median      3Q     Max 
-13.143  -3.915  -0.561   3.333  16.461 

Coefficients:
              Estimate Std. Error t value            Pr(>|t|)    
(Intercept)  76.159901   0.926835   82.17 <0.0000000000000002 ***
netTRUE     -10.441932   0.264906  -39.42 <0.0000000000000002 ***
income       -0.075144   0.001035  -72.58 <0.0000000000000002 ***
temperature   1.005855   0.030953   32.50 <0.0000000000000002 ***
health        0.148362   0.010668   13.91 <0.0000000000000002 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 5.244 on 1747 degrees of freedom
Multiple R-squared:  0.8852,    Adjusted R-squared:  0.8849 
F-statistic:  3366 on 4 and 1747 DF,  p-value: < 0.00000000000000022

Según estos resultados, el uso de un mosquitero causa una disminución de 10.44 puntos en el riesgo de paludismo, en promedio. Nótese que debido a que hemos ajustado los confusores, ahora podemos usar justificadamente el lenguaje causal en lugar de hablar simplemente de asociaciones.

Podemos hablar de causalidad después de una regresión lineal cuando se han ajustado los confusores porque estos se incorporan explícitamente como variables adicionales en el modelo. Al incluirlos, se controla estadísticamente su influencia tanto sobre la variable independiente como sobre la dependiente, lo que permite “bloquear” los caminos espurios que podrían generar asociaciones engañosas. Este ajuste se realiza al estimar el efecto de la variable de interés mientras se mantienen constantes los confusores, lo que simula un escenario de intervención. Así, el coeficiente estimado refleja el efecto directo de la variable explicativa sobre la respuesta, justificando el uso de lenguaje causal.

10 Emparejamiento (Matching)

El principal problema de utilizar datos de observación en lugar de datos experimentales es que las personas que utilizaron un mosquitero lo hicieron sin ser asignadas a un grupo de tratamiento. Las características individuales llevaron a las personas a auto-seleccionarse para el tratamiento, lo que hace que las personas que usaron mosquiteros sean esencialmente diferentes de las que no lo hicieron. Gracias al DAG, conocemos muchos de los factores que causaron que la gente eligiera usar los mosquiteros: ingresos, salud y temperatura.

El objetivo del código es emparejar las observaciones del dataset mosquito_nets entre los grupos tratados (personas que recibieron una red antimosquitos, representado por net) y no tratados, basándose en las covariables income, temperature, y health. Este emparejamiento busca reducir el sesgo al balancear las covariables entre los dos grupos para estimar correctamente el efecto del tratamiento (net).

Nota: recordemos que una covariable es una característica o variable que puede influir tanto en el tratamiento como en el resultado que se quiere estudiar.

El emparejamiento en el código que emplea MatchIt se basa en la idea de crear grupos tratados y no tratados que sean comparables entre sí. Esto se logra al seleccionar observaciones de ambos grupos que sean similares en términos de las covariables que podrían confundir la relación entre el tratamiento (net) y el resultado. Este proceso emula un experimento controlado aleatorizado en un estudio observacional

Ejemplo de MatchIt: un ejemplo sencillo para ilustrar el emparejamiento con covariables:

Imagina que estás evaluando si entregar una red antimosquitos (tratamiento) mejora la salud de las personas (resultado). Sin embargo, las personas con mayor ingreso (income) o que viven en zonas con más mosquitos (temperature) podrían tener más probabilidades de recibir la red y también de tener mejor o peor salud, independientemente del tratamiento. Entonces, si comparas directamente a tratados y no tratados, podrías confundir el efecto real de la red con el efecto de esas otras variables.

Para evitar eso, el emparejamiento busca, por ejemplo, una persona tratada con ingreso de $500, temperatura de 28°C y salud regular, y la empareja con una persona no tratada con características muy similares. Así, al comparar sus resultados de salud, puedes atribuir con mayor confianza cualquier diferencia al efecto de la red, porque las otras condiciones están equilibradas. Esto es lo que se busca al emular un experimento aleatorizado en un estudio observacional.

Continuado con el analisis

La similitud se mide utilizando la distancia de Mahalanobis, que calcula la distancia multidimensional entre observaciones considerando las correlaciones y escalas de las covariables. Observaciones más cercanas en esta métrica son consideradas más similares.

# Indicamos semilla de reproducibilidad
set.seed(123)

# Indicamos modelo de emparejamiento
matched <- matchit(net ~ income + temperature + health, # Definir covariables
                   data = mosquito_nets, # Definir conjunto de datos
                   method = "nearest",  # Emperejamiento por vecinos mas cercanos
                   distance = "mahalanobis",  # Tipo de distancia multidimensional
                   replace = TRUE # permitir reutilizacion para emparejamiento
                   )

summary(matched)                  

Call:
matchit(formula = net ~ income + temperature + health, data = mosquito_nets, 
    method = "nearest", distance = "mahalanobis", replace = TRUE)

Summary of Balance for All Data:

Summary of Balance for Matched Data:

Sample Sizes:

Este código realiza un emparejamiento entre personas que recibieron una red antimosquitos y las que no, utilizando el paquete MatchIt en R. Lo hace comparando las covariables income, temperature y health, que podrían influir tanto en la asignación del tratamiento como en el resultado. Usa el método de vecinos más cercanos (method = “nearest”) con la distancia de Mahalanobis, que considera la correlación entre variables, para encontrar pares lo más similares posible. Además, permite que una misma persona no tratada sea reutilizada en varios emparejamientos (replace = TRUE), lo que mejora la calidad del emparejamiento cuando hay pocos controles disponibles. El objetivo es crear grupos comparables y reducir el sesgo para estimar el efecto causal del tratamiento.

En el resultado de Sample Sizestenemos que:

• Había 1071 controles y 681 tratados antes del emparejamiento.

• Se utilizaron 436 observaciones del grupo control para emparejar a las 681 del grupo tratado.

• Dado que replace = TRUE, algunas observaciones del grupo control se reutilizaron para emparejar a varios tratados.

• 635 observaciones del grupo control no pudieron ser emparejadas. Discarded:

• No se descartaron observaciones en este análisis.

El valor Matched (ESS), o Tamaño Efectivo de la Muestra (por sus siglas en ingles), representa una estimación ajustada del número de unidades de control emparejadas, teniendo en cuenta que algunas de ellas pueden haber sido reutilizadas varias veces durante el emparejamiento (porque se usó replace = TRUE). Aunque el número real de controles emparejados fue 436, el ESS fue de 320.94, lo que indica que la reutilización de controles reduce la diversidad efectiva del grupo de comparación, y por tanto, el poder estadístico del análisis. Es una forma más precisa de reflejar cuánta “información independiente” se tiene realmente en el grupo de control emparejado

El emparejamiento fue exitoso, logrando balancear las covariables entre los dos grupos, lo que permitirá estimar el efecto del tratamiento (net) con menor sesgo.

Ahora,podemos crear una nueva base de datos basada en este matching con la función match.data(). Observa cómo ahora hay sólo 1117 filas en lugar de 1752, ya que hemos descartado las observaciones no emparejadas. También ahora hay una nueva columna llamada weights. La función matchit() asigna a los pares de observaciones emparejadas diferentes pesos dependiendo de lo cerca o lejos que estén los emparejamientos en un intento de controlar la variación de la distancia. Podemos utilizar estos pesos en nuestro modelo de regresión para mejorar nuestra estimación del efecto causal

mosquito_nets_matched <- match.data(matched) 
glimpse(mosquito_nets_matched) 

Rows: 1,117
Columns: 11
$ id           <dbl> 1, 3, 4, 5, 6, 7, 8, 11, 12, 13, 14, 16, 17, 20, 22, 23, 24, 25, 27, 28, 31, 33, 35, 36, 40, 43, 44, 4…
$ net          <lgl> TRUE, FALSE, TRUE, FALSE, FALSE, TRUE, FALSE, FALSE, TRUE, TRUE, FALSE, TRUE, TRUE, TRUE, FALSE, FALSE…
$ net_num      <dbl> 1, 0, 1, 0, 0, 1, 0, 0, 1, 1, 0, 1, 1, 1, 0, 0, 1, 0, 1, 1, 1, 1, 0, 1, 1, 1, 1, 1, 0, 1, 1, 1, 0, 0, …
$ malaria_risk <dbl> 33, 80, 34, 44, 25, 19, 35, 30, 14, 31, 34, 54, 23, 65, 59, 41, 21, 25, 23, 26, 50, 29, 66, 27, 21, 21…
$ income       <dbl> 781, 502, 671, 728, 1050, 1146, 1093, 889, 1133, 807, 950, 594, 936, 431, 655, 896, 1033, 1115, 967, 8…
$ health       <dbl> 56, 15, 20, 17, 48, 65, 75, 75, 62, 42, 53, 29, 47, 6, 30, 45, 64, 62, 71, 37, 37, 53, 22, 67, 40, 53,…
$ household    <dbl> 2, 3, 5, 5, 1, 3, 5, 6, 3, 4, 3, 5, 2, 6, 4, 3, 7, 3, 3, 3, 6, 4, 1, 3, 2, 3, 2, 2, 2, 2, 5, 1, 1, 2, …
$ eligible     <lgl> FALSE, FALSE, TRUE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, TRUE, FALSE, TRUE, FALSE, …
$ temperature  <dbl> 21.1, 25.6, 21.3, 19.2, 25.3, 27.4, 29.8, 17.4, 18.9, 27.6, 24.0, 25.5, 26.1, 26.8, 23.6, 25.9, 23.4, …
$ resistance   <dbl> 59, 65, 46, 54, 34, 45, 65, 33, 39, 37, 53, 55, 32, 38, 43, 53, 51, 33, 52, 71, 59, 52, 84, 31, 63, 53…
$ weights      <dbl> 1.0000000, 0.6402349, 1.0000000, 0.6402349, 1.2804699, 1.0000000, 1.9207048, 0.6402349, 1.0000000, 1.0…

Por ultimo,podemos hacer una regresión usando los datos emparejados:

# Indicamos modelo lineal
model_matched <- lm(malaria_risk ~ net, data = mosquito_nets_matched,
weights = weights)

# Vemos resumen estadistico
summary(model_matched)

Call:
lm(formula = malaria_risk ~ net, data = mosquito_nets_matched, 
    weights = weights)

Weighted Residuals:
    Min      1Q  Median      3Q     Max 
-38.295  -7.605  -2.365   5.395  55.395 

Coefficients:
            Estimate Std. Error t value            Pr(>|t|)    
(Intercept)  36.0896     0.5971   60.45 <0.0000000000000002 ***
netTRUE     -10.4846     0.7647  -13.71 <0.0000000000000002 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 12.47 on 1115 degrees of freedom
Multiple R-squared:  0.1443,    Adjusted R-squared:  0.1435 
F-statistic:   188 on 1 and 1115 DF,  p-value: < 0.00000000000000022

Según estos resultados, el uso de un mosquitero causa una disminución de 10.48 puntos en el riesgo de malaria, en promedio. Una vez más, podemos utilizar el lenguaje causal ahora porque hemos ajustado los confusores al hacer la comparación, identificando así el camino causal entre el uso de mosquiteros y el riesgo de malaria.

11 Ponderación de probabilidad inversa

Una desventaja del uso de matching es que desechamos mucha información: las observaciones no coincidentes se descartan, y el tamaño de nuestra muestra puede reducirse significativamente. El emparejamiento también tiende a ser muy riguroso, ya que cada observación tratada debe ser emparejada con una (o más) observaciones no tratadas. Aunque que algunas observaciones no tratadas están en realidad muy cerca de las observaciones tratadas, aún así se descartan porque fueron superadas por observaciones que tienen una distancia ligeramente menor.

En lugar de desechar datos potencialmente útiles, podemos utilizar otros métodos para crear coincidencias que sean menos rigurosas pero más informativas. Un método común en epidemiología y bioestadística es la ponderación de probabilidad inversa (IPW, en inglés). En la IPW, a cada observación se le asigna un peso basado en lo bien que su asignación real al tratamiento coincide con la probabilidad prevista del tratamiento, y esos pesos se utilizan luego en un modelo de regresión para estimar el efecto causal del tratamiento en el resultado.

Primero, podemos ajustar los factores de confusión de nuestro mosquitero y el ejemplo de riesgo de malaria usando una ponderación de probabilidad inversa.

Para ello, usamos una regresión logística para predecir la propensión a usar un mosquitero usando los ingresos, la temperatura y la salud. Luego usamos los puntajes de propensión para calcular los pesos de probabilidad inversa usando esta fórmula:

\[ \frac{\text{Mosquitero}}{\text{Propensión}} + \frac{1 - \text{Mosquitero}}{1 - \text{Propensión}} \]

# Creamos regresion logistica
model_mosquito_net <- glm(net ~ income + temperature + health,
                          data = mosquito_nets, family = binomial(link = "logit"))



# Agregamos ipw
mosquito_nets_ipw <- mosquito_nets %>% # Indicamos el df
  mutate(propensity = predict(model_mosquito_net, type = "response"), #Aplicamos el modelo logistico, con las probabilidades
         ipw = (net_num / propensity) + (1 - net_num) / (1 - propensity)) # creamos columna ipw


# Vemos resultados
print(mosquito_nets_ipw)

Este código ajusta un modelo de regresión logística para estimar la probabilidad de que una persona reciba un mosquitero (net) en función de sus ingresos, temperatura y salud, y luego utiliza esas probabilidades (propensiones) para calcular pesos de ponderación por la inversa de la probabilidad de tratamiento (IPW). Estos pesos permiten corregir posibles sesgos de selección al comparar grupos tratados y no tratados, asignando mayor peso a las observaciones que recibieron un tratamiento inesperado según sus características. Así, se busca simular un escenario más parecido a un experimento aleatorio y estimar de forma más justa el efecto del tratamiento.

Segundo, ahora que tenemos pesos de probabilidad inversa, podemos usarlos en una regresión:

# Indicamos modelo lineal
model_ipw <- lm(malaria_risk ~ net, data = mosquito_nets_ipw,
weights = ipw)

# Vemos resultados
summary(model_ipw)

Call:
lm(formula = malaria_risk ~ net, data = mosquito_nets_ipw, weights = ipw)

Weighted Residuals:
    Min      1Q  Median      3Q     Max 
-45.705 -14.622  -4.924  10.791 162.642 

Coefficients:
            Estimate Std. Error t value            Pr(>|t|)    
(Intercept)  39.6788     0.4684   84.71 <0.0000000000000002 ***
netTRUE     -10.1312     0.6583  -15.39 <0.0000000000000002 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 19.54 on 1750 degrees of freedom
Multiple R-squared:  0.1192,    Adjusted R-squared:  0.1187 
F-statistic: 236.8 on 1 and 1750 DF,  p-value: < 0.00000000000000022

Según este modelo de IPW, el uso de un mosquitero causa una disminución de 10.13 puntos en el riesgo de malaria, en promedio. Una vez más, podemos utilizar con seguridad el lenguaje causal porque hemos identificado el camino causal entre el uso de mosquiteros y el riesgo de malaria teniendo en cuenta los factores de confusión en los pesos de probabilidad inversos.

12 Comparando todos los métodos

Ahora que hemos ejecutado varios modelos de regresión que se ajustan a los factores de confusión de diferentes maneras, podemos comparar los resultados todos juntos. La ingenua estimación de -16 parece ser definitivamente una sobreestimación: después de ajustar con la regresión, el matching y la ponderación de probabilidad inversa, el efecto causal del uso de un mosquitero en el riesgo de malaria es consistentemente de alrededor de -10. Asumiendo que nuestro DAG está correcto, encontramos con éxito un efecto causal a partir de datos no experimentales y observacionales.

texreg::screenreg(list(model_naive, model_regression,
                       model_matched, model_ipw),
                  custom.model.names = c("Naive", "Regression",
                                         "Matching", "IPW"))

===============================================================
             Naive        Regression   Matching     IPW        
---------------------------------------------------------------
(Intercept)    41.94 ***    76.16 ***    36.09 ***    39.68 ***
               (0.40)       (0.93)       (0.60)       (0.47)   
netTRUE       -16.33 ***   -10.44 ***   -10.48 ***   -10.13 ***
               (0.65)       (0.26)       (0.76)       (0.66)   
income                      -0.08 ***                          
                            (0.00)                             
temperature                  1.01 ***                          
                            (0.03)                             
health                       0.15 ***                          
                            (0.01)                             
---------------------------------------------------------------
R^2             0.27         0.89         0.14         0.12    
Adj. R^2        0.27         0.88         0.14         0.12    
Num. obs.    1752         1752         1117         1752       
===============================================================
*** p < 0.001; ** p < 0.01; * p < 0.05