Intervalo para la media muestral bajo normalidad (\(\sigma\) conocida)

# Parámetros poblacionales
n <- 30; sigma <- 10; alpha <- 0.05

# Generamos la muestra
x <- rnorm(n, mean = 52, sd = sigma)

# IC
mean(x) - qnorm(alpha/2, lower.tail = FALSE)*sigma/sqrt(n)

[1] 47.54094

mean(x) + qnorm(alpha/2, lower.tail = FALSE)*sigma/sqrt(n)

[1] 54.69772

Visualización del intervalo

Test de Hipótesis para la media bajo normalidad (\(\sigma\) conocida)

# Parámetros poblacionales
n <- 30; sigma <- 10; mu0 <- 50; alpha <- 0.05

# Generamos la muestra
x <- rnorm(n, mean = 52, sd = sigma)

# Calculamos el estadístico del contraste
z <- abs(mean(x) - mu0) / (sigma / sqrt(n))

#Rechazo si
z > qnorm(0.05/2, lower.tail = FALSE) 

[1] FALSE

# Obtenemos el pvalor
2 * (1 - pnorm(z))  # valor p

[1] 0.06005003

Ejercicio de simulación

Compruebe por simulación que el Intervalo de Confianza al \(95 \%\) continene efectivamente el parámetro poblacional \(\mu\) en esa proporción.

Cobertura empírica: 0.955 

Intervalo para la media muestral bajo normalidad (\(\sigma\) desconocida)

• \(X_1, \dots, X_n \sim \mathcal{N}(\mu, \sigma^2)\), \(\sigma\) desconocida.

\[ \bar{x} \pm t_{n-1, 1 - \alpha/2} \cdot \frac{s}{\sqrt{n}} \]

Test de hipótesis para la media muestral bajo normalidad (\(\sigma\) desconocida)

t.test(x)

    One Sample t-test

data:  x
t = 25.962, df = 29, p-value < 2.2e-16
alternative hypothesis: true mean is not equal to 0
95 percent confidence interval:
 49.22386 57.64248
sample estimates:
mean of x 
 53.43317 

Intervalo de confianza para la proporción

• \(X \sim \text{Binomial}(n, p)\), estimador \(\hat{p} = X/n\).

\[ \hat{p} \pm z_{1 - \alpha/2} \cdot \sqrt{ \frac{\hat{p}(1 - \hat{p})}{n} } \]

Test de hipótesis para la proporción

x <- rbinom(1000, size = 1, prob = 0.8)
mean(x) + c(-1,1) * qnorm(0.975) * sqrt(mean(x)*(1-mean(x))/n)

[1] 0.6380391 0.9539609

t.test(x, mu = 0.5)

    One Sample t-test

data:  x
t = 23.217, df = 999, p-value < 2.2e-16
alternative hypothesis: true mean is not equal to 0.5
95 percent confidence interval:
 0.7709814 0.8210186
sample estimates:
mean of x 
    0.796 

prop.test(sum(x), length(x), p = 0.5, correct = FALSE)

    1-sample proportions test without continuity correction

data:  sum(x) out of length(x), null probability 0.5
X-squared = 350.46, df = 1, p-value < 2.2e-16
alternative hypothesis: true p is not equal to 0.5
95 percent confidence interval:
 0.7699136 0.8198210
sample estimates:
    p 
0.796 

Intevalo de confianza para la varianza bajo normalidad

• \(X_1, \dots, X_n \sim \mathcal{N}(\mu, \sigma^2)\).

\[ \left[ \frac{(n-1)s^2}{\chi^2_{1 - \alpha/2}}, \frac{(n-1)s^2}{\chi^2_{\alpha/2}} \right] \]

Test de hipótesis para la varianza bajo normalidad

x <- rnorm(20, sd = 2)
s2 <- var(x)
n <- length(x)
sigma0_2 <- 4 

# Estadístico
chi2 <- (n - 1) * s2 / sigma0_2

# Valor p bilateral
2 * min(pchisq(chi2, df = n - 1), 1 - pchisq(chi2, df = n - 1))

[1] 0.04707332

# valor p H1 sigma>sigma_0
1 - pchisq(chi2, df = n-1)

[1] 0.02353666

# valor p H1 sigma<sigma_0
pchisq(chi2, df = n-1)

[1] 0.9764633

Diferencia de medias (independientes)

x <- rnorm(30, mean = 100, sd = 10)
y <- rnorm(30, mean = 95, sd = 10)
t.test(x, y, var.equal = TRUE)

    Two Sample t-test

data:  x and y
t = 2.9128, df = 58, p-value = 0.005078
alternative hypothesis: true difference in means is not equal to 0
95 percent confidence interval:
  2.415902 13.031272
sample estimates:
mean of x mean of y 
100.32036  92.59678 

Diferencia de medias (pareadas)

before <- rnorm(20, mean = 100)
after  <- before + rnorm(20, mean = -5)
t.test(before, after, paired = TRUE)

    Paired t-test

data:  before and after
t = 25.32, df = 19, p-value = 4.21e-16
alternative hypothesis: true mean difference is not equal to 0
95 percent confidence interval:
 4.733021 5.586008
sample estimates:
mean difference 
       5.159514 

Diferencia de proporciones

x <- c(40, 60)
n <- c(100, 100)
prop.test(x, n)

    2-sample test for equality of proportions with continuity correction

data:  x out of n
X-squared = 7.22, df = 1, p-value = 0.00721
alternative hypothesis: two.sided
95 percent confidence interval:
 -0.34579029 -0.05420971
sample estimates:
prop 1 prop 2 
   0.4    0.6 

smokers  <- c(83, 90, 129, 70)
patients <- c(86, 93, 136, 82)
prop.test(smokers, patients)

    4-sample test for equality of proportions without continuity correction

data:  smokers out of patients
X-squared = 12.6, df = 3, p-value = 0.005585
alternative hypothesis: two.sided
sample estimates:
   prop 1    prop 2    prop 3    prop 4 
0.9651163 0.9677419 0.9485294 0.8536585 

Comparación de varianzas (dos muestras)

# Simulación de dos muestras
x <- rnorm(25, sd = 2)
y <- rnorm(30, sd = 3)

# Test F de igualdad de varianzas
var.test(x, y)

    F test to compare two variances

data:  x and y
F = 0.33037, num df = 24, denom df = 29, p-value = 0.007167
alternative hypothesis: true ratio of variances is not equal to 1
95 percent confidence interval:
 0.1533751 0.7325787
sample estimates:
ratio of variances 
         0.3303709 

ANOVA

• Objetivo: comparar medias de \(k\) poblaciones
• Hipótesis: \[ H_0: \mu_1 = \mu_2 = \dots = \mu_k \quad vs \quad H_1: \text{al menos una media difiere} \]
• Supuestos:
  • Independencia de las observaciones.
  • Normalidad en cada grupo.
  • Igualdad de varianzas \(\sigma^2\).

Modelo

\[ Y_{ij} = \mu_i + \varepsilon_{ij}, \quad \varepsilon_{ij} \sim \mathcal{N}(0, \sigma^2) \]

• \(i = 1, \dots, k\): grupo
• \(j = 1, \dots, n_i\): observación dentro del grupo
• \(N\): número total de observaciones

Se basa en la descomposición de la variabilidad suma de cuadrados total (SCT):

\[ \text{SCT} = \sum_{i=1}^k \sum_{j=1}^{n_i} (Y_{ij} - \bar{Y})^2 \]

\[ \text{SCT} = \underbrace{\sum_{i=1}^k n_i (\bar{Y}_i - \bar{Y})^2}_{\text{SCB}} + \underbrace{\sum_{i=1}^k \sum_{j=1}^{n_i} (Y_{ij} - \bar{Y}_i)^2}_{\text{SCD}} \]

Es la suma de la Suma de Cuadrados entre (between) Grupos (SCB) y la Suma de Cuadrados dentro de los Grupos (SCD).

Estadístico de contraste

• Media cuadrática entre grupos: \[\text{MCE} = \frac{\text{SCB}}{k - 1}\]
• Media cuadrática dentro de grupos: \[\text{MCI} = \frac{\text{SCD}}{N - k}\]

Estadístico:

\[F = \frac{\text{MCE}}{\text{MCI}} \sim F_{k-1, N-k} \quad \text{(bajo $H_0$)}\]

Regla: Rechazar \(H_0\) si: \[ F_{\text{obs}} > F_{k-1, N-k}^{(1-\alpha)}\]

ANOVA en R

group <- factor(rep(1:3, each = 20))
valores <- c(rnorm(20, 5), rnorm(20, 5), rnorm(20, 6))
anova_model <- aov(valores ~ group)
summary(anova_model)

            Df Sum Sq Mean Sq F value  Pr(>F)   
group        2  12.56   6.282   5.602 0.00601 **
Residuals   57  63.92   1.121                   
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Test de Chi-cuadrado (χ²)

• Evaluar si existe asociación entre dos variables categóricas (test de independencia).
• Evaluar si una distribución observada se ajusta a una teórica (bondad de ajuste).

Estadístico

\[X^2 = \sum_{i=1}^{k} \frac{(O_i - E_i)^2}{E_i}\]

• \(O_i\): frecuencia observada en la categoría \(i\)
  
• \(E_i\): frecuencia esperada bajo la hipótesis nula
  
• \(k\): número de categorías o celdas

Distribución bajo \(H_0\)

\[X^2 \overset{H_0}{\sim} \chi^2_{\nu}\]

• \(\nu = (r - 1)(c - 1)\) en una tabla de contingencia \(r \times c\)

• \(\nu = k - 1\) (test de bondad de ajuste sin parámetros estimados)

Test de Chi-cuadrado en R: Independencia

# Ejemplo: tabla 3x2
# filas: Grupo A, Grupo B, Grupo C
# columnas: Sí, No

tabla <- matrix(c(20, 15, 30,   # columna 1 (Sí) por fila
                  10, 25, 15),  # columna 2 (No)
                nrow = 3, byrow = FALSE)
rownames(tabla) <- c("Grupo1", "Grupo2", "Grupo3")
colnames(tabla) <- c("Si", "No")

tabla

       Si No
Grupo1 20 10
Grupo2 15 25
Grupo3 30 15

res <- chisq.test(tabla)
res

    Pearson's Chi-squared test

data:  tabla
X-squared = 9.0304, df = 2, p-value = 0.01094

# Mostrar tablas de frecuencias esperadas
res$expected

             Si       No
Grupo1 16.95652 13.04348
Grupo2 22.60870 17.39130
Grupo3 25.43478 19.56522

Test de Chi-cuadrado en R: Bondad de Ajuste

# Datos: frecuencias observadas
observed <- c(50, 30, 20)

# Frecuencias esperadas (proporciones teóricas)
expected_prop <- c(0.5, 0.3, 0.2)
expected <- sum(observed) * expected_prop # en valores absolutos

# Test de bondad de ajuste
chisq.test(x = observed, p = expected_prop) 

    Chi-squared test for given probabilities

data:  observed
X-squared = 0, df = 2, p-value = 1

Test de Kolmogorov-Smirnov (K-S)

Objetivo:

• Evaluar si una muestra proviene de una distribución teórica específica (una muestra)

• Evaluar si dos muestras provienen de la misma distribución (dos muestras).

Estadístico (una muestra)

\[D_n = \sup_x \left| F_n(x) - F(x) \right|\]

• \(F_n(x)\): función de distribución empírica (ECDF)
  
• \(F(x)\): función de distribución teórica bajo \(H_0\)

Estadístico (dos muestras)

\[D_{n,m} = \sup_x \left| F_n(x) - G_m(x) \right|\]

• \(F_n(x)\): ECDF de la primera muestra (tamaño \(n\))
  
• \(G_m(x)\): ECDF de la segunda muestra (tamaño \(m\))

Distribución bajo \(H_0\)

• Para una muestra:

\[\sqrt{n} D_n \xrightarrow{d} \sup_{t \in [0,1]} |B(t)|\]

donde \(B(t)\) es el puente browniano (Kolmogorov distribution)

• Para dos muestras:

\[ \sqrt{\frac{nm}{n + m}} D_{n,m} \xrightarrow{d} \sup_{t \in [0,1]} |B(t)| \]

Test de Kolmogorov-Smirnov en R (una muestra)

• Contrasta si una muestra sigue una distribución teórica.

Hipótesis:

\[H_0: F(x) = F_0(x)\]

\[ H_1$: F(x) \ne F_0(x)\]

x <- runif(100)
ks.test(x, "punif", min = 0, max = 1)

    Asymptotic one-sample Kolmogorov-Smirnov test

data:  x
D = 0.080825, p-value = 0.5308
alternative hypothesis: two-sided

Test de Kolmogorov-Smirnov en R (dos muestras)

• Contrasta si dos muestras provienen de la misma distribución.

Hipótesis:

\[H_0: F_1(x) = F_2(x)\]

\[H_1: F_1(x) \ne F_2(x)\]

x <- rnorm(100)
y <- rnorm(100, mean = 0.5)
ks.test(x, y)

    Asymptotic two-sample Kolmogorov-Smirnov test

data:  x and y
D = 0.28, p-value = 0.0007873
alternative hypothesis: two-sided

Resumen Test de hipótesis una muestra

Resumen Test de hipótesis varias muestra

Máxima Verosímilud

La función de verosimilitud o likelihood de un vector aleatorio \(X_1,\dots, X_n\) cuya distribución depende de un parámetro o vector de parámetros \(\theta\) se define como la función de densidad o de probabilidad conjunta

\[L(\theta) = f_{X_1, \dots, X_n}(x_1, \dots, x_n; \theta).\] con dominio en \(\mathbb{R}^k\).

Bajo muestras aleatorias independiente e identicamente distribuídas, la función de verosimilitud se reduce a

\[L(\theta) = \prod_{i=1}^n f_X(x; \theta).\]

Los estimadores de máxima verosimilitud vienen de encontrar el valor de \(\hat{\theta}\)que maximiza la verosimilitud. Esto es

\[\hat{\theta}_{\text{MV}} = \arg\max_{\theta \in \Theta} L(\theta).\]

Verosimilitud en R

Dada una muestra de una Poisson:

muestra <- rexp(100, rate = 3)

Supongamos que desconocemos su parámetro y queremos obtener una estimación por máxima verosimilitud.

Para ello, vamos a evaluar la función de verosimilitud en diferentes valores del parámetro y nos quedaremos con aquel que la maximice.

grid_rate <- seq(1, 5, by = 0.1)

verosimilitud <- sapply(grid_rate, function(x) prod(dexp(muestra, rate = x)))

print(paste0("Estimador Máximo Verosímil: ", grid_rate[which.max(verosimilitud)]))

[1] "Estimador Máximo Verosímil: 3.2"

Podemos graficar la verosimilutud para cada valor del parámetro:

plot(grid_rate, verosimilitud, type = "l")
points(grid_rate[which.max(verosimilitud)], max(verosimilitud), col = "red", pch = 16)

Motivación

• En muchos problemas estadísticos, la distribución muestral de un estimador es desconocida o difícil de derivar analíticamente o difícil de aproximar asintóticamente.

El bootstrap es un enfoque computacional desarrollado por Bradley Efron (1979) que permite aproximar la distribución muestral de un estadístico a partir de la muestra observada, sin necesidad de suposiciones fuertes sobre la población.

Si la muestra representa bien a la población, ¡podemos simular nuevas muestras re-muestreando con reemplazo!

Fundamento teórico

Sea \(X = (X_1, X_2, \dots, X_n)\) una muestra i.i.d. de una población desconocida con función de distribución \(F\), y sea \(T = T(X_1, \dots, X_n)\) un estadístico cualquiera.

El bootstrap clásico consiste en:

1. Construir una distribución empírica \(\hat{F}_n\) asignando masa \(1/n\) a cada observación de la muestra.

2. Generar \(B\) muestras de tamaño \(n\) con reemplazo de \(\hat{F}_n: X^{*b} = (X_1^{*b}, \dots, X_n^{*b})\), para \(b = 1, \dots, B\).

3. Calcular el estadístico en cada muestra: \(\hat{T^{*b}}\).

4. Utilizar \(\hat{\{T^{*1}}, \dots, \hat{T^{*B}\}}\) como aproximación a la distribución muestral de \(T\).

Ejemplo Bootstrap para la Mediana

# Muestra original
x <- c(-0.38, 0.65, 0.5, 1.47, -1.4, 1.32, -0.54, -1.5, -0.65, 0.26)
t_hat <- median(x)

# Bootstrap
B <- 1000
t_boot <- vector(mode = "numeric", B)*NA

for(b in 1:B){
  x_star <- sample(x, size = length(x), replace = TRUE)
  t_boot[b] <- median(x_star)
}

IC Bootstrap 95%:

round(quantile(t_boot, probs = c(alpha/2, 1-alpha/2)), digits = 2)

 2.5% 97.5% 
-0.97  0.79 

Aproximación a la distribución Bootstrapp e IC

# Muestra los Resultados
hist(t_boot, freq = FALSE, 
     main = "Distribución Bootstrap de la Mediana",
     xlab = "Mediana Bootstrap")
abline(v = t_hat, col = "red", lwd = 2, lty = 2)
abline(v = quantile(t_boot, probs = c(alpha/2, 1-alpha/2)), col = "blue", lwd = 2, lty = 2)