Raking, Introdución al
I.
Ni que decirse tiene que a partir de las probabilidades conjuntas pueden construirse las marginales: se integra (o suma) y ya.
II.
El problema inverso es irresoluble: es imposible reconstruir las conjuntas a partir de las marginales. Las conjuntas, condicionadas a las marginales, pueden tener muchos grados de libertad.
Sin embargo, a petición de los usuarios finales, los comerciales de la estadística se han comprometido históricamente a resolver ese problema de manera científica. Así que los curritos de la estadística, supongo que muy a su pesar, han tenido que desarrollar cosas como las cópulas —esas sí que son verdaderas weapons of math destruction— y el raking, que es lo que nos ocupa hoy.
III.
Cuando uno usa cópulas puede definir ciertos parámetros, i.e., integrar información adicional a la presente en las mismas marginales, a la conjunta resultante y tener mucha fe. Nunca he visto, de hecho, una aplicación sólida de las cópulas y sí una muy calamitosa: la conocemos como crisis financiera de 2008 (véase esto).
El llamado raking, por su parte, no tiene, no parece tener, o quienes hablan de él no le dan importancia a, las hipótesis adicionales a partir de las cuales intentan reconstruir las conjuntas a partir de las marginales. (Probablemente, la culpa sea mía por informarme sobre el asunto en páginas web publicadas por vaya uno a saber qué diletante en lugar de consultar libros sesudos; no obstante, en lo que sigue, voy a usar código pro.)
IV.
Voy a comenzar, a diferencia de otros, que usan datos reales en los que la verdad no es conocida, datos artificiales de los que sabremos todo:
library(weights)
library(mvtnorm)
library(plyr)
library(anesrake)
n_pop <- 1e5
# random covariance matrix
m <- matrix(runif(9) - .3, 3, 3)
my_cov <- m %*% t(m)
# random mvn
muestra <- rmvnorm(n_pop, rep(0, 3), my_cov)
pop <- data.frame(
id = 1:n_pop,
sexo = cut(
muestra[,1],
breaks = c(-Inf, quantile(muestra[,1], .52), Inf)),
edad = cut(
muestra[,2],
breaks = c(-Inf, quantile(muestra[,2], c(.35, .72)), Inf)),
zona = cut(
muestra[,3],
breaks = c(-Inf, quantile(muestra[,3], c(.20, .40, .60)), Inf))
)
levels(pop$sexo) <- c("h", "m")
levels(pop$edad) <- c("25", "50", "75")
levels(pop$zona) <- c("z1", "z2", "z3", "z4")
conjunta_pop <- ddply(
pop,
.(sexo, edad, zona),
summarize,
prop_real = length(sexo) / nrow(pop))
marginales_pop <- list(
sexo = wpct(pop$sexo, pop$prop_real),
edad = wpct(pop$edad, pop$prop_real),
zona = wpct(pop$zona, pop$prop_real)
)
En el código anterior he creado una población con n_pop
individuos categorizados de acuerdo con las variables sexo
, edad
y zona
y he calculado las correspondientes proporciones conjuntas y marginales.
Con
muestra <- ddply(pop, .(sexo, edad, zona), head, 100)
muestra$id <- 1:nrow(muestra)
conjunta_muestra <- ddply(muestra, .(sexo, edad, zona),
summarize,
prop_muestra = length(sexo) / nrow(muestra))
marginales_muestra <- list(
sexo = wpct(muestra$sexo, muestra$prop_real),
edad = wpct(muestra$edad, muestra$prop_real),
zona = wpct(muestra$zona, muestra$prop_real)
)
creo una muestra de esa población. Se trata de una muestra de unos 100 sujetos por cada combinación de sexo
, edad
y zona
. Luego, calculo las correspondientes conjuntas y marginales de la muestra.
Obviamente, van a ser distintas. En la población, por diseño, las tres variables están correlacionadas y en la muestra son (aproximadamente) independientes. De hecho, la muestra está prácticamente equilibrada.
V.
La promesa del raking consiste en unos pesos que asignar a los sujetos de la muestra de tal forma que las proporciones ponderadas marginales coincidan con las de la población.
Esto se hace, obviamente, cuando:
- Las distribuciones conjuntas de la población son desconocidas.
- Las distribuciones marginales de la población son conocidas.
- Se espera que al asignar pesos a las observaciones de la muestra, al equipararse las distribuciones marginales reponderadas a las de la población, las distribuciones conjuntas reponderadas también acaben pareciéndose a las de la población.
Esto último es un pequeño acto de fe que no puede comprobarse experimentalente siempre con datos reales, aunque sí con datos simulados como los que he construido hoy aquí.
Efectivamente, si uno hace
raking <- anesrake(
marginales_pop,
muestra,
muestra$id,
cap = 5,
choosemethod = "total",
type = "pctlim",
pctlim = 0.05
)
muestra$pesos_raking <- raking$weightvec
seudoconjunta <- ddply(
muestra,
.(sexo, edad, zona),
summarize,
prop_raking = sum(pesos_raking) / sum(muestra$pesos_raking))
marginales_raking <- list(
sexo = wpct(muestra$sexo, muestra$pesos_raking),
edad = wpct(muestra$edad, muestra$pesos_raking),
zona = wpct(muestra$zona, muestra$pesos_raking)
)
verá cómo las marginales_pop
y marginales_raking
coinciden. De hecho,
marginales_pop
# $sexo
# h m
# 0.52 0.48
#
# $edad
# 25 50 75
# 0.35 0.37 0.28
#
# $zona
# z1 z2 z3 z4
# 0.2 0.2 0.2 0.4
y
marginales_raking
# $sexo
# h m
# 0.5 0.5
#
# $edad
# 25 50 75
# 0.35 0.37 0.28
#
# $zona
# z1 z2 z3 z4
# 0.2 0.2 0.2 0.4
VI.
Pero la gran pregunta es: ¿se parecerán las conjuntas? El lector interesado podrá responder a la pregunta por sí mismo ejecutando
res <- merge(conjunta_pop, seudoconjunta)
res <- merge(res, conjunta_muestra)
para fabricar una tabla, res
, en la que tendrá a su disposición las tres proporciones conjuntas: las de la población, las de la muestra y las ponderada vía raking.
VII.
He invertido 30 minutos de mi vida en el raking y no emitir juicios con pretensiones de ningún tipo. La gente que lleva 30 años con eso, entiendo, sabrá muy bien dónde, cuándo y en qué variante aplicar esta técnica para alcanzar resultados que no sean cuestionables como los que he presentado hoy. Pero los invito a divulgar más y divulgar mejor sobre esta técnica. No vaya a ser que caiga en manos de una chavalada optimista que nos arme otra marimorena.