Media, mediana y desviación estándar en R

Media y mediana

Los dos resúmenes más usados de toda la estadística son una llamada de función cada uno:

mpg <- mtcars$mpg

print(mean(mpg))
print(median(mpg))

mean() suma todo y divide entre el número de valores. median() ordena los valores y elige el del medio (o promedia los dos centrales cuando el recuento es par). Para los 32 coches de mtcars, el consumo medio es de unas 20.1 mpg y la mediana es 19.2 - valores cercanos, lo que indica que los datos no están muy sesgados. Cuando discrepan mucho, eso también es información; lo veremos al final.

Algo que sorprende a todo el mundo: si el vector contiene aunque sea un solo NA, ambas funciones devuelven NA:

x <- c(4, 8, NA, 6)

print(mean(x))
print(mean(x, na.rm = TRUE))

Es deliberado - R se niega a fingir en silencio que el valor faltante no existe. na.rm = TRUE significa "calcula sobre los valores que tienes". Casi todas las funciones de resumen de este artículo lo aceptan. La historia completa de cómo se propaga NA está en valores faltantes.

Desviación estándar y varianza

sd() mide la dispersión: a grandes rasgos, a qué distancia de la media se sitúa un valor típico, en las mismas unidades que los datos. var() es su cuadrado:

mpg <- mtcars$mpg

print(sd(mpg))
print(var(mpg))
print(sd(mpg)^2)

Una desviación estándar de unas 6 mpg significa que los coches suelen situarse a unas 6 mpg del promedio de 20.1. Como sd() está en las unidades propias de los datos, es la que se reporta; var() aparece sobre todo dentro de otras fórmulas.

Aquí viene el detalle que importa para los trabajos de clase: sd() y var() calculan el estadístico muestral - dividen la suma de desviaciones al cuadrado entre n − 1, no entre n:

x <- c(2, 4, 4, 4, 5, 5, 7, 9)
n <- length(x)

sample_var <- sum((x - mean(x))^2) / (n - 1)
population_var <- sum((x - mean(x))^2) / n

print(sample_var)      # what var(x) gives you
print(var(x))
print(population_var)  # divide by n instead

¿Por qué n − 1? Como estimaste la media a partir de los mismos datos, las desviaciones alrededor de esa media estimada son sistemáticamente un poco demasiado pequeñas; dividir entre n − 1 lo corrige. Dado que tus datos casi siempre son una muestra de algo más grande, la versión con n − 1 es la que quieres. Si de verdad tienes la población completa (cada estudiante de la clase, cada producto del catálogo), multiplica: var(x) * (n - 1) / n.

Error estándar de la media

La desviación estándar y el error estándar se confunden constantemente, así que mantenlos separados: la sd describe los datos, el SE describe tu estimación de la media. R no tiene una se() incorporada, pero la fórmula es una línea:

mpg <- mtcars$mpg

se <- sd(mpg) / sqrt(length(mpg))
print(se)

El error estándar se reduce a medida que crece la muestra - recoge cuatro veces más datos y el SE se reduce a la mitad - porque una muestra más grande fija la media con más precisión. La sd no se reduce con el tamaño de la muestra; los coches son tan variados como son sin importar cuántos midas. El SE es el bloque de construcción de los intervalos de confianza, que es donde se gana el sueldo.

summary() - la visión general en una llamada

summary() te da el resumen de cinco números más la media de una sola vez, y funciona sobre data frames completos:

print(summary(mtcars$mpg))

print(summary(mtcars[, c("mpg", "wt", "hp")]))

Llamado sobre un data frame, resume cada columna - las columnas numéricas reciben mín/cuartiles/media/máx, y los factores reciben recuentos por nivel. Es lo primero que hay que ejecutar sobre cualquier conjunto de datos recién cargado: los valores imposibles (una edad negativa, un máximo de 9999) saltan a la vista de inmediato.

Cuantiles, rango e IQR

quantile() generaliza la mediana a cualquier punto de corte:

mpg <- mtcars$mpg

print(quantile(mpg))
print(quantile(mpg, probs = c(0.1, 0.5, 0.9)))

print(range(mpg))
print(IQR(mpg))
print(max(mpg) - min(mpg))

Sin argumentos devuelve el mínimo, los cuartiles y el máximo. Pasa probs = para puntos de corte concretos - los percentiles 10 y 90 de arriba delimitan dónde vive el grueso de los datos. IQR() (la distancia entre los percentiles 25 y 75) es una medida de dispersión que, a diferencia de sd(), no se deja arrastrar por los valores atípicos. range() devuelve el mínimo y el máximo como una pareja.

La moda: mode() en R NO hace esto

Esta atrapa a todo el mundo exactamente una vez. R tiene una función llamada mode(), y no tiene nada que ver con la estadística - reporta el tipo de almacenamiento de un objeto:

x <- c(3, 7, 7, 2, 7, 5, 3)

print(mode(x))   # "numeric" - not what you wanted

counts <- table(x)
print(counts)
print(names(which.max(counts)))

El idiom que hay que recordar: table(x) cuenta cuántas veces aparece cada valor, which.max() encuentra el recuento más grande y names() extrae el valor en sí. Fíjate en que vuelve como una cadena de caracteres (los nombres de una tabla siempre lo son); envuélvelo en as.numeric() si necesitas calcular con él. Si dos valores empatan, which.max() devuelve en silencio solo el primero - revisa la tabla tú mismo cuando los empates sean plausibles.

Media vs mediana: cuál reportar

La media usa todos los valores, lo cual es su fortaleza y su debilidad - un solo valor extremo la arrastra. A la mediana solo le importa el centro, así que los valores atípicos apenas la tocan:

salaries <- c(42000, 45000, 47000, 48000, 51000, 54000, 58000)
print(mean(salaries))
print(median(salaries))

# The CEO joins the dataset:
salaries <- c(salaries, 900000)
print(mean(salaries))
print(median(salaries))

Un solo valor añadido empuja la media de unos 49,300 a más de 155,000 - un número que no describe a nadie en los datos - mientras que la mediana solo se mueve de 48,000 a 49,500. Por eso los ingresos, los precios de la vivienda y las estancias hospitalarias se reportan como medianas: los datos sesgados con una cola larga hacen que la media resulte engañosa. Para datos aproximadamente simétricos las dos coinciden y la media está bien (y es estadísticamente más eficiente). Un histograma rápido te dice en qué situación estás - y comparar tu media con tu mediana es en sí mismo una comprobación de sesgo de una línea.

Lo que te llevas

• mean(x) y median(x) - añade na.rm = TRUE cuando haya valores faltantes.
• sd(x) y var(x) calculan el estadístico muestral (el denominador n − 1) - que es lo que quieres.
• El error estándar es sd(x) / sqrt(length(x)) - mide qué tan bien conoces la media, no cuán dispersos están los datos.
• summary() sobre un data frame recién cargado es la comprobación de cordura más rápida de R.
• La moda es names(which.max(table(x))) - mode() en R trata sobre tipos de almacenamiento.
• Datos sesgados o con valores atípicos: reporta la mediana. Datos simétricos: la media está bien.

A continuación: medir cómo dos variables se mueven juntas - correlación con cor() y cor.test().