7.3 Proceso de Benchmark

  1. Del censo extraer el total de personas por DAM2
total_pp <- readRDS(file = "Recursos/Día2/Sesion3/Data/total_personas_dam2.rds")

N_dam_pp <- total_pp %>%   group_by(region) %>%  
            mutate(region_pp = sum(pp_dam2) ) 

tba(N_dam_pp %>% slice(1:20))
region nombre_region dam dam2 id_municipio pp_dam2 region_pp
01 Región Cibao Norte 09 00901 010901 179829 1516957
01 Región Cibao Norte 09 00902 010902 6911 1516957
01 Región Cibao Norte 09 00903 010903 37378 1516957
01 Región Cibao Norte 09 00904 010904 7820 1516957
01 Región Cibao Norte 18 01801 011801 158756 1516957
01 Región Cibao Norte 18 01802 011802 18868 1516957
01 Región Cibao Norte 18 01803 011803 6333 1516957
01 Región Cibao Norte 18 01804 011804 22058 1516957
01 Región Cibao Norte 18 01805 011805 12639 1516957
01 Región Cibao Norte 18 01806 011806 16464 1516957
01 Región Cibao Norte 18 01807 011807 49593 1516957
01 Región Cibao Norte 18 01808 011808 17169 1516957
01 Región Cibao Norte 18 01809 011809 19717 1516957
01 Región Cibao Norte 25 02501 012501 691262 1516957
01 Región Cibao Norte 25 02502 012502 42092 1516957
01 Región Cibao Norte 25 02503 012503 16993 1516957
01 Región Cibao Norte 25 02504 012504 25539 1516957
01 Región Cibao Norte 25 02505 012505 38628 1516957
01 Región Cibao Norte 25 02506 012506 51695 1516957
01 Región Cibao Norte 25 02507 012507 37349 1516957
02 Región Cibao Sur 13 01301 021301 248089 710821
02 Región Cibao Sur 13 01302 021302 59052 710821
02 Región Cibao Sur 13 01303 021303 56803 710821
02 Región Cibao Sur 13 01304 021304 30261 710821
02 Región Cibao Sur 24 02401 022401 76554 710821
02 Región Cibao Sur 24 02402 022402 13759 710821
02 Región Cibao Sur 24 02403 022403 22117 710821
02 Región Cibao Sur 24 02404 022404 38962 710821
02 Región Cibao Sur 28 02801 022801 125338 710821
02 Región Cibao Sur 28 02802 022802 18952 710821
02 Región Cibao Sur 28 02803 022803 20934 710821
03 Región Cibao Nordeste 06 00601 030601 188118 624186
03 Región Cibao Nordeste 06 00602 030602 14062 624186
03 Región Cibao Nordeste 06 00603 030603 15709 624186
03 Región Cibao Nordeste 06 00604 030604 17864 624186
03 Región Cibao Nordeste 06 00605 030605 33663 624186
03 Región Cibao Nordeste 06 00606 030606 14661 624186
03 Región Cibao Nordeste 06 00607 030607 5497 624186
03 Región Cibao Nordeste 14 01401 031401 76993 624186
03 Región Cibao Nordeste 14 01402 031402 24524 624186
03 Región Cibao Nordeste 14 01403 031403 24240 624186
03 Región Cibao Nordeste 14 01404 031404 15168 624186
03 Región Cibao Nordeste 19 01901 031901 39557 624186
03 Región Cibao Nordeste 19 01902 031902 27765 624186
03 Región Cibao Nordeste 19 01903 031903 24871 624186
03 Región Cibao Nordeste 20 02001 032001 58156 624186
03 Región Cibao Nordeste 20 02002 032002 24509 624186
03 Región Cibao Nordeste 20 02003 032003 18829 624186
04 Región Cibao Noroeste 05 00501 040501 28071 394068
04 Región Cibao Noroeste 05 00502 040502 15624 394068
04 Región Cibao Noroeste 05 00503 040503 6951 394068
04 Región Cibao Noroeste 05 00504 040504 7274 394068
04 Región Cibao Noroeste 05 00505 040505 6035 394068
04 Región Cibao Noroeste 15 01501 041501 24644 394068
04 Región Cibao Noroeste 15 01502 041502 14921 394068
04 Región Cibao Noroeste 15 01503 041503 35923 394068
04 Región Cibao Noroeste 15 01504 041504 10559 394068
04 Región Cibao Noroeste 15 01505 041505 9136 394068
04 Región Cibao Noroeste 15 01506 041506 14424 394068
04 Región Cibao Noroeste 26 02601 042601 34540 394068
04 Región Cibao Noroeste 26 02602 042602 11183 394068
04 Región Cibao Noroeste 26 02603 042603 11753 394068
04 Región Cibao Noroeste 27 02701 042701 76863 394068
04 Región Cibao Noroeste 27 02702 042702 62205 394068
04 Región Cibao Noroeste 27 02703 042703 23962 394068
05 Región Valdesia 02 00201 050201 91345 1028129
05 Región Valdesia 02 00202 050202 11243 1028129
05 Región Valdesia 02 00203 050203 17620 1028129
05 Región Valdesia 02 00204 050204 20041 1028129
05 Región Valdesia 02 00205 050205 15257 1028129
05 Región Valdesia 02 00206 050206 19020 1028129
05 Región Valdesia 02 00207 050207 11235 1028129
05 Región Valdesia 02 00208 050208 17647 1028129
05 Región Valdesia 02 00209 050209 5263 1028129
05 Región Valdesia 02 00210 050210 5640 1028129
05 Región Valdesia 17 01701 051701 157316 1028129
05 Región Valdesia 17 01702 051702 27028 1028129
05 Región Valdesia 21 02101 052101 232769 1028129
05 Región Valdesia 21 02102 052102 15466 1028129
05 Región Valdesia 21 02103 052103 124193 1028129
05 Región Valdesia 21 02104 052104 31057 1028129
05 Región Valdesia 21 02105 052105 84312 1028129
05 Región Valdesia 21 02106 052106 42325 1028129
05 Región Valdesia 21 02107 052107 30268 1028129
05 Región Valdesia 21 02108 052108 9540 1028129
06 Región Enriquillo 03 00301 060301 36511 368594
06 Región Enriquillo 03 00302 060302 15702 368594
06 Región Enriquillo 03 00303 060303 26772 368594
06 Región Enriquillo 03 00304 060304 10619 368594
06 Región Enriquillo 03 00305 060305 7709 368594
06 Región Enriquillo 04 00401 060401 83619 368594
06 Región Enriquillo 04 00402 060402 14823 368594
06 Región Enriquillo 04 00403 060403 13164 368594
06 Región Enriquillo 04 00404 060404 15390 368594
06 Región Enriquillo 04 00405 060405 21605 368594
06 Región Enriquillo 04 00406 060406 3970 368594
06 Región Enriquillo 04 00407 060407 9112 368594
06 Región Enriquillo 04 00408 060408 8042 368594
06 Región Enriquillo 04 00409 060409 4703 368594
06 Región Enriquillo 04 00410 060410 8186 368594
06 Región Enriquillo 04 00411 060411 4491 368594
06 Región Enriquillo 10 01001 061001 16510 368594
06 Región Enriquillo 10 01002 061002 12029 368594
06 Región Enriquillo 10 01003 061003 8310 368594
06 Región Enriquillo 10 01004 061004 5668 368594
07 Región EL Valle 07 00701 070701 25924 295362
07 Región EL Valle 07 00702 070702 6533 295362
07 Región EL Valle 07 00703 070703 8344 295362
07 Región EL Valle 07 00704 070704 10587 295362
07 Región EL Valle 07 00705 070705 7281 295362
07 Región EL Valle 07 00706 070706 4360 295362
07 Región EL Valle 22 02201 072201 132177 295362
07 Región EL Valle 22 02202 072202 9685 295362
07 Región EL Valle 22 02203 072203 20843 295362
07 Región EL Valle 22 02204 072204 13062 295362
07 Región EL Valle 22 02205 072205 44163 295362
07 Región EL Valle 22 02206 072206 12403 295362
08 Región Yuma 08 00801 080801 66867 606323
08 Región Yuma 08 00802 080802 20813 606323
08 Región Yuma 11 01101 081101 251243 606323
08 Región Yuma 11 01102 081102 21967 606323
08 Región Yuma 12 01201 081201 139671 606323
08 Región Yuma 12 01202 081202 16558 606323
08 Región Yuma 12 01203 081203 89204 606323
09 Región Higuamo 23 02301 092301 195307 561431
09 Región Higuamo 23 02302 092302 22573 561431
09 Región Higuamo 23 02303 092303 8901 561431
09 Región Higuamo 23 02304 092304 30051 561431
09 Región Higuamo 23 02305 092305 19034 561431
09 Región Higuamo 23 02306 092306 14592 561431
09 Región Higuamo 29 02901 092901 46723 561431
09 Región Higuamo 29 02902 092902 31889 561431
09 Región Higuamo 29 02903 092903 31096 561431
09 Región Higuamo 29 02904 092904 55348 561431
09 Región Higuamo 29 02905 092905 20900 561431
09 Región Higuamo 30 03001 093001 61517 561431
09 Región Higuamo 30 03002 093002 16272 561431
09 Región Higuamo 30 03003 093003 7228 561431
10 Región Ozama 01 00101 100101 965040 3339410
10 Región Ozama 32 03201 103201 948885 3339410
10 Región Ozama 32 03202 103202 363321 3339410
10 Región Ozama 32 03203 103203 529390 3339410
10 Región Ozama 32 03204 103204 142019 3339410
10 Región Ozama 32 03205 103205 43963 3339410
10 Región Ozama 32 03206 103206 272776 3339410
10 Región Ozama 32 03207 103207 74016 3339410
  1. Obtener las estimaciones directa por región o el nivel de agregación en el cual la encuesta es representativa.

En este código, se lee un archivo RDS de una encuesta y se utilizan las funciones transmute() y paste0() para seleccionar y transformar las variables de interés.

En primer lugar, se crea una variable dam que corresponde al identificador de la división administrativa mayor de la encuesta. A continuación, se utiliza la columna dam_ee para crear una variable dam, se selecciona la variable dam2 que corresponde al identificador de la división administrativa municipal de segundo nivel (subdivisión del departamento) de la encuesta.

Luego, se crea una variable wkx que corresponde al peso de la observación en la encuesta, y una variable upm que corresponde al identificador del segmento muestral en la encuesta.

La variable estrato se crea utilizando la función paste0(), que concatena los valores de dam y area_ee (una variable que indica el área geográfica en la que se encuentra la vivienda de la encuesta).

Finalmente, se crea una variable pobreza que toma el valor 1 si el ingreso de la vivienda es menor que un umbral lp, y 0 en caso contrario.

encuesta <- readRDS("Recursos/Día2/Sesion3/Data/encuestaDOM21N1.rds") %>% 
  transmute(
    dam = haven::as_factor(dam_ee,levels = "values"),
    dam = str_pad(dam,width = 2,pad = "0"),
    dam2,
    wkx = `_fep`, 
    upm = `_upm`,
    estrato =`_estrato`,
    pobreza = ifelse(ingcorte < lp, 1 , 0)) %>% 
  inner_join(N_dam_pp %>% select(region,dam2) )

El código está realizando un análisis de datos de encuestas utilizando el paquete survey de R. Primero, se crea un objeto diseno de diseño de encuestas usando la función as_survey_design() del paquete srvyr, que incluye los identificadores de la unidad primaria de muestreo (upm), los pesos (wkx), las estratos (estrato) y los datos de la encuesta (encuesta). Posteriormente, se agrupa el objeto diseno por la variable “Agregado” y se calcula la media de la variable pobreza con un intervalo de confianza para toda la población utilizando la función survey_mean(). El resultado se guarda en el objeto directoDam y se muestra en una tabla.

library(survey)
library(srvyr)
options(survey.lonely.psu = "adjust")

diseno <-
  as_survey_design(
    ids = upm,
    weights = wkx,
    strata = estrato,
    nest = TRUE,
    .data = encuesta
  )

directoDam <- diseno %>% 
   group_by(region) %>% 
  summarise(
    theta_dir = survey_mean(pobreza, vartype = c("ci"))
    )
tba(directoDam)
region theta_dir theta_dir_low theta_dir_upp
01 0.1696 0.1472 0.1919
02 0.1618 0.1230 0.2007
03 0.1808 0.1371 0.2245
04 0.1714 0.1269 0.2160
05 0.2629 0.2287 0.2972
06 0.4221 0.3741 0.4702
07 0.2957 0.2316 0.3599
08 0.1998 0.1605 0.2392
09 0.2258 0.1817 0.2700
10 0.2409 0.2213 0.2604
  1. Realizar el consolidando información obtenida en 1 y 2.
temp <- estimacionesPre %>%
  inner_join(N_dam_pp %>% select(region,dam2,pp_dam2,region_pp) ) %>% 
  inner_join(directoDam )

tba(temp %>% slice(1:10))
dam2 theta_pred region pp_dam2 region_pp theta_dir theta_dir_low theta_dir_upp
00101 0.2204 10 965040 3339410 0.2409 0.2213 0.2604
00201 0.2234 05 91345 1028129 0.2629 0.2287 0.2972
00206 0.3097 05 19020 1028129 0.2629 0.2287 0.2972
00301 0.4158 06 36511 368594 0.4221 0.3741 0.4702
00302 0.4343 06 15702 368594 0.4221 0.3741 0.4702
00303 0.4579 06 26772 368594 0.4221 0.3741 0.4702
00304 0.4922 06 10619 368594 0.4221 0.3741 0.4702
00401 0.3379 06 83619 368594 0.4221 0.3741 0.4702
00402 0.2213 06 14823 368594 0.4221 0.3741 0.4702
00403 0.4380 06 13164 368594 0.4221 0.3741 0.4702
  1. Con la información organizada realizar el calculo de los pesos para el Benchmark
R_dam2 <- temp %>% group_by(region) %>% 
  summarise(
  R_dam_RB = unique(theta_dir) / sum((pp_dam2  / region_pp ) * theta_pred)
) 

tba(R_dam2)
region R_dam_RB
01 1.0526
02 0.9339
03 1.1303
04 1.0137
05 1.0035
06 1.0894
07 0.9750
08 1.0018
09 0.8851
10 0.9941

calculando los pesos para cada dominio.

pesos <- temp %>% 
  mutate(W_i = pp_dam2  / region_pp) %>% 
  select(dam2, W_i)
tba(pesos %>% slice(1:10))
dam2 W_i
00101 0.2890
00201 0.0888
00206 0.0185
00301 0.0991
00302 0.0426
00303 0.0726
00304 0.0288
00401 0.2269
00402 0.0402
00403 0.0357
  1. Realizar la estimación FH Benchmark

En este proceso, se realiza la adición de una nueva columna denominada R_dam_RB, que es obtenida a partir de un objeto denominado R_dam2. Posteriormente, se agrega una nueva columna denominada theta_pred_RBench, la cual es igual a la multiplicación de R_dam_RB y theta_pred. Finalmente, se hace un left_join con el dataframe pesos, y se seleccionan únicamente las columnas dam, dam2, W_i, theta_pred y theta_pred_RBench para ser presentadas en una tabla (tba) que muestra únicamente las primeras 10 filas.

estimacionesBench <- estimacionesPre %>% 
  inner_join(N_dam_pp %>% select(region,dam2) )%>%
  left_join(R_dam2, by = c("region")) %>%
  mutate(theta_pred_RBench = R_dam_RB * theta_pred) %>%
  left_join(pesos) %>% 
  select(region, dam2, W_i, theta_pred, theta_pred_RBench)  

  tba(estimacionesBench %>% slice(1:10))
region dam2 W_i theta_pred theta_pred_RBench
10 00101 0.2890 0.2204 0.2191
05 00201 0.0888 0.2234 0.2242
05 00206 0.0185 0.3097 0.3108
06 00301 0.0991 0.4158 0.4529
06 00302 0.0426 0.4343 0.4731
06 00303 0.0726 0.4579 0.4988
06 00304 0.0288 0.4922 0.5362
06 00401 0.2269 0.3379 0.3681
06 00402 0.0402 0.2213 0.2411
06 00403 0.0357 0.4380 0.4771
  1. Validación: Estimación FH con Benchmark
estimacionesBench %>% group_by(region) %>%
  summarise(theta_reg_RB = sum(W_i * theta_pred_RBench)) %>%
  left_join(directoDam, by = "region") %>% 
  tba()
region theta_reg_RB theta_dir theta_dir_low theta_dir_upp
01 0.1696 0.1696 0.1472 0.1919
02 0.1618 0.1618 0.1230 0.2007
03 0.1808 0.1808 0.1371 0.2245
04 0.1714 0.1714 0.1269 0.2160
05 0.2629 0.2629 0.2287 0.2972
06 0.4221 0.4221 0.3741 0.4702
07 0.2957 0.2957 0.2316 0.3599
08 0.1998 0.1998 0.1605 0.2392
09 0.2258 0.2258 0.1817 0.2700
10 0.2409 0.2409 0.2213 0.2604