10.2 Función Generalizada de Varianza

La Función Generalizada de Varianza (FGV) es una técnica estadística utilizada para suavizar las estimaciones de las varianzas directas de los estimadores. Esta técnica busca estimar la varianza suavizada del estimador directo a través de un modelo log-lineal que involucra un vector de covariables auxiliares. La GVF es particularmente útil para modelar las varianzas de los estimadores directos, ya que permite lidiar con la naturaleza positiva de este parámetro. Además, esta técnica ha sido ampliamente utilizada en la literatura para estimar la varianza de los estimadores directos en diferentes contextos, incluyendo la estimación de ingreso per-cápita en los Estados Unidos, cifras oficiales del mercado de trabajo en Canadá y las tasas de pobreza comunal en la región. En este sentido, la GVF se plantea en términos de una relación log-lineal con un vector de covariables auxiliares que puede variar dependiendo del contexto en que se aplique.

El proceso continua con la selección de las dam que posean una varianza estimada mayor que cero, un deff mayor que 1 y 2 o más UPMs. Para los dominios que superan estas condiciones se realiza la transformación \(\log(\hat{\sigma}^2_d)\), además se realiza la selección de las columnas identificador del municipio (id_dominio), la estimación directa del indicador (Rd), El número de personas en el dominio (n) y la varianza estimada del para la estimación directa Rd_var,siendo esta la que transforma mediante la función log().

indicador_dom1 <- indicador_dom %>% 
  filter(Rd_var>0 & Rd_deff>=1 & n_upm >= 2) 

baseFGV <-  indicador_dom1 %>%  
  dplyr::select(dam2 , Rd, n, Rd_var) %>%
  mutate(ln_sigma2 = log(Rd_var))

10.2.1 Gráficas exploratorias

El código muestra la creación de cuatro gráficos usando la librería ggplot2 y el uso de los datos baseFGV. Estos gráficos tienen como objetivo explorar la relación entre el logaritmo de la varianza y diferentes transformaciones de la n y Rd.

El primer gráfico (p1) representa la relación entre la estimación directa y el logaritmo de la varianza. El segundo gráfico (p2) representa la relación entre el tamaño de muestra y el logaritmo de la varianza. El tercer gráfico (p3) representa la relación entre \(n_d \times Rd\) y el logaritmo de la varianza. Finalmente, el cuarto gráfico (p4) representa la relación entre la raíz cuadrada de la estimación directa y el logaritmo de la varianza.

library(patchwork)
p1 <- ggplot(baseFGV, aes(x = Rd, y = ln_sigma2)) +
  geom_point() +
  geom_smooth(method = "loess") +
  xlab("Formal")

p2 <- ggplot(baseFGV, aes(x = n, y = ln_sigma2)) + 
  geom_point() +
  geom_smooth(method = "loess") + 
  xlab("Tamaño de muestra")

p3 <- ggplot(baseFGV, 
             aes(x = Rd * n, y = ln_sigma2)) + 
  geom_point() +
  geom_smooth(method = "loess") + 
  xlab("Número de Formales")

p4 <- ggplot(baseFGV, 
             aes(x = sqrt(Rd), y = ln_sigma2)) + 
  geom_point() +
  geom_smooth(method = "loess") + 
  xlab("Raiz cuadrada de tasa de formalidad")


(p1 | p2) / (p3 | p4)

rm('p1','p2','p3','p4')

10.2.2 Ajustando el modelo log-lineal de la varianza

El código ajusta un modelo de regresión lineal múltiple (utilizando la función lm()), donde ln_sigma2 es la variable respuesta y las variables predictoras son Rd, n, y varias transformaciones de éstas. El objetivo de este modelo es estimar la función generalizada de varianza (FGV) para los dominios observados.

library(gtsummary)
FGV1 <- lm(ln_sigma2 ~ 1 + Rd + 
             n + I(n ^ 2) + I(Rd * n) +
             I(sqrt(Rd)) + I(sqrt(n)) + 
             I(sqrt(Rd * n)) ,
           data = baseFGV)

tbl_regression(FGV1) %>% 
  add_glance_table(include = c(r.squared, adj.r.squared))
Characteristic Beta 95% CI1 p-value
Rd -22 -50, 5.8 0.12
n -0.01 -0.03, 0.01 0.3
I(n^2) 0.00 0.00, 0.00 0.5
I(Rd * n) 0.02 -0.01, 0.05 0.2
I(sqrt(Rd)) 45 -4.7, 94 0.076
I(sqrt(n)) 0.87 -0.48, 2.2 0.2
I(sqrt(Rd * n)) -1.4 -3.1, 0.34 0.11
0.683
Adjusted R² 0.656
1 CI = Confidence Interval

Después de tener la estimación del modelo se debe obtener el valor de la constante \(\Delta\) para lo cual se usa el siguiente código.

delta.hat = sum(baseFGV$Rd_var) / sum(exp(fitted.values(FGV1)))

De donde se obtiene que \(\Delta = 1.2364739\). Final es posible obtener la varianza suavizada ejecutando el siguiente comando.

baseFGV <-
  baseFGV %>% mutate(hat_var = delta.hat * exp(fitted.values(FGV1)))

10.2.3 Validaciones sobre el modelo

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

varianza suavizada Vs varianza estimada

ggplot(baseFGV, 
       aes(x = Rd_var, y = hat_var)) + 
  geom_point() +
  geom_smooth(method = "loess")

Este código está realizando una Consolidación de los dominios observados y no observados para lo cual hace una unión izquierda (left_join()) entre: indicador_dom y baseFGV de la cual selecciona las columnas de id_dominio y hat_var. El argumento by = id_dominio especifica que la unión debe realizarse mediante la columna id_dominio.

Luego, se utiliza la función mutate() para crear dos nuevas variables. La primera variable Rd_var se asigna el valor de Rd_var de baseFGV si hat_var no es un valor nulo (NA), de lo contrario se le asigna un valor NA_real_ (NA pero de tipo numérico). De manera similar, se crea la variable Rd_deff con el valor de Rd_deff de baseFGV si hat_var no es nulo, de lo contrario se le asigna un valor NA_real_.

base_sae <- left_join(indicador_dom,
                      baseFGV %>% select(dam2, hat_var), 
                      by = "dam2") %>%
  mutate(
    Rd_var = ifelse(is.na(hat_var), NA_real_, Rd_var),
    Rd_deff = ifelse(is.na(hat_var), NA_real_, Rd_deff)
  )

Ahora, se debe estimar deff_FGV y n_eff_FGV a parir de la varianza suvizada (hat_var).

base_FH <- base_sae %>%
  mutate(
    Rd_deff = ifelse(is.nan(Rd_deff), 1, Rd_deff),
    deff_FGV = ifelse(Rd_var == 0 ,
      1,
      hat_var / (Rd_var / Rd_deff) #Fórmula del nuevo DEFF
    ),
   # Criterio MDS para regularizar el DeffFGV
    deff_FGV = ifelse(deff_FGV <= 1, NA_real_, deff_FGV), #Deff estimado
    n_eff_FGV = n / deff_FGV, #Número efectivo de personas encuestadas
   # Si no se estimó varianza para ese municipio, también excluir
   # la estimación directa de este municipio, esto es relevante para el modelo FH  
    hat_var = ifelse(deff_FGV <= 1, NA_real_, hat_var), 
    Rd = ifelse(is.na(hat_var), NA_real_, Rd) 
  )
tba(head(base_FH %>% select(dam2,n,n_upm,Rd, Rd_var,hat_var:n_eff_FGV), 10))
dam2 n n_upm Rd Rd_var hat_var deff_FGV n_eff_FGV
0101 2951 126 0.4147 0.0005 0.0005 5.9862 492.9681
3201 2840 108 0.4233 0.0003 0.0005 6.0488 469.5114
2501 3057 87 0.4108 0.0004 0.0005 6.2443 489.5639
3203 1944 59 0.4858 0.0004 0.0006 5.0816 382.5577
3202 1046 42 0.4221 0.0006 0.0011 4.5844 228.1662
1101 1198 38 0.3788 0.0012 0.0009 4.7384 252.8274
3206 836 32 0.3968 0.0006 0.0015 5.3642 155.8470
0901 743 20 0.5236 0.0028 0.0013 3.9091 190.0700
1301 738 20 0.4899 0.0011 0.0014 4.2704 172.8174
2101 505 20 0.4522 0.0015 0.0025 5.0969 99.0789

10.2.4 Otras validaciones sobre el resultado del modelo.

Continuando con el proceso de validación se construye el siguiente gráfico de dispersión con la variable de la varianza del estimador directo en el eje y y la varianza FGV en el eje x, para los municipios que tienen valores válidos para ambas variables. La línea de regresión lineal se ajusta a los puntos usando el método de mínimos cuadrados.

La visualización del gráfico permite evaluar si la FGV está capturando adecuadamente la variabilidad de la variable de interés (en este caso, la variable de varianza del estimador directo). Si la FGV captura la variabilidad, se espera que los puntos estén relativamente cerca de la línea de regresión, lo que indicaría que la FGV explica una gran parte de la variabilidad de la varianza del estimador directo. Por otro lado, si la FGV no captura la variabilidad, los puntos estarán más dispersos y alejados de la línea de regresión.

nDom <- sum(!is.na(base_FH$hat_var))
temp_FH <- base_FH %>% filter(!is.na(hat_var))

ggplot(temp_FH %>% arrange(n), aes(x = hat_var, y = Rd_var)) + 
  geom_point() + 
  geom_smooth(method = "lm", col = 2) + 
  labs(x = "FGV", y = "VaRdirEst") +
  ylab("Varianza del Estimador Directo")

Ahora, se realiza la comparación de la variabilidad de la varianza del estimador directo frente a la varianza suavizada a medida que el tamaño de muestra aumenta. El eje x representa el tamaño de la muestra y el eje y representa las varianzas. La línea azul representa la varianza FGV, mientras que la línea roja representa la varianza del estimador directo. En el gráfica es posible notar que la varianza FGV tiene una menos volatilidad que la varianza directa.

ggplot(temp_FH %>% 
         arrange(n), aes(x = 1:nDom)) +
  geom_line(aes(y = Rd_var, color = "VarDirEst")) +
  geom_line(aes(y = hat_var, color = "FGV")) +
  labs(y = "Varianzas", x = "Tamaño muestral", color = " ") +
  scale_x_continuous(breaks = seq(1, nDom, by = 10),
                     labels = temp_FH$n[order(temp_FH$n)][seq(1, nDom, by = 10)]) +
  scale_color_manual(values = c("FGV" = "Blue", "VarDirEst" = "Red"))

Siguiendo en la misma línea, se realiza la comparación del efectivo directo (n_eff_DIR) y el efectivo FGV (n_eff_DIR). El código que se muestra a continuación produce un gráfico que compara el tamaño de muestra efectivo obtenido a través de la estimación del DEFF con el tamaño de muestra directo. En el eje x se muestra el tamaño de muestra directo (n) y en el eje y se muestra el tamaño de muestra efectivo, calculado a través de la fórmula n/DEFF para la estimación directa (en rojo) y para la FGV (en azul).

Se puede observar que, en general, el tamaño de muestra efectivo estimado a través de la FGV es menos variable que el estimado a través de la estimación directa, lo que indica que la FGV reduce la varianza de la estimación. Además, se puede observar que para algunos dominios, el tamaño de muestra efectivo estimado a través de la FGV es menor que el tamaño de muestra directo, lo que podría deberse a la estimación de la varianza a través de la FGV. En general, este gráfico es útil para comparar la eficiencia de la estimación a través de la FGV y la estimación directa para cada dominio.

ggplot(temp_FH %>%
         arrange(n), aes(x = 1:nDom)) +
  geom_line(aes(y =  n / Rd_deff, color = "n_eff_DIR")) +
  geom_line(aes(y = n_eff_FGV, color = "n_eff_FGV")) +
  labs(y = "Tamaño de muestra efectivo",
       x = "Tamaño muestral", color = " ") +
  scale_x_continuous(breaks = seq(1, nDom, by = 10),
                     labels = temp_FH$n[order(temp_FH$n)][seq(1, nDom, by = 10)]) +
  scale_color_manual(values = c("n_eff_FGV" = "Blue", "n_eff_DIR" = "red"))

Por último, guardamos la base resultante.

saveRDS(object = base_FH, "Recursos/Día3/Sesion3/Data/base_FH.Rds")