Estadística descriptiva con NumPy

Medidas de tendencia central

Las medidas de tendencia central describen el valor representativo de un conjunto de datos. NumPy proporciona funciones directas para calcular las más habituales: media y mediana. Para la moda se puede recurrir a una combinación de funciones auxiliares.

Media aritmética con np.mean

La media es la suma de todos los valores dividida por el número de elementos. Se calcula con np.mean o con el método .mean() del array.

import numpy as np

# Temperaturas diarias en grados Celsius durante una semana
temperaturas = np.array([22.1, 24.5, 19.8, 23.0, 25.2, 21.7, 20.3])

media = np.mean(temperaturas)
print(media)  # 22.37...

Mediana con np.median

La mediana es el valor que ocupa la posición central cuando los datos están ordenados. Es más robusta frente a valores atípicos que la media.

salarios = np.array([1800, 2100, 2300, 2500, 15000])

print(np.mean(salarios))    # 4740.0 (distorsionada por 15000)
print(np.median(salarios))  # 2300.0 (refleja mejor el centro)

Moda con np.unique y np.bincount

NumPy no incluye una función de moda directa, pero se puede obtener combinando np.unique o np.bincount. La función np.bincount cuenta las apariciones de cada entero no negativo en un array.

# Puntuaciones de una encuesta (valores enteros de 1 a 5)
puntuaciones = np.array([3, 5, 3, 4, 3, 2, 5, 4, 3, 1])

conteo = np.bincount(puntuaciones)
moda = np.argmax(conteo)
print(moda)  # 3 (aparece 4 veces)

Para datos no enteros se puede usar np.unique con return_counts:

alturas = np.array([1.70, 1.75, 1.70, 1.80, 1.75, 1.70])

valores, conteos = np.unique(alturas, return_counts=True)
moda = valores[np.argmax(conteos)]
print(moda)  # 1.70

La media es sensible a valores extremos, mientras que la mediana los ignora. Elegir una u otra depende de la distribución de los datos y del objetivo del análisis.

Medidas de dispersión

Las medidas de dispersión cuantifican la variabilidad de los datos respecto a su centro. Complementan a las medidas de tendencia central y son imprescindibles para interpretar correctamente un dataset.

Varianza con np.var

La varianza mide la dispersión promedio de los datos respecto a la media. Por defecto np.var calcula la varianza poblacional (divide entre $N$). Para obtener la varianza muestral (divide entre $N-1$) se utiliza el parámetro ddof=1.

import numpy as np

notas = np.array([7.5, 8.0, 6.5, 9.0, 7.0, 8.5])

var_poblacional = np.var(notas)
var_muestral = np.var(notas, ddof=1)

print(f"Poblacional: {var_poblacional:.4f}")  # 0.7222
print(f"Muestral:    {var_muestral:.4f}")     # 0.8667

Desviación estándar con np.std

La desviación estándar es la raíz cuadrada de la varianza y tiene la ventaja de expresarse en las mismas unidades que los datos originales. Admite el mismo parámetro ddof.

print(f"Desviación poblacional: {np.std(notas):.4f}")        # 0.8498
print(f"Desviación muestral:    {np.std(notas, ddof=1):.4f}")  # 0.9309

A partir de NumPy 2.0, tanto np.var como np.std aceptan el parámetro mean para pasar una media precalculada, lo que evita recalcularla internamente y mejora el rendimiento cuando la media ya se ha obtenido previamente.

media_notas = np.mean(notas)
var_con_media = np.var(notas, ddof=1, mean=media_notas)
std_con_media = np.std(notas, ddof=1, mean=media_notas)

print(f"Varianza (con mean=): {var_con_media:.4f}")   # 0.8667
print(f"Desv. std (con mean=): {std_con_media:.4f}")   # 0.9309

El parámetro ddof (delta degrees of freedom) indica cuántos grados de libertad se restan del denominador. Con ddof=0 se obtiene la versión poblacional y con ddof=1 la muestral, que es la habitual en inferencia estadística.

Rango

El rango es la diferencia entre el valor máximo y el mínimo. En NumPy 2 el método ptp fue eliminado del ndarray, por lo que se calcula de forma directa con np.max y np.min.

datos = np.array([12, 45, 23, 67, 34, 89, 11])

rango = np.max(datos) - np.min(datos)
print(rango)  # 78 (89 - 11)

Rango intercuartílico (IQR)

El IQR mide la dispersión del 50% central de los datos y es resistente a valores atípicos. Se calcula como la diferencia entre el percentil 75 y el percentil 25.

edades = np.array([18, 22, 25, 27, 30, 32, 35, 40, 45, 60, 75])

q1 = np.percentile(edades, 25)
q3 = np.percentile(edades, 75)
iqr = q3 - q1

print(f"Q1: {q1}, Q3: {q3}, IQR: {iqr}")  # Q1: 26.0, Q3: 42.5, IQR: 16.5

Percentiles y cuantiles

Los percentiles dividen los datos ordenados en 100 partes iguales, mientras que los cuantiles los dividen en fracciones entre 0 y 1. Ambos son fundamentales para comprender la distribución de un dataset.

Percentiles con np.percentile

import numpy as np

tiempos_respuesta = np.array([120, 150, 180, 200, 210, 230, 250, 300, 350, 500])

p50 = np.percentile(tiempos_respuesta, 50)   # Mediana
p90 = np.percentile(tiempos_respuesta, 90)   # Percentil 90
p99 = np.percentile(tiempos_respuesta, 99)   # Percentil 99

print(f"P50: {p50} ms, P90: {p90} ms, P99: {p99} ms")

Se pueden calcular varios percentiles a la vez pasando una lista:

percentiles = np.percentile(tiempos_respuesta, [10, 25, 50, 75, 90])
print(percentiles)  # Array con los 5 valores

Cuantiles con np.quantile

La función np.quantile trabaja con fracciones entre 0 y 1. Es equivalente a np.percentile pero con una escala normalizada.

# Cuartiles (divide en 4 partes)
cuartiles = np.quantile(tiempos_respuesta, [0.25, 0.5, 0.75])
print(cuartiles)

Los percentiles P90, P95 y P99 son especialmente útiles en sistemas de producción para medir latencias y tiempos de respuesta, ya que reflejan el comportamiento del servicio en los casos peor atendidos.

Correlación y covarianza

Cuando se trabaja con múltiples variables, la correlación y la covarianza permiten cuantificar cómo se relacionan entre sí.

Covarianza con np.cov

La covarianza mide la variación conjunta de dos variables. Un valor positivo indica que tienden a crecer juntas, un valor negativo indica relación inversa.

import numpy as np

# Horas de estudio y calificaciones de 8 alumnos
horas_estudio = np.array([2, 3, 5, 6, 7, 8, 9, 10])
calificaciones = np.array([4.5, 5.0, 6.5, 7.0, 7.5, 8.0, 9.0, 9.5])

# Matriz de covarianza 2x2
cov_matrix = np.cov(horas_estudio, calificaciones)
print(cov_matrix)
# [[8.2857  4.8571]
#  [4.8571  2.9643]]

La diagonal de la matriz contiene las varianzas de cada variable. Los elementos fuera de la diagonal son las covarianzas cruzadas.

Coeficiente de correlación de Pearson con np.corrcoef

El coeficiente de correlación normaliza la covarianza al rango $[-1, 1]$, lo que facilita la interpretación.

corr_matrix = np.corrcoef(horas_estudio, calificaciones)
print(corr_matrix)
# [[1.     0.9803]
#  [0.9803 1.    ]]

r = corr_matrix[0, 1]
print(f"Correlación: {r:.4f}")  # 0.9803 (correlación fuerte positiva)

flowchart LR
  subgraph Interpretación[Coeficiente de correlación r]
    A["r cercano a +1<br>Relación lineal<br>positiva fuerte"]
    B["r cercano a 0<br>Sin relación<br>lineal clara"]
    C["r cercano a -1<br>Relación lineal<br>negativa fuerte"]
  end

Ejemplo con múltiples variables

Cuando se tienen más de dos variables, np.corrcoef genera una matriz cuadrada con todas las correlaciones por pares.

# Tres indicadores económicos de 6 países
pib = np.array([45000, 38000, 52000, 41000, 60000, 35000])
esperanza_vida = np.array([82, 79, 84, 80, 85, 77])
educacion = np.array([12.5, 11.0, 13.2, 11.8, 14.0, 10.5])

datos = np.array([pib, esperanza_vida, educacion])
corr = np.corrcoef(datos)
print(np.round(corr, 3))

Histogramas para análisis exploratorio

NumPy permite calcular histogramas sin necesidad de bibliotecas de visualización. Esto es útil para analizar la distribución de frecuencias de los datos de forma programática.

Histograma unidimensional con np.histogram

La función np.histogram devuelve dos arrays: los conteos por intervalo y los bordes de los intervalos (bins).

import numpy as np

# Pesos en kg de 50 personas (datos simulados)
rng = np.random.default_rng(42)
pesos = rng.normal(loc=70, scale=12, size=50)

conteos, bordes = np.histogram(pesos, bins=5)
print("Conteos:", conteos)
print("Bordes:", bordes)

Se pueden especificar los bordes de forma explícita:

# Intervalos personalizados
bordes_custom = [40, 55, 65, 75, 85, 100]
conteos, _ = np.histogram(pesos, bins=bordes_custom)

for i in range(len(conteos)):
    print(f"  [{bordes_custom[i]}-{bordes_custom[i+1]}): {conteos[i]}")

Histograma bidimensional con np.histogram2d

Cuando se quiere analizar la distribución conjunta de dos variables, np.histogram2d genera una matriz de conteos en una cuadrícula 2D.

# Altura (cm) y peso (kg) de 200 personas
rng = np.random.default_rng(0)
alturas = rng.normal(170, 10, 200)
pesos_2d = rng.normal(70, 12, 200)

conteos_2d, bordes_x, bordes_y = np.histogram2d(alturas, pesos_2d, bins=5)
print("Forma de la matriz:", conteos_2d.shape)  # (5, 5)
print(conteos_2d)

La función np.histogram solo calcula conteos y bordes, no genera gráficos. Para visualizar los resultados se combinan con bibliotecas como Matplotlib o Plotly.

Operaciones estadísticas sobre ejes

Todas las funciones estadísticas de NumPy aceptan el parámetro axis, que permite calcular estadísticas a lo largo de una dimensión concreta de un array multidimensional. Este mecanismo es fundamental cuando se trabaja con datasets organizados en filas y columnas.

import numpy as np

# Ventas trimestrales de 4 productos (filas: productos, columnas: trimestres)
ventas = np.array([
    [1200, 1350, 1100, 1500],  # Producto A
    [800,  950,  870,  1020],  # Producto B
    [2100, 2300, 2050, 2400],  # Producto C
    [600,  720,  680,  750]    # Producto D
])

# Media por producto (a lo largo de columnas, axis=1)
media_producto = np.mean(ventas, axis=1)
print("Media por producto:", media_producto)

# Media por trimestre (a lo largo de filas, axis=0)
media_trimestre = np.mean(ventas, axis=0)
print("Media por trimestre:", media_trimestre)

El parámetro axis=0 colapsa las filas (calcula sobre cada columna) y axis=1 colapsa las columnas (calcula sobre cada fila). Este mismo parámetro funciona de forma consistente en np.mean, np.std, np.var, np.median, np.percentile, np.min, np.max y np.sum.

# Desviación estándar por producto
std_producto = np.std(ventas, axis=1, ddof=1)
print("Desviación por producto:", np.round(std_producto, 1))

# Percentil 75 por trimestre
p75_trimestre = np.percentile(ventas, 75, axis=0)
print("P75 por trimestre:", p75_trimestre)

Caso práctico: análisis de sensor

Un array tridimensional puede representar lecturas de múltiples sensores a lo largo de varios días y horas. Las operaciones sobre ejes permiten resumir la información sin bucles explícitos.

# Forma: (3 sensores, 7 días, 24 horas)
rng = np.random.default_rng(99)
lecturas = rng.normal(loc=25, scale=3, size=(3, 7, 24))

# Temperatura media diaria por sensor (colapsar eje de horas)
media_diaria = np.mean(lecturas, axis=2)
print("Forma resultado:", media_diaria.shape)  # (3, 7)

# Temperatura máxima absoluta por sensor (colapsar días y horas)
max_sensor = np.max(lecturas, axis=(1, 2))
print("Máxima por sensor:", np.round(max_sensor, 2))

El parámetro axis acepta una tupla de ejes para colapsar varias dimensiones a la vez, como axis=(1, 2). Omitir axis o pasar axis=None calcula la estadística sobre todos los elementos del array.