Aquí quedó pendiente hablar de datos y métodos. Los primeros proceden de El Mundo. Solicité a Marta Ley, una coautora, los datos pero, antes de que contestase que sí (¡gracias!), me di cuenta de que podía obtenerlos solito: basta con capturar la llamada que el javascript local hace al servidor.
¿Métodos? Mejorables: se suaviza la intención de voto (con loess) y se estima la diferencia con un modelo de efectos mixtos, i.e.,
modelo<- lmer(delta ~ 1 + (1 | medio), data = misdatos)
¿Caveats? Veo dos: el primero, que loess suaviza teniendo en cuenta también observaciones futuras. Los autores de las encuestas no ven la verdad: solo los resultados de las encuestas previas. Debería haber usado como referencia la mejor predicción basada en observaciones pasadas. El segundo, que los porcentajes de los distintos partidos suman un total. Los sesgos no son independientes y yo los modelo como tales.
Y termino con el código completo:
library(rjson) library(plyr) raw <- fromJSON(readLines("https://spreadsheets.google.com/feeds/list/1vyVTJPr7ZpuQI4m17cekWl485cQ-Zh6O9Yb6zXkPpYI/od6/public/values?alt=json")) dat <- raw$feed$entry res <- ldply(dat, unlist) res[, "id.$t"] <- res[, "updated.$t"] <- NULL res$category.scheme <- res$category.term <- res$title.type <- NULL res$`content.$t` <- res$link.href <- res$link.rel <- res$content.type <- NULL res[, "title.$t"] <- res$link.type <- NULL colnames(res) <- make.names(colnames(res)) res$gsx.casa..t <- NULL res$fecha <- as.Date(res$gsx.fechaok..t, format = "%d/%m/%Y") res$medio <- res$gsx.empresaymedio..t res$margen.error <- as.numeric(gsub(",", ".", res$gsx.margendeerror..t)) res$tamano <- as.numeric(gsub("\\.", "", res$gsx.tamañomuestra..t)) res <- res[res$tamano < 1e6,] hist(res$tamano) res$int.pp <- as.numeric(gsub(",", ".", res$gsx.pp..t)) res$int.psoe <- as.numeric(gsub(",", ".", res$gsx.psoe..t)) res$int.cs <- as.numeric(gsub(",", ".", res$gsx.cs..t)) res$int.podemos <- as.numeric(gsub(",", ".", res$gsx.podemos..t)) res$int.iu <- as.numeric(gsub(",", ".", res$gsx.iu..t)) res <- res[, -grep("^gsx", colnames(res))] library(ggplot2) library(reshape2) tmp <- melt(res, id.vars = c("fecha", "medio", "margen.error", "tamano")) ggplot(tmp, aes(x = fecha, y = value)) + geom_smooth() + geom_point() + facet_wrap(~ variable) library(lme4) library(lattice) ## pp tmp <- res tmp$pred.pp <- predict(loess(int.pp ~ as.numeric(fecha), data = res)) tmp$delta.pp <- tmp$int.pp - tmp$pred.pp mod.pp <- lmer(delta.pp ~ 1 + (1 | medio), data = tmp) png("/tmp/sesgo_encuestas_pp.png", width = 600, height = 500) dotplot(ranef(mod.pp, condVar = TRUE)) dev.off() ## psoe tmp <- res tmp$pred.psoe <- predict(loess(int.psoe ~ as.numeric(fecha), data = res)) tmp$delta.psoe <- tmp$int.psoe - tmp$pred.psoe mod.psoe <- lmer(delta.psoe ~ 1 + (1 | medio), data = tmp) png("/tmp/sesgo_encuestas_psoe.png", width = 600, height = 500) dotplot(ranef(mod.psoe, condVar = TRUE)) dev.off() ## podemos tmp <- res tmp$pred.podemos <- predict(loess(int.podemos ~ as.numeric(fecha), data = res)) tmp$delta.podemos <- tmp$int.podemos - tmp$pred.podemos mod.podemos <- lmer(delta.podemos ~ 1 + (1 | medio), data = tmp) png("/tmp/sesgo_encuestas_podemos.png", width = 600, height = 500) dotplot(ranef(mod.podemos, condVar = TRUE)) dev.off() ## ciudadanos tmp <- res tmp$pred.cs <- predict(loess(int.cs ~ as.numeric(fecha), data = res)) tmp$delta.cs <- tmp$int.cs - tmp$pred.cs mod.cs <- lmer(delta.cs ~ 1 + (1 | medio), data = tmp) png("/tmp/sesgo_encuestas_ciudadanos.png", width = 600, height = 500) dotplot(ranef(mod.cs, condVar = TRUE)) dev.off() ## iu tmp <- res tmp$pred.iu <- predict(loess(int.iu ~ as.numeric(fecha), data = res)) tmp$delta.iu <- tmp$int.iu - tmp$pred.iu mod.iu <- lmer(delta.iu ~ 1 + (1 | medio), data = tmp) png("/tmp/sesgo_encuestas_iu.png", width = 600, height = 500) dotplot(ranef(mod.iu, condVar = TRUE)) dev.off()
Como siempre, interesante. Ya barruntaba yo que había un lmer por algún lado en cuánto vi los dotplots del anterior post.
Saludos.
Buenísima entrada! Tengo muchísimas dudas, pero, por no apabullar, soltaré la que más me corroe: ¿es válido utilizar este tipo de modelo cuando en alguno de los niveles, como por ejemplo, “TNS Demoscopia” o “20 Minutos ( A+M)”, sólo tenemos una observación?
Saludos!
Sí, tiene sentido. Es por lo que uso modelos mixtos. Esencialmente, lo que sucede es que hay una “transferencia de información” sobre la variabilidad de las estimaciones de aquellos medios con varias a aquellas que solo tienen una. De todos modos, la incertidumbre por la existencia de pocas observaciones se manifiesta en la magnitud del intervalo de confianza: a menos observaciones, menor precisión. ¡Incluso podríamos estimar la desviación de un medio del que no tuviésemos información! Obviamente, no sabríamos si es positiva o negativa, pero tendríamos una medida “razonable” de la desviación “a priori” (¡antes de ver datos!) de un medio.
Muchas gracias, Carlos.
Lo de “estimar la desviación de un medio del que no tuviésemos información” me ha dejado descolocado. Voy a buscar información para ahondar en esto de los modelos multinivel, porque, por lo que veo, mi conocimiento de estos modelos es muy limitado.
Saludos.
José Ignacio, un libro que me ayudó mucho a entender esto de los modelos multinivel/mixtos es el de Andrew Gelman. Es bastante asequible, con ejemplos en R (y también en bugs), y lo de la “estimación en categorías de las que no tenemos información” viene muy bien explicado. http://www.stat.columbia.edu/~gelman/arm/
Saludos
Un par de ejemplos que te puedan aclarar un poco. Uno es del libro de Gelman y otro es mío utilizando la EPA.
http://rpubs.com/joscani/mixed_model_example
http://rpubs.com/joscani/modmixtosepa
Muchas gracias por las referencias, Jose Luis.
Acabo de encargar el libro. Saludos.