Bibliotecas

library(readxl)
library(tidyverse)
library(janitor)
library(ggpubr)
library(qqplotr)
library(car)
library(nortest)
library(infer)

Medias muestrales

Población simulada

set.seed(2022)
poblacion <- rnorm(n = 1000000, mean = 92.5, sd = 6.8)

Muestra aleatoria

Obtenemos una muestra (un experimento):

set.seed(2022)
muestra <- sample(x = poblacion, size = 100, replace = FALSE)

Calculamos los estadísticos:

mean(muestra)

## [1] 92.7529

sd(muestra)

## [1] 6.311283

Varias muestras (experimentos)

Convertimos el vector “poblacion” a tipo base de datos (dataframe) con la función enframe():

df_poblacion <-
  poblacion %>%
  enframe(value = "peso_gramos") %>%
  select(-name)

df_poblacion %>% head()

## # A tibble: 6 x 1
##   peso_gramos
##         <dbl>
## 1        98.6
## 2        84.5
## 3        86.4
## 4        82.7
## 5        90.2
## 6        72.8

Simulación 1

# N = 10
set.seed(2022)
simulacion_50 <-
  df_poblacion %>%
  rep_sample_n(size = 10,
               replace = FALSE,
               reps = 10) %>%
  group_by(replicate) %>%
  summarise(promedio_peso_gramos = mean(peso_gramos)) %>% 
  ungroup() %>% 
  mutate(N = 10)

# N 100
set.seed(2022)
simulacion_100 <-
  df_poblacion %>%
  rep_sample_n(size = 100,
               replace = FALSE,
               reps = 10) %>%
  group_by(replicate) %>%
  summarise(promedio_peso_gramos = mean(peso_gramos)) %>% 
  ungroup() %>% 
  mutate(N = 100)

# N = 1000
set.seed(2022)
simulacion_1000 <-
  df_poblacion %>%
  rep_sample_n(size = 1000,
               replace = FALSE,
               reps = 10) %>%
  group_by(replicate) %>%
  summarise(promedio_peso_gramos = mean(peso_gramos)) %>% 
  ungroup() %>% 
  mutate(N = 1000)

# Total de datos
datos_total <- 
  bind_rows(simulacion_50, simulacion_100, simulacion_1000) 

datos_total %>%
  ggplot(aes(x = promedio_peso_gramos)) +
  facet_wrap( ~ N) +
  geom_density(fill = "blue",
               alpha = 0.5,
               color = "blue") +
  geom_vline(xintercept = 92.5,
             lty = 2,
             color = "red")

Simulación 2

set.seed(2022)
prueba <-
  df_poblacion %>%
  rep_sample_n(size = 10,
               replace = FALSE,
               reps = 200) %>%
  group_by(replicate) %>%
  summarise(promedio_peso_gramos = mean(peso_gramos)) %>%
  ungroup()

prueba %>%
  ggplot(aes(x = promedio_peso_gramos)) +
  geom_density(fill = "blue",
               alpha = 0.5,
               color = "blue") +
  geom_vline(xintercept = 92.5,
             lty = 2,
             color = "red") +
  geom_vline(
    xintercept = mean(prueba$promedio_peso_gramos),
    size = 0.8,
    color = "black"
  )

Inferencia

Datos ratones

Fuente de datos

ratones <- read_excel("experimento_ratones.xlsx")
ratones

## # A tibble: 16 x 4
##    Block Treatment Strain    GST
##    <chr> <chr>     <chr>   <dbl>
##  1 A     Control   NIH       444
##  2 A     Treated   NIH       614
##  3 A     Control   BALB/C    423
##  4 A     Treated   BALB/C    625
##  5 A     Control   A/J       408
##  6 A     Treated   A/J       856
##  7 A     Control   129/Ola   447
##  8 A     Treated   129/Ola   719
##  9 B     Control   NIH       764
## 10 B     Treated   NIH       831
## 11 B     Control   BALB/C    586
## 12 B     Treated   BALB/C    782
## 13 B     Control   A/J       609
## 14 B     Treated   A/J      1002
## 15 B     Control   129/Ola   606
## 16 B     Treated   129/Ola   766

Datos créditos

Fuente de datos

creditos <- read_rds("creditos_colombia_genero.Rds")
creditos

## # A tibble: 461,269 x 2
##    departamento_inversion genero
##    <chr>                  <chr> 
##  1 BOYACÁ                 M     
##  2 CAUCA                  H     
##  3 ANTIOQUIA              H     
##  4 VALLE DEL CAUCA        M     
##  5 BOLÍVAR                H     
##  6 CUNDINAMARCA           H     
##  7 ANTIOQUIA              H     
##  8 TOLIMA                 M     
##  9 TOLIMA                 M     
## 10 TOLIMA                 M     
## # ... with 461,259 more rows

Una población

Una media

Ejemplo: en este caso vamos a hacer inferencia sobre la media de GST.

Solución manual

Vamos a suponer que queremos responder a la siguinete hipótesis para la variable GST:

\[H_0: \mu = 500 \\ H_1: \mu \neq 500\]

El nivel de significancia será del 5% (\(\alpha = 0.05\))
Comprobamos si la variable de interés se distribuye de forma normal:

ggqqplot(ratones$GST)

Calculamos el estadístico de prueba:

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

x_barra <- mean(ratones$GST)
mu_referencia <- 500
desviacion_muestral <- sd(ratones$GST)
raiz_n <- sqrt(nrow(ratones))

\[T = \frac{655.125 - 500}{173.3378/4} = 3.579716\]

(x_barra - mu_referencia) / (desviacion_muestral / raiz_n)

## [1] 3.579716

Podemos obtener los límites critícos con R:

qt(p = 0.025, df = nrow(ratones) - 1, lower.tail = TRUE)

## [1] -2.13145

qt(p = 0.025, df = nrow(ratones) - 1, lower.tail = FALSE)

## [1] 2.13145

Límite inferior del intervalo de confianza:

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

x_barra - (2.13145 * (desviacion_muestral / raiz_n))

## [1] 562.7598

Límite superior del intervalo de confianza:

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

x_barra + (2.13145 * (desviacion_muestral / raiz_n))

## [1] 747.4902

Calcular el área que deja un valor de 3.579716 a la izquierda:

pt(q = -3.579716, df = nrow(ratones) - 1, lower.tail = TRUE)

## [1] 0.001368596

Calcular el área que deja un valor de 3.579716 a la derecha:

pt(q = 3.579716, df = nrow(ratones) - 1, lower.tail = FALSE)

## [1] 0.001368596

El valor p es la suma de las dor áreas anteriores:

0.001368596 + 0.001368596

## [1] 0.002737192

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 GST
- 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 500

t.test(x = ratones$GST,
       alternative = "two.sided",
       conf.level = 0.95,
       mu = 500)

## 
##  One Sample t-test
## 
## data:  ratones$GST
## t = 3.5797, df = 15, p-value = 0.002737
## alternative hypothesis: true mean is not equal to 500
## 95 percent confidence interval:
##  562.7598 747.4902
## sample estimates:
## mean of x 
##   655.125

Una proporción

Ejemplo: En este caso estamos interesados en la proporción de mujeres que se hicieron acreedoras de créditos agropecuarios en el año 2021.

Juepo de hipótesis

\[H_0: p = 0.5 \\ H_1: p \neq 0.5\]

Nivel de significancia

En este caso usamos un nivel de significancia del 5%

Calcular el estadístico

Calculamos la proporción de mujeres que se hicieron acreedoras a créditos agropecuarios:

prop.table(table(creditos$genero))

## 
##        H        M 
## 0.618123 0.381877

Solución con R

Usamos la función prop.test() con los siguientes argumentos:
- x: el número de éxitos para el evento de interés. En este caso son mujers que recibieron créditos (n = 176148)
- n: total de ensayos (observaciones). En este caso equivale a 461269.
- p: proporción de referencia. En este caso es 0.5
- alternative: tipo de hipótesis alternativa. En este caso es “two.sided”
- conf.level: nivel de confianza. En este caso es 1 - 0.05 = 0.95

prop.test(
  x = 176148,
  n = 461269,
  p = 0.5,
  alternative = "two.sided",
  conf.level = 0.95
)

## 
##  1-sample proportions test with continuity correction
## 
## data:  176148 out of 461269, null probability 0.5
## X-squared = 25744, df = 1, p-value < 2.2e-16
## alternative hypothesis: true p is not equal to 0.5
## 95 percent confidence interval:
##  0.3804748 0.3832811
## sample estimates:
##        p 
## 0.381877

Dos poblaciones

Dos varianzas

Ejemplo: en este caso estamos interesado en responder si las varianzas de GST discrepan estadísticamente entre bloques (A vs B).

Hipótesis

\[H_0: \sigma^2_{A} / \sigma^2_{B} = 1 \\ H_1: \sigma^2_{A} / \sigma^2_{B} \neq 1\]

Nivel de significancia

En este caso vamos a usar un nivel de significancia del 5%.

Analizando la normalidad

Gráfico:

ggqqplot(data = ratones$GST)

Prueba de shapiro wilk:

shapiro.test(x = ratones$GST)

## 
##  Shapiro-Wilk normality test
## 
## data:  ratones$GST
## W = 0.94891, p-value = 0.4726

Solución con R

Vamos a usar la función var.test() con los siguientes argumentos:
- formula: y ~ x. En este caso “y” es la variable GST y el “x” es el bloque
- ratio: es el resultado del cociente de las dos varianzas. En este caso asumimos en la hipótesis nula el valor de “1”.
- alternative: tipo de prueba. En este es bilateral (“two.sided”)
- conf.level: nivel de confianza. En este caso es 0.95 (1 - 0.5)

var.test(ratones$GST ~ ratones$Block,
         ratio = 1,
         alternative = "two.sided",
         conf.level = 0.95)

## 
##  F test to compare two variances
## 
## data:  ratones$GST by ratones$Block
## F = 1.3537, num df = 7, denom df = 7, p-value = 0.6996
## alternative hypothesis: true ratio of variances is not equal to 1
## 95 percent confidence interval:
##  0.271016 6.761609
## sample estimates:
## ratio of variances 
##             1.3537

Prueba de Bartlett

Esta prueba es recomendada para comparar varianzas cuando la normalidad se cumple.

bartlett.test(ratones$GST ~ ratones$Block)

## 
##  Bartlett test of homogeneity of variances
## 
## data:  ratones$GST by ratones$Block
## Bartlett's K-squared = 0.14923, df = 1, p-value = 0.6993

Prueba de Levene

Esta prueba es recomendada para comparar varianzas cuando la normalidad no se cumple:

leveneTest(ratones$GST ~ ratones$Block)

## Levene's Test for Homogeneity of Variance (center = median)
##       Df F value Pr(>F)
## group  1  0.6138 0.4464
##       14

Dos medias

Ejemplo: en este caso vamos a comparar la variable GST para los dos tratamientos (control vs tratamiento).

Hipótesis

\[H_0: \mu_{control} = \mu_{tratamiento} \\ H1: \mu_{control} \neq \mu_{tratamiento}\]

El juego de hipótesis anterior es equivalente al siguiente:

\[H_0: \mu_{control} - \mu_{tratamiento} = 0 \\ H1: \mu_{control} - \mu_{tratamiento} \neq 0\]

Nivel de significancia

El nivel de significancia será del 5%.

Normalidad

Para la prueba t-student sobre dos poblaciones (medias) se debe garantizar que las variable (en cada grupo) se distribuye de forma normal.

grupo_control <- ratones %>% filter(Treatment == "Control")
grupo_tratamiento <- ratones %>% filter(Treatment == "Treated")

ggqqplot(grupo_control$GST, title = "Control")

ggqqplot(grupo_tratamiento$GST, title = "Tratamiento")

Pruebas de shapiro wilk:

shapiro.test(grupo_control$GST)

## 
##  Shapiro-Wilk normality test
## 
## data:  grupo_control$GST
## W = 0.8742, p-value = 0.1656

shapiro.test(grupo_tratamiento$GST)

## 
##  Shapiro-Wilk normality test
## 
## data:  grupo_tratamiento$GST
## W = 0.95484, p-value = 0.7598

Igualdad de varianzas

Primero debemos verificar si existe la homogeneidad en las varianzas, es decir, si las varianzas entre tratamiento y control no son diferentes.

var.test(ratones$GST ~ ratones$Treatment)

## 
##  F test to compare two variances
## 
## data:  ratones$GST by ratones$Treatment
## F = 0.97645, num df = 7, denom df = 7, p-value = 0.9757
## alternative hypothesis: true ratio of variances is not equal to 1
## 95 percent confidence interval:
##  0.1954891 4.8772803
## sample estimates:
## ratio of variances 
##          0.9764502

Solución con R

Usamos la función t.test() con los siguientes argumentos:
- formula: y ~ x. En este caso “y” es GST y “x” es el grupo (Control, Tratamiento)
- alternative: tipo de hipótesis alternativa. En este caso es bilateral.
- conf.level: nivel de confianza. En este caso es del 95%
- var.equal: toma valores TRUE o FALSE para cuando las varianzas son iguales o diferentes, respectivamente.

t.test(ratones$GST ~ ratones$Treatment,
       alternative = "two.sided",
       conf.level = 0.95,
       var.equal = TRUE)

## 
##  Two Sample t-test
## 
## data:  ratones$GST by ratones$Treatment
## t = -3.7781, df = 14, p-value = 0.002037
## alternative hypothesis: true difference in means between group Control and group Treated is not equal to 0
## 95 percent confidence interval:
##  -373.8939 -103.1061
## sample estimates:
## mean in group Control mean in group Treated 
##               535.875               774.375

Dos proporciones

Ejemplo: en este caso vamos a comprar la proporción de hombres y mujeres que se recibieron créditos agropecuarios en el año 2021.

Hipótesis

\[H_0: p_{hombres} = p_{mujeres} \\ H_1: p_{hombres} \neq p_{mujeres}\]

O de forma equivalente:

\[H_0: p_{hombres} - p_{mujeres} = 0 \\ H_1: p_{hombres} - p_{mujeres} \neq 0\]

Nivel de significancia

En este caso usaremos un nivel de significancia del 5% (0.05)

Calculando frecuencias

table(creditos$genero)

## 
##      H      M 
## 285121 176148

Solución con R

Usamos la función prop.test() con los siguientes argumentos:
- x: el número de éxitos para el evento de interés. En este caso son mujeres (n = 176148 ) y hombres (n = 285121) que recibieron créditos agropecuarios.
- n: total de ensayos (observaciones). En este caso equivale a 461269.
- alternative: tipo de hipótesis alternativa. En este caso es “two.sided”
- conf.level: nivel de confianza. En este caso es 1 - 0.05 = 0.95

prop.test(x = c(176148, 285121),
          n = c(461269, 461269),
          alternative = "two.sided",
          conf.level = 0.95)

## 
##  2-sample test for equality of proportions with continuity correction
## 
## data:  c(176148, 285121) out of c(461269, 461269)
## X-squared = 51488, df = 1, p-value < 2.2e-16
## alternative hypothesis: two.sided
## 95 percent confidence interval:
##  -0.2382311 -0.2342611
## sample estimates:
##   prop 1   prop 2 
## 0.381877 0.618123

Inferencia Estadística

Edimer David Jaramillo

2022-06-07

Bibliotecas

Medias muestrales

Población simulada

Muestra aleatoria

Varias muestras (experimentos)

Simulación 1

Simulación 2

Inferencia

Datos ratones

Datos créditos

Una población

Una media

Solución manual

Solución con R

Una proporción

Juepo de hipótesis

Nivel de significancia

Calcular el estadístico

Solución con R

Dos poblaciones

Dos varianzas

Hipótesis

Nivel de significancia

Analizando la normalidad

Solución con R

Prueba de Bartlett

Prueba de Levene

Dos medias

Hipótesis

Nivel de significancia

Normalidad

Igualdad de varianzas

Solución con R

Dos proporciones

Hipótesis

Nivel de significancia

Calculando frecuencias

Solución con R