Inferencia sobre una población

Autor/a

Edimer David Jaramillo

Fecha de publicación

22 de mayo de 2025

Bibliotecas

Código
library(tidyverse)
library(readxl)
library(infer)

theme_set(theme_bw())

Datos

Código
datos <- read_csv("../datos/Encuesta_Motociclistas.csv") |>
  select(
    municipio,
    sexo = hombre,
    nivel_educativo,
    herramienta_trabajo,
    experiencia,
    cilin_grupo,
    gasto_anual,
    licencia_moto,
    licencia_curso,
  ) |>
  mutate(
    municipio = str_to_title(municipio),
    sexo = if_else(sexo == 1, "Hombre", "Mujer"),
    herramienta_trabajo = if_else(herramienta_trabajo == 0, "No", "Sí"),
    licencia_moto = if_else(licencia_moto == 0, "No", "Sí"),
    licencia_curso = if_else(licencia_curso == 0, "No", "Sí"),
    nivel_educativo = factor(
      nivel_educativo,
      c(
        "Primaria o menos",
        "Secundaria",
        "Técnica / Tecnológica",
        "Universitaria o postgrado"
      )
    ),
    experiencia = factor(experiencia, levels = c("0-2", "3-5", "6-10", "11-20", ">20")),
    cilin_grupo = factor(cilin_grupo, levels = c("<100", "100-125", "126-200", ">200"))
  ) |> 
  filter(!is.na(gasto_anual))

datos

Escenarios de inferencia

Hipótesis estadística

Tipos de pruebas

Tipos de errores

Inferencia sobre una población

Inferencia sobre la media (\(\mu\))

Normalidad de la variable

Código
ggpubr::ggqqplot(datos$gasto_anual)

Juego de hipótesis

\[H_0: \mu = 1.000.000\]

\[H_1: \mu \neq 1.000.000\]

Nivel de significancia

  • Vamos a utilizar un nivel de significancia del 5% (\(\alpha = 0.05\))

Calcular el estadístico observado

  • En este caso calculamos el estadístico T:

\[T = \frac{\bar{X} - \mu}{S/\sqrt{n}}\]

Código
x_barra <- mean(datos$gasto_anual, na.rm = TRUE)
mu_referencia <- 1e+06
desviacion_muestral <- sd(datos$gasto_anual, na.rm = TRUE)
raiz_n <- sqrt(nrow(datos))

\[T = \frac{1039744 - 1000000}{713080/157.31} = 8.768062\]

Código
(x_barra - mu_referencia) / (desviacion_muestral / raiz_n)
[1] 8.768062

Región de rechazo

Código
qt(p = 0.025, df = nrow(datos) - 1, lower.tail = TRUE)
[1] -1.96006
Código
qt(p = 0.025, df = nrow(datos) - 1, lower.tail = FALSE)
[1] 1.96006

  • Conclusión: como el estadístico calculado (8.768062) está dentro de las regiones de rechazo, existe evidencia para rechazar la hipótesis nula, es decir, que se rechaza que el promedio del gasto anual sea igual a un millón de pesos. Esto se concluye con un nivel de significancia del 5%.

Intervalo de confianza para \(\mu\)

  • Límite inferior del intervalo de confianza:

\[\bar{X} - t_{\alpha/2, n-1} \times \frac{s}{\sqrt{n}}\]

Código
x_barra - (1.96006 * (desviacion_muestral / raiz_n))
[1] 1030859
  • Límite superior del intervalo de confianza:

\[\bar{X} + t_{\alpha/2, n-1} \times \frac{s}{\sqrt{n}}\]

Código
x_barra + (1.96006 * (desviacion_muestral / raiz_n))
[1] 1048629
  • Conclusión: como el valor de referencia (1.000.000) no está dentro del intervalo de confianza, existe evidencia para rechazar la hipótesis nula, es decir, que se rechaza que el promedio del gasto_anual sea igual a un millón de pesos. Esto se concluye con un nivel de significancia del 5%. Como el intervalo de confianza está a la derecha del valor de referencia (1.000.000) podemos afirmar que la media del gasto anual de los motociclistas en colombia es mayor a 1.000.000

Valor p

  • Calcular el área que deja un valor de 8.768062 a la izquierda:
Código
pt(q = -8.768062, df = nrow(datos) - 1, lower.tail = TRUE)
[1] 9.661257e-19
  • Calcular el área que deja un valor de 8.768062 a la derecha:
Código
pt(q = 8.768062, df = nrow(datos) - 1, lower.tail = FALSE)
[1] 9.661257e-19
  • El valor p es la suma de las dor áreas anteriores:
Código
9.661257e-19 + 9.661257e-19
[1] 1.932251e-18
  • Conclusión: como el valor p (1.932251e-18) es menor que el nivel de significancia (0.05) existe evidencia para rechazar la hipótesis nula.

Solución con R

  • Utilizamos la función t.test() con los siguientes argumentos:
    • x: la variable sobre la cual estamos haciendo inferencia. En este caso el gasto_anual
    • alternative: tipo de hipótesis alternativa. En este es una prueba bilateral usamos “two.sided”
    • conf.level: nivel de confianza (1 - nivel de significancia = 1 - 0.05 = 0.95)
    • mu: valor promedio de referencia. En este caso es 1.000.000
Código
t.test(x = datos$gasto_anual,
       alternative = "two.sided",
       conf.level = 0.95,
       mu = 1e+06)

    One Sample t-test

data:  datos$gasto_anual
t = 8.7681, df = 24747, p-value < 2.2e-16
alternative hypothesis: true mean is not equal to 1e+06
95 percent confidence interval:
 1030859 1048629
sample estimates:
mean of x 
  1039744 
  • Presentando resultados de la prueba anterior, de forma ordenada (tidy):
Código
prueba_t1 <- t.test(
  x = datos$gasto_anual,
  alternative = "two.sided",
  conf.level = 0.95,
  mu = 1e+06
)

library(broom)
prueba_t1 |> tidy()

Alternativa no parámetrica 1: Test de Wilcoxon

Código
wilcox.test(
  x = datos$gasto_anual,
  alternative = "two.sided",
  conf.int = TRUE,
  conf.level = 0.95,
  mu = 1e+06
)

    Wilcoxon signed rank test with continuity correction

data:  datos$gasto_anual
V = 114442553, p-value < 2.2e-16
alternative hypothesis: true location is not equal to 1e+06
95 percent confidence interval:
 925000 945000
sample estimates:
(pseudo)median 
        935000 
  • Conclusión: aunque la hipótesis nula se rechaza, el intervalo de confianza con esta prueba no paramétrica nos muestra que el valor real del gasto anual es inferior a un millón de pesos, estaría entre 925000 y 945000.

Alternativa no parámetrica 2: Bootstrapping

Boootstrapping

Proceso con biblioteca infer

    1. Especificar las variables con la función specify()
    1. Generar réplicas (remuestreo) con la función generate()
    1. Calcular la estadísticas de resumen con la función calculate()
    1. Visualice los resultados con la función visualize()
    1. (opcional) Construir intervalos de confianza con la función get_confidence_interval(). Nota: para mejorar la visualización de los intervalos de confianza, se puede utilizar la función shade_confidence_interval()

Código
set.seed(2025)
bootstrap_gasto_anual <- 
  datos |> 
  specify(response = gasto_anual) |> 
  generate(reps = 1000, type = "bootstrap") |> 
  calculate(stat = "mean")

bootstrap_gasto_anual
  • Graficamos la distribución bootstrap:
Código
bootstrap_gasto_anual |> 
  visualize()

Intervalos de confianza

{fig-align=“center” width = “70%”}

Método de percentiles

  • Intervalo de confianza con percentiles (95%):
Código
ic_promedio_percentil <-
  bootstrap_gasto_anual |> 
  get_confidence_interval(level = 0.95, type = "percentile")
ic_promedio_percentil
  • Graficamos el intervalo de confianza obtenido con el método de percentiles:
Código
bootstrap_gasto_anual |> 
  visualize() +
  shade_confidence_interval(endpoints = ic_promedio_percentil) +
  geom_vline(
    xintercept = x_barra,
    color = "red",
    lty = 2,
    size = 1.5
  ) +
  geom_vline(
    xintercept = mean(bootstrap_gasto_anual$stat),
    color = "black",
    lty = 2,
    size = 1.5
  )

Método de error estándar

  • Intervalo de confianza con método del error estándar:
Código
ic_promedio_error_est <-
  bootstrap_gasto_anual |> 
  get_confidence_interval(type = "se", point_estimate = x_barra)
ic_promedio_error_est
  • Graficamos el intervalo de confianza obtenido con el método de percentiles:
Código
bootstrap_gasto_anual |> 
  visualize() +
  shade_confidence_interval(endpoints = ic_promedio_error_est) +
  geom_vline(
    xintercept = x_barra,
    color = "red",
    lty = 2,
    size = 1.5
  ) +
  geom_vline(
    xintercept = mean(bootstrap_gasto_anual$stat),
    color = "black",
    lty = 2,
    size = 1.5
  )

:::

Inferencia sobre un proporción (\(p\))

Hipótesis

  • ¿Proporción de motociclistas que obtuvieron su licencia a través de un curso?

\[H_0: p = 50\%\]

\[H_1: p \neq 50\%\]

Proporción muestral

Código
datos2 <- 
  datos |> 
  filter(!is.na(licencia_curso))

datos2$licencia_curso |> 
  table() |> 
  prop.table()

       No        Sí 
0.2835426 0.7164574 

Remuestreo - Bootstrapping

  • ¿Cuál es la verdadera proporción de motociclistas que obtuvieron su licencia a través de un curso?
Código
set.seed(2025)
remuestreo_licencia_curso <- datos2 |>
  specify(response = licencia_curso, success = "Sí") |>
  generate(reps = 1000, type = "bootstrap") |>
  calculate(stat = "prop")

remuestreo_licencia_curso
  • Graficamos la distribución del remuestreo:
Código
remuestreo_licencia_curso |> 
  visualize()

  • Calculamos el intervalo de confianza por el método de percentiles:
Código
ic_perc_licencia_curso <-
  remuestreo_licencia_curso |>
  get_confidence_interval(level = 0.95, type = "percentile")
ic_perc_licencia_curso
  • Graficamos los intervalos de confianza con el método de percentiles:
Código
remuestreo_licencia_curso |>
  visualize() +
  shade_confidence_interval(endpoints = ic_perc_licencia_curso)

Modelo matemático (paramétrico)

\[\hat{p}-Z_{\alpha/2}\sqrt{\frac{\hat{p}(1-\hat{p})}{n}} < p < \hat{p}+Z_{\alpha/2}\sqrt{\frac{\hat{p}(1-\hat{p})}{n}}\]

  • Conteo:
Código
datos2 |> 
  count(licencia_curso)
  • Ahora aplicamos la prueba estadística para proporciones:
Código
total_muestra <- datos2 |> nrow()
total_licencia_curso <- 13809

prop.test(
  x = total_licencia_curso,
  n = total_muestra,
  conf.level = 0.95,
  p = 0.5
)

    1-sample proportions test with continuity correction

data:  total_licencia_curso out of total_muestra, null probability 0.5
X-squared = 3611.4, df = 1, p-value < 2.2e-16
alternative hypothesis: true p is not equal to 0.5
95 percent confidence interval:
 0.7100256 0.7228026
sample estimates:
        p 
0.7164574 

Referencia