Tutorial R: Pruebas de hipótesis en R

Test t para comparación de medias

La prueba t (o test t) es una de las herramientas estadísticas más utilizadas para comparar medias entre grupos. Esta prueba nos permite determinar si existe una diferencia estadísticamente significativa entre las medias de dos conjuntos de datos.

Fundamentos del test t

El test t se basa en la distribución t de Student y es especialmente útil cuando trabajamos con muestras pequeñas. La idea fundamental es contrastar dos hipótesis estadísticas:

• Hipótesis nula (H₀): No hay diferencia entre las medias de los grupos (μ₁ = μ₂)
• Hipótesis alternativa (H₁): Existe una diferencia entre las medias (μ₁ ≠ μ₂)

En R, podemos implementar diferentes variantes del test t según nuestras necesidades:

• Test t para muestras independientes
• Test t para muestras pareadas
• Test t de una muestra

Test t para muestras independientes

Este tipo de test se utiliza cuando queremos comparar las medias de dos grupos independientes, como por ejemplo, comparar el rendimiento entre un grupo de control y un grupo experimental.

Vamos a crear un ejemplo sencillo para ilustrar su uso:

# Creamos dos vectores con datos de ejemplo
grupo_control <- c(5.2, 5.8, 6.1, 5.9, 5.5, 6.2, 5.7, 6.0)
grupo_tratamiento <- c(6.5, 7.1, 6.8, 7.2, 6.9, 7.0, 7.3, 6.7)

# Visualizamos los datos con un diagrama de cajas
boxplot(grupo_control, grupo_tratamiento, 
        names = c("Control", "Tratamiento"),
        col = c("lightblue", "lightgreen"),
        main = "Comparación entre grupos")

Para realizar el test t, utilizamos la función t.test():

# Realizamos el test t para muestras independientes
resultado <- t.test(grupo_tratamiento, grupo_control, 
                    var.equal = FALSE)  # Asumimos varianzas desiguales

# Mostramos los resultados
resultado

La salida nos mostrará:

• El valor t calculado
• Los grados de libertad
• El valor p (p-value)
• El intervalo de confianza para la diferencia de medias
• Las medias de cada grupo

Interpretación de resultados

La interpretación del test t se basa principalmente en el valor p:

• Si p < 0.05 (nivel de significancia común), rechazamos la hipótesis nula y concluimos que existe una diferencia significativa entre las medias.
• Si p ≥ 0.05, no podemos rechazar la hipótesis nula, lo que sugiere que no hay evidencia suficiente para afirmar que las medias son diferentes.

En nuestro ejemplo, probablemente obtendríamos un valor p muy pequeño, indicando que el grupo de tratamiento tiene una media significativamente mayor que el grupo de control.

Test t para muestras pareadas

El test t pareado se utiliza cuando tenemos mediciones repetidas en los mismos sujetos, como mediciones antes y después de un tratamiento.

# Datos de ejemplo: mediciones antes y después
antes <- c(68, 72, 74, 69, 71, 73, 70, 75)
despues <- c(65, 70, 71, 65, 69, 71, 66, 72)

# Realizamos el test t pareado
resultado_pareado <- t.test(antes, despues, paired = TRUE)

# Mostramos los resultados
resultado_pareado

La diferencia clave aquí es el parámetro paired = TRUE, que indica a R que debe tratar los datos como muestras pareadas.

Test t de una muestra

Este test compara la media de una muestra con un valor de referencia específico.

# Datos de ejemplo
muestra <- c(25.2, 24.8, 26.1, 25.6, 25.9, 25.3, 24.9, 25.7, 25.8, 26.2)

# Valor de referencia (por ejemplo, la media poblacional conocida)
valor_referencia <- 24.5

# Realizamos el test t de una muestra
resultado_una_muestra <- t.test(muestra, mu = valor_referencia)

# Mostramos los resultados
resultado_una_muestra

El parámetro mu especifica el valor de referencia contra el cual queremos comparar nuestra muestra.

Visualización de resultados

Una buena práctica es visualizar los datos junto con los resultados estadísticos:

# Para el ejemplo de muestras independientes
library(ggplot2)

# Preparamos los datos para ggplot
datos <- data.frame(
  valor = c(grupo_control, grupo_tratamiento),
  grupo = factor(rep(c("Control", "Tratamiento"), 
                    c(length(grupo_control), length(grupo_tratamiento))))
)

# Creamos el gráfico
ggplot(datos, aes(x = grupo, y = valor, fill = grupo)) +
  geom_boxplot(alpha = 0.7) +
  geom_jitter(width = 0.2, alpha = 0.5) +
  labs(title = "Comparación entre grupos",
       subtitle = paste("Valor p =", round(resultado$p.value, 4)),
       y = "Valor medido") +
  theme_minimal()

Consideraciones importantes

Al aplicar el test t, debemos tener en cuenta algunas consideraciones:

• Normalidad: El test t asume que los datos siguen una distribución normal. Para muestras grandes (n > 30), el test es bastante robusto ante desviaciones de la normalidad.

• Varianzas iguales o diferentes: Por defecto, R utiliza la variante de Welch del test t, que no asume varianzas iguales. Si queremos asumir varianzas iguales, debemos especificar var.equal = TRUE.

• Tamaño de la muestra: Para muestras muy pequeñas, los resultados del test t pueden no ser confiables.

# Ejemplo con varianzas iguales
resultado_var_iguales <- t.test(grupo_tratamiento, grupo_control, 
                               var.equal = TRUE)

# Comparamos con el resultado anterior (varianzas desiguales)
resultado_var_iguales
resultado  # Resultado con varianzas desiguales

Pruebas unilaterales y bilaterales

Por defecto, el test t realiza una prueba bilateral (two-tailed), que evalúa si las medias son diferentes en cualquier dirección. Sin embargo, a veces queremos probar específicamente si una media es mayor o menor que la otra:

# Test t unilateral: ¿Es la media del grupo de tratamiento mayor?
resultado_unilateral <- t.test(grupo_tratamiento, grupo_control, 
                              alternative = "greater")

# Mostramos los resultados
resultado_unilateral

Las opciones para el parámetro alternative son:

• "two.sided" (bilateral, por defecto)
• "less" (unilateral, primera muestra menor que la segunda)
• "greater" (unilateral, primera muestra mayor que la segunda)

Ejemplo práctico con datos reales

Veamos un ejemplo más realista utilizando un conjunto de datos incorporado en R:

# Cargamos el conjunto de datos de peso de pollos según dieta
data(chickwts)

# Extraemos dos dietas para comparar
dieta1 <- chickwts$weight[chickwts$feed == "sunflower"]
dieta2 <- chickwts$weight[chickwts$feed == "linseed"]

# Realizamos el test t
resultado_dietas <- t.test(dieta1, dieta2)

# Visualizamos los datos y resultados
boxplot(dieta1, dieta2, 
        names = c("Girasol", "Linaza"),
        col = c("yellow", "brown"),
        main = "Peso de pollos según dieta",
        ylab = "Peso (g)")

# Añadimos el valor p al gráfico
mtext(paste("Valor p =", round(resultado_dietas$p.value, 4)), 
      side = 3)

Este ejemplo muestra cómo podemos aplicar el test t para comparar el efecto de diferentes dietas en el peso de los pollos, utilizando datos reales disponibles en R.

Chi-cuadrado para variables categóricas

La prueba de chi-cuadrado es una herramienta estadística fundamental para analizar la relación entre variables categóricas. A diferencia del test t que vimos anteriormente (que se usa para variables numéricas), chi-cuadrado nos permite trabajar con datos que se clasifican en categorías, como género, tipo de tratamiento, o respuestas de tipo "sí/no".

Esta prueba evalúa si existe una asociación significativa entre dos variables categóricas o si una variable categórica sigue una distribución específica. En R, podemos implementar fácilmente diferentes variantes de esta prueba.

Prueba de independencia

La prueba de independencia de chi-cuadrado nos permite determinar si dos variables categóricas están relacionadas. Las hipótesis que evaluamos son:

• Hipótesis nula (H₀): Las variables son independientes (no existe relación)
• Hipótesis alternativa (H₁): Las variables están relacionadas (existe asociación)

Veamos un ejemplo sencillo donde analizamos si existe relación entre el género y la preferencia por un producto:

# Creamos una tabla de contingencia
datos <- matrix(c(45, 30, 25, 40), nrow = 2, 
               dimnames = list(Género = c("Masculino", "Femenino"),
                              Preferencia = c("Producto A", "Producto B")))

# Visualizamos la tabla
datos

Para realizar la prueba chi-cuadrado, utilizamos la función chisq.test():

# Realizamos la prueba chi-cuadrado
resultado <- chisq.test(datos)

# Mostramos los resultados
resultado

La salida nos proporcionará:

• El estadístico chi-cuadrado (X-squared)
• Los grados de libertad (df)
• El valor p (p-value)

Interpretación de resultados

La interpretación se basa principalmente en el valor p:

• Si p < 0.05, rechazamos la hipótesis nula, lo que sugiere que existe una asociación significativa entre las variables.
• Si p ≥ 0.05, no podemos rechazar la hipótesis nula, lo que indica que no hay evidencia suficiente para afirmar que las variables están relacionadas.

Podemos examinar más a fondo los resultados para entender la naturaleza de la asociación:

# Calculamos los valores esperados bajo la hipótesis de independencia
resultado$expected

# Calculamos los residuos estandarizados
resultado$residuals

# Residuos estandarizados ajustados (más útiles para interpretar)
resultado$stdres

Los residuos estandarizados nos indican qué celdas contribuyen más a la significancia del estadístico chi-cuadrado. Valores absolutos mayores a 1.96 sugieren una desviación significativa de lo esperado bajo independencia.

Visualización de tablas de contingencia

Es útil visualizar los datos para complementar el análisis estadístico:

# Gráfico de mosaico
mosaicplot(datos, 
          main = "Relación entre Género y Preferencia",
          color = c("lightblue", "lightgreen"),
          shade = TRUE)

# Usando ggplot2 para una visualización más personalizada
library(ggplot2)
library(reshape2)

# Convertimos la matriz a formato largo
datos_largo <- melt(datos)
colnames(datos_largo) <- c("Género", "Preferencia", "Frecuencia")

# Creamos el gráfico de barras
ggplot(datos_largo, aes(x = Género, y = Frecuencia, fill = Preferencia)) +
  geom_bar(stat = "identity", position = "dodge") +
  labs(title = "Preferencia por producto según género",
       subtitle = paste("Chi-cuadrado p-valor =", round(resultado$p.value, 4))) +
  theme_minimal()

Prueba de bondad de ajuste

Otra aplicación común de chi-cuadrado es la prueba de bondad de ajuste, que evalúa si una variable categórica sigue una distribución específica.

Supongamos que queremos verificar si los nacimientos en los días de la semana siguen una distribución uniforme:

# Datos observados de nacimientos por día de la semana
nacimientos <- c(160, 140, 135, 142, 148, 110, 165)
names(nacimientos) <- c("Lunes", "Martes", "Miércoles", "Jueves", 
                       "Viernes", "Sábado", "Domingo")

# Realizamos la prueba de bondad de ajuste
# Bajo H0, esperamos igual proporción para cada día
resultado_bondad <- chisq.test(nacimientos)

# Mostramos los resultados
resultado_bondad

También podemos especificar una distribución esperada diferente:

# Supongamos que esperamos menos nacimientos en fin de semana
proporciones_esperadas <- c(0.15, 0.15, 0.15, 0.15, 0.15, 0.125, 0.125)

# Realizamos la prueba con proporciones específicas
resultado_prop <- chisq.test(nacimientos, p = proporciones_esperadas)

# Mostramos los resultados
resultado_prop

Corrección de continuidad y pruebas exactas

Cuando trabajamos con muestras pequeñas o tablas con frecuencias bajas, la aproximación chi-cuadrado puede no ser adecuada. R ofrece alternativas:

# Tabla de contingencia pequeña
tabla_pequeña <- matrix(c(5, 1, 2, 8), nrow = 2)

# Chi-cuadrado con corrección de Yates (corrección de continuidad)
chisq.test(tabla_pequeña, correct = TRUE)  # Por defecto es TRUE

# Test exacto de Fisher (alternativa para tablas pequeñas)
fisher.test(tabla_pequeña)

La corrección de Yates ajusta el estadístico chi-cuadrado para tablas 2x2, mientras que el test exacto de Fisher calcula la probabilidad exacta de observar la tabla dada bajo la hipótesis nula.

Ejemplo práctico con datos reales

Veamos un ejemplo utilizando un conjunto de datos incorporado en R:

# Cargamos los datos de supervivencia del Titanic
data(Titanic)

# Convertimos a tabla de contingencia de dos vías (Clase vs Supervivencia)
tabla_titanic <- apply(Titanic, c(1, 4), sum)
tabla_titanic

# Realizamos la prueba chi-cuadrado
resultado_titanic <- chisq.test(tabla_titanic)
resultado_titanic

# Visualizamos los datos
barplot(t(tabla_titanic), 
       beside = TRUE, 
       col = c("darkred", "darkgreen"),
       legend.text = c("No sobrevivió", "Sobrevivió"),
       main = "Supervivencia en el Titanic por clase",
       xlab = "Clase", ylab = "Número de pasajeros")

# Añadimos el valor p al gráfico
mtext(paste("Chi-cuadrado p-valor =", round(resultado_titanic$p.value, 6)), 
     side = 3, line = 0.5)

Este análisis nos permitiría determinar si la clase social (primera, segunda, tercera clase o tripulación) estaba relacionada con la supervivencia en el desastre del Titanic.

Análisis de residuos

Para entender mejor qué categorías contribuyen más a la significancia, podemos analizar los residuos estandarizados:

# Calculamos los residuos estandarizados
residuos <- resultado_titanic$stdres

# Visualizamos los residuos
library(corrplot)
corrplot(residuos, is.corr = FALSE, 
        method = "color", 
        addCoef.col = "black",
        tl.col = "black", 
        cl.lim = c(-5, 5),
        title = "Residuos estandarizados")

Los valores positivos altos indican que hay más casos de los esperados bajo independencia, mientras que los valores negativos altos indican menos casos de los esperados.

Tablas de contingencia multidimensionales

Chi-cuadrado también puede aplicarse a tablas con más de dos dimensiones, aunque la interpretación se vuelve más compleja:

# Analizamos la relación entre clase, género y supervivencia
tabla_3d <- apply(Titanic, c(1, 2, 4), sum)

# Realizamos pruebas separadas para cada género
chisq.test(tabla_3d[, "Male", ])
chisq.test(tabla_3d[, "Female", ])

# Visualizamos con gráficos de mosaico
mosaicplot(tabla_3d, 
          main = "Supervivencia por clase y género",
          color = TRUE,
          shade = TRUE)

Este tipo de análisis nos permite identificar patrones complejos de asociación entre múltiples variables categóricas, como en este caso donde podríamos observar que la relación entre clase y supervivencia podría ser diferente para hombres y mujeres.

Pruebas de normalidad y homocedasticidad

Antes de aplicar muchas pruebas estadísticas paramétricas (como el test t que vimos anteriormente), es necesario verificar que nuestros datos cumplan ciertos supuestos estadísticos. Dos de los más importantes son la normalidad (que los datos sigan una distribución normal) y la homocedasticidad (que los grupos comparados tengan varianzas similares).

R ofrece diversas herramientas para evaluar estos supuestos de forma sencilla y visual. Veamos cómo implementarlas e interpretarlas correctamente.

Pruebas de normalidad

La normalidad es un requisito fundamental para muchos análisis estadísticos. Existen varias formas de evaluar si nuestros datos siguen una distribución normal:

Métodos gráficos

Los métodos visuales son una primera aproximación intuitiva:

# Creamos datos de ejemplo
set.seed(123)  # Para reproducibilidad
datos_normales <- rnorm(100, mean = 50, sd = 10)
datos_sesgados <- rexp(100, rate = 0.1)  # Distribución no normal

# Histograma con curva normal superpuesta
hist(datos_normales, freq = FALSE, 
     main = "Histograma con curva normal",
     col = "lightblue", border = "white")
curve(dnorm(x, mean = mean(datos_normales), sd = sd(datos_normales)), 
      add = TRUE, col = "darkblue", lwd = 2)

El gráfico Q-Q (quantile-quantile) es especialmente útil para evaluar normalidad:

# Gráfico Q-Q para datos normales
qqnorm(datos_normales, main = "Gráfico Q-Q para datos normales")
qqline(datos_normales, col = "red")

# Gráfico Q-Q para datos no normales
qqnorm(datos_sesgados, main = "Gráfico Q-Q para datos no normales")
qqline(datos_sesgados, col = "red")

En un gráfico Q-Q, si los puntos se alinean aproximadamente sobre la línea diagonal, sugiere que los datos siguen una distribución normal. Desviaciones de esta línea, especialmente en los extremos, indican alejamiento de la normalidad.

Pruebas estadísticas de normalidad

Para complementar el análisis visual, podemos utilizar pruebas estadísticas formales:

Test de Shapiro-Wilk

Es una de las pruebas más potentes para evaluar normalidad, especialmente para muestras pequeñas (n < 50):

# Test de Shapiro-Wilk
shapiro.test(datos_normales)
shapiro.test(datos_sesgados)

Interpretación:

• Hipótesis nula: Los datos provienen de una distribución normal
• Hipótesis alternativa: Los datos no provienen de una distribución normal
• Si p < 0.05, rechazamos la normalidad
• Si p ≥ 0.05, no podemos rechazar que los datos sean normales

Test de Kolmogorov-Smirnov

Útil para muestras más grandes:

# Test de Kolmogorov-Smirnov
ks.test(datos_normales, "pnorm", mean = mean(datos_normales), 
        sd = sd(datos_normales))

Test de Anderson-Darling

Disponible en el paquete nortest, es más sensible a las desviaciones en las colas de la distribución:

# Instalamos y cargamos el paquete si es necesario
if (!requireNamespace("nortest", quietly = TRUE)) {
  install.packages("nortest")
}
library(nortest)

# Test de Anderson-Darling
ad.test(datos_normales)

Pruebas de homocedasticidad

La homocedasticidad (igualdad de varianzas) es otro supuesto importante cuando comparamos grupos. Veamos cómo evaluarla:

Visualización de la homocedasticidad

Los diagramas de caja son útiles para comparar visualmente la dispersión entre grupos:

# Creamos datos de ejemplo para dos grupos
grupo1 <- rnorm(50, mean = 10, sd = 2)
grupo2 <- rnorm(50, mean = 12, sd = 2)  # Misma varianza
grupo3 <- rnorm(50, mean = 14, sd = 4)  # Mayor varianza

# Combinamos los datos
datos_combinados <- data.frame(
  valor = c(grupo1, grupo2, grupo3),
  grupo = factor(rep(c("Grupo 1", "Grupo 2", "Grupo 3"), 
                    each = 50))
)

# Creamos el diagrama de cajas
boxplot(valor ~ grupo, data = datos_combinados,
        col = c("lightblue", "lightgreen", "lightpink"),
        main = "Comparación de variabilidad entre grupos")

Prueba de Levene

La prueba de Levene es robusta y no requiere normalidad en los datos:

# Instalamos y cargamos el paquete car si es necesario
if (!requireNamespace("car", quietly = TRUE)) {
  install.packages("car")
}
library(car)

# Prueba de Levene
leveneTest(valor ~ grupo, data = datos_combinados)

Interpretación:

• Hipótesis nula: Las varianzas son iguales entre los grupos
• Hipótesis alternativa: Al menos un grupo tiene varianza diferente
• Si p < 0.05, rechazamos la homocedasticidad
• Si p ≥ 0.05, no podemos rechazar que las varianzas sean iguales

Prueba de Bartlett

La prueba de Bartlett es otra opción, aunque es más sensible a desviaciones de la normalidad:

# Prueba de Bartlett
bartlett.test(valor ~ grupo, data = datos_combinados)

Esta prueba es más potente que la de Levene cuando los datos son normales, pero menos robusta cuando no lo son.

Prueba F para comparar dos varianzas

Cuando solo comparamos dos grupos, podemos usar la prueba F:

# Prueba F para comparar varianzas de dos grupos
var.test(grupo1, grupo2)

Aplicación práctica con datos reales

Veamos un ejemplo completo utilizando un conjunto de datos incorporado en R:

# Cargamos el conjunto de datos iris
data(iris)

# Verificamos la normalidad para la longitud del sépalo por especie
par(mfrow = c(2, 2))  # Dividimos la ventana gráfica en 2x2

# Histogramas por especie
for (sp in unique(iris$Species)) {
  datos_especie <- iris$Sepal.Length[iris$Species == sp]
  hist(datos_especie, main = paste("Histograma para", sp),
       xlab = "Longitud del sépalo", col = "lightblue")
}

# Gráficos Q-Q por especie
par(mfrow = c(2, 2))
for (sp in unique(iris$Species)) {
  datos_especie <- iris$Sepal.Length[iris$Species == sp]
  qqnorm(datos_especie, main = paste("Q-Q plot para", sp))
  qqline(datos_especie, col = "red")
}

# Pruebas de Shapiro-Wilk por especie
for (sp in unique(iris$Species)) {
  datos_especie <- iris$Sepal.Length[iris$Species == sp]
  cat("Shapiro-Wilk para", sp, ":\n")
  print(shapiro.test(datos_especie))
}

# Verificamos homocedasticidad entre especies
# Diagrama de cajas
par(mfrow = c(1, 1))
boxplot(Sepal.Length ~ Species, data = iris,
        col = c("lightblue", "lightgreen", "lightpink"),
        main = "Longitud del sépalo por especie")

# Prueba de Levene
leveneTest(Sepal.Length ~ Species, data = iris)

Transformaciones para corregir violaciones de supuestos

Cuando los datos no cumplen con los supuestos, podemos aplicar transformaciones para intentar corregir el problema:

# Datos con distribución sesgada
datos_sesgados <- rexp(100, 0.1)

# Verificamos normalidad antes de transformar
shapiro.test(datos_sesgados)

# Aplicamos transformación logarítmica
datos_log <- log(datos_sesgados)

# Verificamos normalidad después de transformar
shapiro.test(datos_log)

# Comparamos gráficamente
par(mfrow = c(1, 2))
hist(datos_sesgados, main = "Datos originales", col = "lightblue")
hist(datos_log, main = "Datos transformados (log)", col = "lightgreen")

# Comparamos gráficos Q-Q
par(mfrow = c(1, 2))
qqnorm(datos_sesgados, main = "Q-Q original")
qqline(datos_sesgados, col = "red")
qqnorm(datos_log, main = "Q-Q transformado (log)")
qqline(datos_log, col = "red")

Las transformaciones comunes incluyen:

• Logarítmica: útil para datos con asimetría positiva
• Raíz cuadrada: menos drástica que la logarítmica
• Box-Cox: familia de transformaciones que incluye las anteriores

Alternativas no paramétricas

Cuando las transformaciones no resuelven los problemas o no son apropiadas, podemos recurrir a pruebas no paramétricas que no requieren normalidad:

# En lugar de t.test, podemos usar wilcox.test (prueba de Wilcoxon)
wilcox.test(Sepal.Length ~ Species, 
           data = iris[iris$Species %in% c("setosa", "versicolor"),])

# En lugar de ANOVA, podemos usar kruskal.test (Kruskal-Wallis)
kruskal.test(Sepal.Length ~ Species, data = iris)

Recomendaciones prácticas

Al verificar los supuestos estadísticos, ten en cuenta estas recomendaciones:

• Combina siempre métodos visuales con pruebas formales
• Para muestras grandes (n > 50), las pruebas de normalidad pueden ser demasiado sensibles
• El teorema del límite central sugiere que para muestras grandes, muchas pruebas son robustas ante desviaciones moderadas de la normalidad
• La homocedasticidad suele ser más importante que la normalidad para pruebas como ANOVA
• Cuando tengas dudas, opta por métodos no paramétricos o correcciones como la de Welch para el test t

# Ejemplo de test t con corrección de Welch (no asume varianzas iguales)
t.test(Sepal.Length ~ Species, 
      data = iris[iris$Species %in% c("setosa", "versicolor"),],
      var.equal = FALSE)  # Esta es la opción por defecto

Estas pruebas de supuestos son fundamentales para garantizar la validez de nuestros análisis estadísticos y evitar conclusiones erróneas basadas en técnicas inapropiadas para nuestros datos.