Estimar la incertidumbre

Autores importantes — como Andrew Gelman — señalan que la estadística se enfrenta a tres retos fundamentales de inferencia:

Ya hemos hablado un poco del primer problema cuando discutimos la semántica de los datos. Y el tercero realmente es un problema de diseño de investigación que sólo en un segundo momento se convierte en un tema estadístico; para ello les sugiero leer a Judea Pearl, en especial The Book of Why y si quieren algo más estadístico, puede ser el de Hernán y Robbins, What if. Por ello este tema escapa a los alcances de la clase.

Así es que ahora nos vamos a concentrar en el segundo problema. En este contexto, la inferencia es la manera en que hacemos juicios sobre parámetros poblacionales a partir de una muestra.

Censos y Muestras

Cuando queremos hacerle una pregunta a la población entonces hay básicamente dos estrategias:

Si nosotros queremos una respuesta muy precisa entonces realizamos un censo. Sin embargo, esto sólo es factible si la población objetivo es muy pequeña (digamos, el salón de clases o la maestría) o si tenemos MUCHO dinero y MUCHO tiempo. Por lo general, cuando hablamos de poblaciones grandes, estamos hablando de muestras.

Alt text

Alt text

Las muestras son muy útiles porque al conocer la opinión de un grupo muy pequeño de la población, entonces podemos asumir que conocemos la opinión de la población. Esto significa que un parámetro poblacional (MIU) esperamos que coincida con una estimación muestral (Media).

La ventaja es que esto es más barato y más rápido que un censo. La desventaja es que hacer una muestra es un procedimiento técnico especializado y necesariamente introducimos incertidumbre en los resultados. En otras palabras, conocemos la estimación muestral pero no sabemos si le atinamos al parámetro poblacional…

Tipos de muestras

Por lo general las muestras se pueden realizar siguiendo dos estrategias con implicaciones cruciales:

  1. Probabilísticas: en las muestras probabilísticas todos los miembros de la población objetivo tienen la misma probabilidad de ser elegidos. En otras palabras, la selección es estrictamente aleatoria o al azar. La aleatoriedad suele ilustrarse con algún mecanismo, como las urnas en la lotería, o computacional con la generación de números aleatorios, como cuando asignamos exposiciones.

Recuerden que para toda muestra se requiere de un marco muestral, esto es, el listado a partir del cual se calcula la probabilidad de inclusión.

Hay tres estrategias para conformar muestras probabilísticas:

  1. Por Conveniencia: En ocasiones no es posible hacer una selección aleatoria. Si el mecanismo de aleatoriedad no es explícito: sospechen. Cuando el investigador u otro elige a quien se tiene disponible entonces se utiliza un método de conveniencia. Si bien las estimaciones provenientes de esta estrategia tienen valor exploratorio, este método de selección no permite generalizar a una población. ¿Qué tan sólidas son las encuestas en línea con millones de observaciones?

Por ello, siempre que se pueda, preferimos aleatorizar.

Una muestra involucra un proceso técnico especializado. No lo intenten en sus casas sin la supervisión de un experto. Dos tips generales y poco científicos para determinar el tamaño de muestra:

¿Qué tipo de muestra es la ENCOVID-19?

¿Qué piensan de que sea n=833?

Tipos de error muestral

Recuerden que todas las muestras tienen error, tienen incertidumbre en sus resultados. Esto es inevitable. El objetivo no es eliminar el error sino identificar sus fuentes y minimizarlo.

Hay dos grandes tipos de error con las muestras:

La recolección de datos debe de ser de alta calidad para que las estimaciones sean confiables

¿Qué tipo de errores sistemáticos puede tener la ENCOVID-19?

Vean la siguiente gráfica que presentamos al publicar la primera ENCOVID

Nos falló por 4%

¿Es aceptable?

Alt text

Alt text

Y ven ahora por edad:

¿Es aceptable la incertidumbre?

Alt text

Alt text

¿De qué se trata la columnade “calibrada”?

Tipos de estimaciones

¿Prefieren pescar con un harpón o con una red?

Cuando consumimos resultados basados en encuestas provenientes de una muestra, con frecuencia solemos ver 2 tipos de estimaciones. La diferencia se refiere al grado de incertidumbre de una aseveración.

El error muestral (+/- 3%) vuelve incierto afirmar que la estimación puntual coincida exactamente con el parámetro poblacional.

En cambio, una estimación de intervalo incluye un rango dentro del cual el parámetro poblacional se encuentra, dada cierta probabilidad o nivel de confianza.

En términos de comunicación es más fácil reportar prevalencias o comparaciones con estimaciones puntuales, pero no es preciso porque dan la falsa impresión de que se dice con certidumbre. Aunque sea más complejo, deberíamos de utilizar intervalos, ya que comunican con mayor precisión la incertidumbre.

Vean la siguiente gráfica que publicamos hace un par de meses

El 20% de octubre es la estimación puntual de octubre

¿Es menor que en junio?

Alt text

Alt text

Intervalos de confianza

El intervalo de confianza (IC) es un rango que expresa el margen de error con el cual se hace una inferencia con datos limitados; inferencia de una muestra a una población.

Veamos las siguientes tres aseveraciones. ¿Cuál tiene menor incertidumbre?

  1. 54% de los mexicanos aprueba la gestión de AMLO. n=1,000

  2. 14% reporta inseguridad alimentaria severa. n= 9,578

  3. 38% se atendió en una farmacia en el ultimo año. n= 250

Claramente la respuesta se asocia con el tamaño de muestra, de ahí la pista en cada renglón. Pero noten cómo la información periodística suele relegar el tamaño muestral, en el mejor de los casos, a una nota al pie. El más preciso es la número 2, con un IC de 13.32%- 14.71%. Le sigue la 1 con un IC de 51%-57%.Y vean la diferencia con el 3, con un IC de 32% - 44%.

¿Por qué son diferentes los IC entre abril y junio en la ENCOVID?

Intervalos de confianza y tamaño de muestra

El tamaño de los intervalos de confianza la decide el investigador por diseño. El tamaño del intervalo refleja el alpha (0.05). Y en este sentido significa la probabilidad con la cual el parámetro poblacional se encuentra dentro del rango.

95% de confianza= 5% de probabilidad de error

Una manera de entender esta proproción es suponer que, si la muestra se volviera a realizar 100 veces, en 95 de esas ocasiones el parámetro poblacional estaría dentro del intervalo de confianza estimado.

Ojo: el intervalo no es una probabilidad.

Vean el siguiente diagrama de Motulsky, en Intuitive Biosatistics. Imaginen que yo quiero saber la prevalencia de depresión. En el diagrama, cada rectángulo es una muestra. Y el final superior e inferior de los rectángulos son los intervalos. Como es un ejemplo teórico, la línea horizontal de referencia muestra el parámetro poblacional.

Alt text

Alt text

Lo que queremos es que nuestro IC (el rectángulo) capture el parámetro poblacional (línea horizontal en .25). Pero en una de cada 20 muestras nuestros resultados no lo van a capturar.

En otras palabras, si 20 universidades hubieran hecho la ENCOVID-19 en el mismo periodo y con la misma población, 19 de ellas reportarían un intervalo que captura el parámetro de la depresión (el promedio censal) y 1 de ellas habría fallado.

¿Si fuera Las Vegas apostarían a que su muestra, con su IC, capturó el parámetro de la depresión?

La mayor parte de las ciencias sociales sí hace la apuesta. Y este es el gran drama del muestreo. Nosotros sólo haremos una muestra, no 20, y nunca sabremos si ganamos o perdimos. Debemos de asumir que ganamos la apuesta y que la muestra que realizamos, efectivamente, captura el parámetro poblacional.

Ustedes se preguntarán, con justa razón: ¿qué tiene de especial el 95%? Yo Podría querer mayor certidumbre.

El 95% es una convención, es una tradición. Pero no es el único nivel de confianza.

Vean esta figura. Aquí las tres cajas quieren capturar el parámetro poblacional. Mientras más alto, mayor certidumbre porque la “red” es más grande. Pero mientras más alto es menos útil. Piensen en la Depresión. Si no tuviera alas tan anchas, llegaríamos a otras inferencias más útiles.

Alt text

Alt text

E imaginen sus perspectivas con un IC al 90%. Son menos atractivas. Si 10 universidades levantan la ENCOVID, 1 de ellas va a fallar. ¿Se la juegan?

Mejor vamos a una apuesta más segura. Si quieren IC cortos al 95%, lo que deben hacer es levantar muestras más grandes. Lo que suele representar más dinero. Pero incluso si son muy ricos, nunca van a llegar al 100% de certidumbre; a menos que hagan un censo.

Vean el diagrama siguiente. Todos los IC están al 95%. Pero esa definción de confianza, en una muestra mas grande, significa que el IC es más corto. Si quieren más precisión, suban la muestra. Y si son pobres y no pueden pagar una muestra grande… Pues su solución es abrir el IC de 95%.

¿Recuerdan la diferencia entre abril y junio en la ENCOVID?

Alt text

Alt text

La decisión del IC es una decisión contextual que depende de la naturaleza y capacidades del estudio. Si pueden, quedarse en 95% es una buena opcion. Y siempre buscamos las muestras más grandes posibles - hasta que el bolsillo aguante!

Muestreo aplicado

Vamos a imaginar esta situación. La Vicerrectora les pregunta cuál es el puntaje de preparatoria promedio de los alumnos de la IBERO. Ustedes tienen dos maneras de dar esta respuesta:

Primero hacen 20 llamadas a alumnos. Veamos qué hay dentro del objeto.

muestra20 <- telefonos_alumnos %>% 
  rep_sample_n(size = 20)

muestra20

Ahora veamos la media y la desviación estandard

muestra20 %>% 
  summarize(n = n(), media= mean(promedio), desv.std= sd(promedio))

Y graficamos. en el ceje de la Y, el count se refiere a número de personas. si la barra llega a 1, significa que una persona dio esa respuesta. si llega a 2, fueron dos personas. No pierdan de vista la escala de los counts porque pronto irá aumentando.

La línea roja representa la media de esta muestra

¿Cómo interpretamos la desviación estandard?

muestra20 %>%
  ggplot(aes(x= promedio))+
  geom_histogram(binwidth = 1, fill= "blue", color="white")+
  geom_vline(xintercept = mean(muestra20$promedio), color= "red")+
  theme_minimal()

¿Le entregamos este resultado a la Vicerrectora? ¿Están confiados de su resultado?

Este es el drama del muestreo. ¿Qué pasa si el vecino repite exactamente el mismo ejercicio?. ¿Obtendremos lo mismo?

Repitamos las 20 llamadas a alumnos y vemos las estadísticas: ¿coinciden las medias?

muestra20b <- telefonos_alumnos %>% 
  rep_sample_n(size = 20)

muestra20b %>% 
  summarize(n = n(), media= mean(promedio), desv.std= sd(promedio))

Veamos la gráfica:

muestra20b %>%
  ggplot(aes(x= promedio))+
  geom_histogram(binwidth = 1, fill= "blue", color="white")+
  geom_vline(xintercept = mean(muestra20$promedio), color= "red")+
  theme_minimal()

¿Entonces a cuál le creemos?

La mejor salida a este problema es repetir la muestra muchas veces. Vamos a suponer que tomamos a 50 alumnos de la maestría y les pedimos que cada uno haga 20 llamadas.

muestra20_50 <- telefonos_alumnos %>% 
  rep_sample_n(size = 20, reps = 50)

muestra20_50

Noten varias cosas. Primero, ya tenemos 1,000 filas, en lugar de 20.

La columna de “replicate” nos dice el número de alumno que entrevistó. Si le van dando a la flecha vemos los resultados que obtuvo cada encuestador.

Una posibilidad es obtener el promedio general, como si fuera una muestra grandota. Pero eso nos lleva casi al mismo lugar. Si otros encuestadores eligen otros 1,000 alumnos distintos, entonces el promedio sería diferente y no sabríamos cuál de los dos elegir.

Vamos a tomar otro camino. Nosotros queremos saber qué tan diferentes son los resultados entre los encuestadores.

medias20_50 <- muestra20_50 %>%
  group_by(replicate) %>% 
  summarize(n = n(), media= mean(promedio), desv.std= sd(promedio))

medias20_50

Veamos detenidamente estos resultados. El primer encuestador (replicate=1), hizo 20 llamdas, y la media de los promedios de prepa es… De igual manera, el segundo entrevistador, con 20 llamadas, obtuvo algo similar, pero no igual. También vemos que la distancia con respecto a la media es diferente.

Vamos a graficar estas medias para entender qué tanto difieren entre sí los encuestadores

medias20_50 %>%
  ggplot(aes(x= media))+
  geom_histogram(binwidth = 1, fill= "coral", color="white")+
  geom_vline(xintercept = mean(medias20_50$media), color= "black")+
  theme_minimal()

Es fundamental destacar varias diferencias entre esta gráfica y la previa.

La desviación estandard de las medias tiene un nombre diferente: error estandard Este es el número clave para saber, en promedio, qué tan diferentes son los resultados de nuestros encuestadores.

medias20_50 %>% 
  summarize(ErrorStd= sd(media), media= mean(media),  n = n())

Sabemos que la media de medias estará alrededor del 78-79 y que, en promedio, cada encuestador se separará por poquito más de 2 puntos de ese valor (un punto arriba y un punto abajo). Nos acercamos a poder calcular los intervalos de confianza.

Pero antes, vamos a entender un poco mejor qué sucede si cambiamos algunas decisiones.

muestra80_100 <- telefonos_alumnos %>% 
  rep_sample_n(size = 80, reps = 100)

muestra80_100

Tenemos ahora 8,000 filas.

muestra80_100 <- muestra80_100 %>%
  group_by(replicate) %>% 
  summarize(n = n(), media= mean(promedio), desv.std= sd(promedio))

muestra80_100

Ahí están las 80 llamadas por los 100 encuestadores

muestra80_100 %>%
  ggplot(aes(x= media))+
  geom_histogram(binwidth = 1, fill= "coral", color="white")+
  geom_vline(xintercept = mean(medias20_50$media), color= "black")+
  theme_minimal()

¿Qué diferencia vemos?

El rango se va haciendo cada vez más corto.

Pero volvamos a nuestra metáfora previa del 95% de confianza: De 100 muestra que realice (n=80), 5 de ellas van a fallar. Según esta gráfica, ¿Cuáles 5 dirían que es probable que estén lejos del parámetro poblacional?

Y ahora nuestros números clave

muestra80_100 %>% 
  summarize(ErrorStd= sd(media), media= mean(media),  n = n())

Varias cosas relevantes. La media es muy similar; no igual, pero similar. Al provenir de muestras más grandes, probablemente es más precisa. Ahora fíjense cómo el error estandar se redujo a la mitad. El error estandard es un indicador de la precisión de nuestro muestreo: más bajo, más preciso.

Vamos por otra ronda… ahora vamos a meter más encuestadores de maestría y que cada uno haga 1000 llamadas. ¿Les parecen 150 alumnos?

muestra1000_150 <- telefonos_alumnos %>% 
  rep_sample_n(size = 1000, reps = 150)

muestra1000_150

150,000 filas!!

Algunos me dirán que esto ya se descontroló. ¿Cómo es posible tener 150,000 respuestas si sólo hay 11,800 alumnos en la IBERO (población)? aquí debemos recordar que cada encuestador hace sus muestras de forma independiente y sin reemplazo. ¿Qué significa esto? Que si el encuestador 1 llamó al número de cuenta 95,999, no le avisó a nadie que ya tiene esta info. Y entonces el encuestador 2 volvió a llamar al número de cuenta 95,999. Tenemos duplicados. Pero esto no nos preocupa mucho porque aquí nos interesa saber qué tanto difieren los encuestadores entre sí. Así que seguimos.

muestra1000_150 <- muestra1000_150 %>%
  group_by(replicate) %>% 
  summarize(n = n(), media= mean(promedio), desv.std= sd(promedio))
muestra1000_150 %>%
  ggplot(aes(x= media))+
  geom_histogram(binwidth = 1, fill= "coral", color="white")+
  geom_vline(xintercept = mean(medias20_50$media), color= "black")+
  theme_minimal()

¿Cómo va nuestro rango?

muestra1000_150 %>% 
  summarize(ErrorStd= sd(media), media= mean(media),  n = n())

Otra vez tenemos una media muy parecida. Es decir, no importa cuántas veces hacemos llamadas aleatorias, la media de medias converge hacia el mismo valor. Y se reduce el error estandard.

Teorema del Límite Central

El hecho de que las distribuciones muestrales tengan una forma de campana tiene consecuencias cruciales. Esto permite extraer inferencias valiosas de su forma.

Analicemos las implicaciones en le siguiente diagrama, del Apéndice A.2 de ModernDive

Alt text

Alt text

En un momento quedarán más claras las implicaciones. Por ahora basta pensar que distribuciones más estrechas, más espigadas, permiten decir que el intervalo en el que se encuentra la proporción poblacional es menor. Mientras más ancha es la distribución muestral, el intervalo es más ancho.

Una muestra

Todo lo anterior ha sido un ejercicio mental. En realidad no vamos a conseguir tantos alumnos y no van a hacer tantas llamadas. Lo MÁS común es hacer una sola muestra. Tenemos sólo una oportunidad y debemos aprovecharla bien.

Quiero mostrarles que si sólo van a hacer una muestra, más les vale que sea grande. Veamos esto:

telefonos_alumnos %>% 
  rep_sample_n(size = 20) %>%
  summarize(n = n(), media= mean(promedio), desv.std= sd(promedio))
telefonos_alumnos %>% 
  rep_sample_n(size = 100) %>%
  summarize(n = n(), media= mean(promedio), desv.std= sd(promedio))
telefonos_alumnos %>% 
  rep_sample_n(size = 1000) %>%
  summarize(n = n(), media= mean(promedio), desv.std= sd(promedio))
telefonos_alumnos %>% 
  rep_sample_n(size = 9000) %>%
  summarize(n = n(), media= mean(promedio), desv.std= sd(promedio))

La moraleja de estos ejericios es esta: una mayor muestra, hace que la variación en los resultados sea menor. En otras palabras, mientras más grande es la muestra, tenemos una menor incertidumbre con los resultados. En otras palabras, estamos más seguros de que la media muestral coincide con el parámetro. Siguiendo el ejemplo: estamos más seguros que los resutados de nuestras medias coinciden con la media del Director de Servicios Escolares.

Error Estándard

Cada una de las gráficas anteriores (naranja) muestra las distribuciones muestrales porque se grafican los resultados obtenidos con muchas muestras. Estas distribuciones ilustran el concepto de variación muestral . Qué tanto cambian las medias de los encuestadores al hacer nuevas muestras. Con muestras grandes cambia poco (baja variación muestral) y por ello obtenemos estimaciones más precisas.

Si en lugar de verlo gráficamente lo queremos ver en un número, en un estadístico, requerimos de la desviación estándard de esas medias, es decir, del error estandard. el problema es que sólo tenemos una muestra.

Afortunadamente, hay una fórmula muy sencilla para ello. Asumamos que esta fue nuestra muestra real, nuestra única bala:

muestra_real <- telefonos_alumnos %>% 
  rep_sample_n(size = 830) 


tamanio_muestra <- 830
media_muestral <- mean(muestra_real$promedio)
desviacion_std_muestral <- sd(muestra_real$promedio)

muestra_real%>%
  summarize(n = n(), media= mean(promedio), desv.std= sd(promedio))

Esta es la fórmula del error estandard:

error_estandard <- desviacion_std_muestral / sqrt(tamanio_muestra)

error_estandard
## [1] 0.3329194

Ahora construimos el margen de error

margen_error <- 1.96*error_estandard

margen_error
## [1] 0.652522

Y recuperamos la media de nuestra muestra final

ci_inf <- media_muestral - margen_error

ci_inf
## [1] 78.17302
ci_sup <- media_muestral + margen_error

ci_sup
## [1] 79.47806

¿Le atinamos?

La Vicerrectora habló con el Director de Servicios Escolares y nos abrió su base de datos para ver qué tan bien lo hicimos.

Primero los estadísticos básicos: ¿cómo nos fue?

telefonos_alumnos %>%
  summarize(n = n(), media= mean(promedio), desv.std= sd(promedio))

Ahora en gráfica. La línea negra es la media poblacional: el parámetro que buscamos encontrar.

telefonos_alumnos %>%
  ggplot(aes(x= promedio))+
  geom_histogram(binwidth = 1, fill= "sea green", color="white")+
  geom_vline(xintercept = mean(telefonos_alumnos$promedio), color= "black") +
  theme_minimal()

Y comparemos con nuestros resultados. La línea roja es la media muestral: nuestra estimación que buscamos le “atine” al parámetro poblacional. Pero nosotros preferimos pescar con una red. Por eso hablamos de intervalos de confianza. Los límites inferior y superior se muestran con líneas naranja.

telefonos_alumnos %>%
  ggplot(aes(x= promedio))+
  geom_histogram(binwidth = 1, fill= "sea green", color="white")+
  geom_vline(xintercept = mean(telefonos_alumnos$promedio), color= "black") +
  geom_vline(xintercept = mean(muestra_real$promedio), color= "red") +
  geom_vline(xintercept = ci_inf, color= "orange") +
  geom_vline(xintercept = ci_sup, color= "orange") +
  theme_minimal()

¿Entonces podemos confiar en un muestreo aleatorio y de buen tamaño?

Construcción de intervalos de confianza

Hay dos métodos para construir los IC a un nivel de confianza al 95%. Ambos son práctimente equivalentes

  • Por Percentiles: A la distribución muestral se le quitan 2.5% en cada cola y se deja todo lo que está en medio. Esos valores en los extremos son el IC.

  • Por Error estándard. Aprovechar que la distribución normal, a 1.96 errores estandard, se obtiene el 95% de la distribuicón. La fórmula es:

  • IC inferior: proporción - (1.96 * error_estándard)

  • IC superior: proporción + (1.96 * error_estándard)

El error estándard se puede obtener por dos vías. Uno es al generar la distribución muestral por métodos computacionales y la otra es por medio de una fórmula (Ver conclusión de cap8 de ModernDive).

Para construirlos noten que en la realidad seguimos el proceso opuesto al de la sección previa. Antes, de muchas muestras obtuvimos el error estándard. Ahora el reto —mucho más realista— es construir el IC con la única muestra que tenemos.

IC de Ansiedad

Vamos a poner esto en práctica… Con nuestra proporción de ansiedad en la ENCOVID.

Utilicemos los resultados reportados como referencia

Alt text

Alt text

Ansiedad: 32.38% [29.25, 35.67]. El error es de 1.64

Vamos a intentar replicar esto

library(readr)
base_covid_abril2020_final <- read_csv("BD/base_covid_abril2020_final.csv")

abril <-base_covid_abril2020_final

abril <- abril %>% 
  rename(
  nervios= "p5_ansiedad",
  control= "p6_ansiedad")

abril <- abril %>% 
  mutate(gad2_total= nervios+ control)

abril <- abril %>% 
  mutate(gad2_level= if_else(gad2_total >= 3, 1, 0))

abril <- abril %>% 
  mutate(gad2_level2= ordered(gad2_level, levels=c("1", "0"), labels= c("Con síntomas", "Sin síntomas")))

abril2 <- abril %>% 
  select(gad2_level, gad2_level2)

prop.table(table(abril2$gad2_level2))*100
## 
## Con síntomas Sin síntomas 
##     32.65306     67.34694

Fíjense cómo la media no es igual; pero es cercano

IC de ansiedad con fórmula

el primer paso es calcular el error estándard

prop_abril <- abril2 %>% 
  summarize(porc_ans = mean(gad2_level))

se_abril <- prop_abril %>% 
  mutate(SE = sqrt(porc_ans * (1 - porc_ans) / 833)) # 833 es el tamaño de muestra de la encovid abril

se_abril

El error estandard es 0.016. Tampoco es exacto, pero muy parecido. Recuerden recorrer dos lugares el punto decimal.

Ahora a construir el margen de error. Este es el famoso +/- 3%

MdE_abril <- se_abril %>% 
  mutate(MdE= 1.96*SE)

MdE_abril

Y ahora sí el intervalo de confianza

ic_ans <- MdE_abril %>% 
  mutate(ci_inf = porc_ans - MdE,
         ci_sup = porc_ans + MdE)
ic_ans

Aquí está: proporción de Ansiedad = 32.65% [29.47, 35.84]

En un comando:

prop_abril2 <- abril2 %>% 
  summarize(porc_ans = mean(gad2_level)) %>%
  mutate(
        SE = sqrt(porc_ans * (1 - porc_ans) / 833),
        MdE= 1.96*SE,
        ci_inf = porc_ans - MdE,
        ci_sup = porc_ans + MdE)

prop_abril2

CI ansiedad con bootstrap

Vamos a intentar ahora el método de bootstrap largo

Hagamos 3,000 muestras de 500 personas. Noten que estamos partiendo de la ENCOVID. No son simulaciones

Primer paso

medias_remuestreo <- abril2 %>% 
  rep_sample_n(size = 500, replace = TRUE, reps = 3000)

medias_remuestreo

Vean cuántas filas tenemos: 1,500,000!

Ahora saquemos las medias para cada una de las muestras (“replicate”)

Sintetizo los dos pasos en uno

medias_remuestreo1 <- abril2 %>% 
  rep_sample_n(size = 500, replace = TRUE, reps = 3000) %>% 
  group_by(replicate) %>% 
  summarize(medias_ans = mean(gad2_level))

medias_remuestreo1

Vean la variabilidad en ansiedad. Tenemos 3,000 estimaciones puntuales de la media de anisedad.

Ahora estimemos la media de medias y su desviación estandard

medias_remuestreo1 %>% 
summarize(media = mean (medias_ans), sd = sd(medias_ans))

La media de medias es 32.72% Bastante cercano al oficial.

medias_remuestreo1 %>%
  ggplot(aes(x= medias_ans))+
  geom_histogram(binwidth = .01, fill= "coral", color="white")+
  geom_vline(xintercept = mean(medias_remuestreo1$medias_ans), color= "black")+
  theme_minimal()

Vamos a sacar los percentiles 2.5 y 97.5. Las colas muy improbables

quantile(medias_remuestreo1$medias_ans, c(.025, .975))
##  2.5% 97.5% 
## 0.288 0.368

Son un poco más anchos, pero no están mal. Si subimos los tamaños de las muestras pueden cerrarse

ggplot(medias_remuestreo1, aes(x = medias_ans)) +
  geom_histogram(binwidth = .01,fill= "coral", color = "white") +
  geom_vline(xintercept = .323, color= "black") + #oficial
  geom_vline(xintercept = c(.292, .357), color= "blue") + #oficial
  geom_vline(xintercept = c(.288, .366), color= "purple") + #estimado
  labs(x = "medias muestrales")+
  theme_minimal()

No se ve mal.

Prueben con otros parámetros. Sobretodo unos más pequeños para que se abra la distribución.

infer

Vamos a hacerlo con el paquete infer

medias_remuestreo_infer <- abril2 %>% 
  specify(response = gad2_level) %>% 
  generate(reps = 3000, type = "bootstrap") %>% 
  calculate(stat = "mean")

medias_remuestreo_infer
visualize(medias_remuestreo_infer)

percentile_ci <- medias_remuestreo_infer %>% 
  get_confidence_interval(level = 0.95, type = "percentile")

percentile_ci

Practicamente igual al oficial

visualize(medias_remuestreo_infer) + 
  shade_confidence_interval(endpoints = percentile_ci, color = "hotpink", fill = "khaki")

Ahora vamos a ver los intervalos con error estándard. Yo le asigno la media de la muestra original, la normalita de 32.6

se_ci <- medias_remuestreo_infer %>% 
  get_confidence_interval(level = 0.95, type = "se", point_estimate = 0.326) 

se_ci

MUY parecidos

visualize(medias_remuestreo_infer) + 
  shade_confidence_interval(endpoints = se_ci, color = "hotpink", fill = "khaki")

Ahora ya tenemos una manera de conocer las estimaciones de incertidumbre.

Vean esto cómo se publicó en Salud Pública de México; es un procedimiento similar, pero repetido mes a mes.

Alt text

Alt text

Rmarkdown

Vamos a ver una manera más interesante de utilizar RSTUDIO. Así es más fácil seguir las notas y el output

Empezamos la exploración de Rmarkdown por Notebooks y poco a poco iremos a cosas ligeramente más sofisticadas.

Antes que otra cosa, instala el paquete; a la derecha: Packages -> Install -> “Rmarkdown”

El motivo no sólo es estético o pedagógico (ya no tenemos que poner el #). El objetivo que perseguimos es la reproducibilidad de la ciencia.

La meta es que en un archivo tengan todo su texto y todo su análisis estadístico juntos.

Así ustedes mismos pueden volver al archivo original y hacer las modificaciones que quieran.

¿Qué nos dice el machote precargado?

This is an R Markdown Notebook. When you execute code within the notebook, the results appear beneath the code.

Notebook

Aquí pongan atención a que es un “Notebook”. Esta es la versión más sencilla de Rmarkdown. Nos permite hacer cosas que, de inicio, no pensamos publicar o compartir. Más adelante veremos que también podemos hacer html, .pdf,.doc y un fascinante etcétera.

En Rmarkdown tenemos tres tipos de “áreas”:

  1. YAML: Parte inicial, entre tres guiones. Ahí ponemos el título y metadatos del documento. Entre otras cosas anotamos ahí que queremos un notebook

  2. Texto: Es la sección en la que estamos ahora. Podemos escribir libremente como en un word en blanco. Es el equivalente a poner un # en el script tradicional.

  3. Código: estos son los comandos que le pedimos a R para ejecutar. Diferenciamos ambos tipos de texto con estos signos: ```{r} Noten que el pedazo de código termina con otras tres comillas. La r ahí dentro indica que usaremos el lenguaje de R; pero podríamos usar otros lenguajes, como Python o Julia. Dentro de esos corchetes podemos meter indicaciones a R, mismas que veremos poco a poco. Hay dos maneras de correr sólo ese “chunk” o pedazo de código:

Try executing this chunk by clicking the Run button within the chunk or by placing your cursor inside it and pressing Cmd+Shift+Enter.

plot(cars)

Noten cómo ahora el output sale inmediatamente después.

Este es el shortcut para insertar código:

Add a new chunk by clicking the Insert Chunk button on the toolbar or by pressing Cmd+Option+I.

22+4
## [1] 26

¿Y para qué tanta molestia si apenas vamos cachando qué onda con el Script?

El Script me interesa porque es una funcionalidad similar a la de otros softwares, como STATA con el dofile o SPSS con syntax

Pero lo bueno, bueno, es lo siguiente:

When you save the notebook, an HTML file containing the code and output will be saved alongside it (click the Preview button or press Cmd+Shift+K to preview the HTML file).

The preview shows you a rendered HTML copy of the contents of the editor. Consequently, unlike Knit, Preview does not run any R code chunks. Instead, the output of the chunk when it was last run in the editor is displayed.

Juntos iremos aprendiendo poco a poco a usar Rmarkdown.

Mientras tanto les dejo varios recursos para que le entrena su ritmo:

En las últimas semanas hubo un atractivo de Curso de Markdown en la IBERO