← Volver al blog

Cómo identificar una distribución estadística usando Python

22 de abril de 2026 7 min de lectura statistics distributions python scipy metallurgy data analysis

En ingeniería, minería y análisis de procesos trabajamos constantemente con datos que no siempre se comportan como una distribución normal. Las recuperaciones metalúrgicas, los tiempos entre fallas de un equipo, los tamaños de partícula o los consumos de reactivo obedecen a familias de distribuciones distintas — lognormal, Weibull, gamma, beta — y asumir normalidad cuando no corresponde puede invalidar intervalos de confianza, pruebas de hipótesis y simulaciones Monte Carlo.

En este artículo mostraré un flujo reproducible en Python para identificar la distribución estadística que mejor describe un conjunto de datos, combinando exploración visual, ajuste por máxima verosimilitud, pruebas de bondad de ajuste y criterios de información.

1. ¿Por qué importa identificar la distribución correcta?

Identificar la familia de distribución subyacente tiene implicancias prácticas:

Pruebas de hipótesis: la prueba t, ANOVA o intervalos clásicos asumen normalidad.
Control de procesos: los límites 3σ pierden sentido si la variable es asimétrica.
Simulación Monte Carlo: muestrear de una distribución equivocada sesga todas las conclusiones.
Confiabilidad: el tiempo a falla de bombas, motores o revestimientos suele ser Weibull, no normal.
Granulometría / química: los tamaños de partícula y las leyes de mineral suelen ser lognormales.

El flujo general es:

Exploración visual  →  Candidatos razonables  →  Ajuste MLE  →
Bondad de ajuste  →  Criterios de información  →  Validación visual

2. Librerías y datos de ejemplo

Usaremos el stack científico estándar:

import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
from scipy import stats
import warnings
warnings.filterwarnings("ignore")

rng = np.random.default_rng(seed=42)

Como ejemplo generaremos una muestra sintética que simule recuperaciones de cobre (%) de pruebas de flotación de laboratorio — típicamente asimétricas y acotadas en [0, 100]:

# Muestra simulada: recuperación Cu en un circuito rougher
datos = rng.lognormal(mean=np.log(85), sigma=0.08, size=200)
datos = np.clip(datos, 0, 100)

En un caso real reemplazarías esta línea por pd.read_excel("tus_datos.xlsx")[columna].values.

3. Exploración inicial

Antes de ajustar nada, mira los datos. Calcula estadísticos descriptivos y grafica.

3.1 Estadísticos descriptivos

desc = pd.Series(datos).describe()
print(desc)
print(f"Asimetría (skewness): {stats.skew(datos):.3f}")
print(f"Curtosis (exceso):    {stats.kurtosis(datos):.3f}")

Interpretación rápida:

Indicador	Distribución sugerida
Skewness ≈ 0, kurtosis ≈ 0	Normal
Skewness > 0, cola derecha	Lognormal, gamma, Weibull
Skewness < 0, cola izquierda	Beta (β > α), Weibull invertida
Datos acotados en [0, 1] o [a, b]	Beta, uniforme
Datos enteros de conteo	Poisson, binomial negativa

3.2 Histograma y KDE

fig, ax = plt.subplots(1, 2, figsize=(12, 4))

ax[0].hist(datos, bins=25, density=True, alpha=0.6, edgecolor="black")
ax[0].set_title("Histograma")
ax[0].set_xlabel("Recuperación Cu (%)")

pd.Series(datos).plot.kde(ax=ax[1])
ax[1].set_title("Densidad estimada (KDE)")
plt.tight_layout()
plt.show()

3.3 Q-Q plot contra la normal

Un primer diagnóstico barato:

stats.probplot(datos, dist="norm", plot=plt)
plt.title("Q-Q plot vs Normal")
plt.show()

Si los puntos caen sobre la recta, hay buen indicio de normalidad. Si se curvan en los extremos, la cola es más pesada o más liviana que lo normal.

4. Candidatos razonables según el dominio

No tiene sentido probar las ~80 distribuciones de SciPy. Se eligen candidatos plausibles según el soporte de la variable y el fenómeno físico:

Variable	Candidatos típicos
Recuperación, ley, % sólidos	Beta, lognormal, gamma
Tiempo a falla	Weibull, exponencial, lognormal
Tamaño de partícula	Lognormal, Weibull (Rosin-Rammler)
Error de medición	Normal, t de Student
Conteos (partículas, eventos)	Poisson, binomial negativa

Para nuestro ejemplo consideraremos cuatro candidatos continuos:

candidatos = {
    "Normal":     stats.norm,
    "Lognormal":  stats.lognorm,
    "Gamma":      stats.gamma,
    "Weibull":    stats.weibull_min,
}

5. Ajuste por máxima verosimilitud (MLE)

SciPy ajusta los parámetros de cualquier distribución continua con .fit(), que resuelve el problema de máxima verosimilitud:

\[\hat{\theta} = \arg\max_{\theta} \prod_{i=1}^{n} f(x_i \mid \theta)\]

ajustes = {}
for nombre, dist in candidatos.items():
    params = dist.fit(datos)
    ajustes[nombre] = {"dist": dist, "params": params}
    print(f"{nombre:12s} → parámetros: {np.round(params, 4)}")

Nota: en SciPy, los dos últimos parámetros siempre son loc y scale. Los parámetros de forma (shape) van primero. Para lognorm el primer parámetro es s (sigma de la lognormal); para weibull_min es c (parámetro de forma k de Weibull).

6. Pruebas de bondad de ajuste

Tres pruebas clásicas, en orden de uso:

6.1 Kolmogorov–Smirnov (KS)

Compara la CDF empírica con la teórica. Útil pero poco sensible a las colas.

for nombre, info in ajustes.items():
    dist, params = info["dist"], info["params"]
    D, p = stats.kstest(datos, dist.name, args=params)
    print(f"{nombre:12s}  D={D:.4f}  p={p:.4f}")

Regla práctica: si p > 0.05, no rechazas la hipótesis de que los datos provienen de esa distribución.

6.2 Anderson–Darling

Más sensible a las colas que KS. SciPy la implementa para normal, exponencial, logística, Gumbel y Weibull:

result = stats.anderson(datos, dist="norm")
print(f"Estadístico AD: {result.statistic:.4f}")
for sl, cv in zip(result.significance_level, result.critical_values):
    print(f"  Nivel {sl:.1f}%  valor crítico {cv:.3f}")

6.3 Shapiro–Wilk (solo normalidad)

Potente con muestras pequeñas a medianas:

W, p = stats.shapiro(datos)
print(f"Shapiro-Wilk: W={W:.4f}  p={p:.4f}")

Advertencia: con n muy grande (por ejemplo n > 2000), casi cualquier prueba rechaza normalidad por diferencias triviales. En ese régimen, privilegia las herramientas gráficas y los criterios de información.

7. Selección del mejor modelo: AIC y BIC

Las pruebas de bondad de ajuste te dicen si un modelo es compatible con los datos, no cuál es el mejor. Para eso se usan criterios de información que penalizan el número de parámetros:

\[\mathrm{AIC} = 2k - 2\ln(\hat{L})\]

\[\mathrm{BIC} = k \ln(n) - 2 \ln(\hat{L})\]

Donde k es el número de parámetros y L la verosimilitud. Menor es mejor. BIC penaliza más fuerte la complejidad y tiende a preferir modelos más simples.

def aic_bic(datos, dist, params):
    logL = np.sum(dist.logpdf(datos, *params))
    k = len(params)
    n = len(datos)
    aic = 2 * k - 2 * logL
    bic = k * np.log(n) - 2 * logL
    return aic, bic, logL

resumen = []
for nombre, info in ajustes.items():
    aic, bic, logL = aic_bic(datos, info["dist"], info["params"])
    resumen.append({"Distribución": nombre, "logL": logL, "AIC": aic, "BIC": bic})

df = pd.DataFrame(resumen).sort_values("AIC")
print(df.to_string(index=False))

El modelo en la fila superior tiene el mejor compromiso ajuste/complejidad.

8. Validación visual: el paso que nunca se debe saltar

Dos gráficos finales confirman (o desmienten) la elección: histograma con densidades superpuestas + Q-Q plot del mejor candidato.

Ajustes candidatos + Q-Q plot

A la izquierda el histograma con los 4 candidatos superpuestos: la curva roja (Lognormal) sigue mucho mejor la forma empírica que la azul (Normal). A la derecha, el Q-Q plot contra la Lognormal con \(R^2 \approx 1\) confirma que los cuantiles observados caen sobre la diagonal teórica en todo el rango.

x = np.linspace(datos.min(), datos.max(), 500)
plt.hist(datos, bins=25, density=True, alpha=0.5, edgecolor="black", label="Datos")
for nombre, info in ajustes.items():
    plt.plot(x, info["dist"].pdf(x, *info["params"]), label=nombre, lw=2)
plt.legend()
plt.xlabel("Recuperación Cu (%)")
plt.ylabel("Densidad")
plt.title("Ajustes candidatos")
plt.show()

# Q-Q plot del mejor
mejor = df.iloc[0]["Distribución"]
dist, params = ajustes[mejor]["dist"], ajustes[mejor]["params"]
stats.probplot(datos, dist=dist, sparams=params, plot=plt)
plt.title(f"Q-Q plot vs {mejor}")
plt.show()

9. Flujo automatizado

Para reutilizar el análisis, una función que toma el arreglo de datos y una lista de candidatos:

def identificar_distribucion(datos, candidatos=None):
    if candidatos is None:
        candidatos = {
            "Normal":    stats.norm,
            "Lognormal": stats.lognorm,
            "Gamma":     stats.gamma,
            "Weibull":   stats.weibull_min,
            "Beta":      stats.beta,
            "Expon":     stats.expon,
        }

    filas = []
    for nombre, dist in candidatos.items():
        try:
            params = dist.fit(datos)
            logL = np.sum(dist.logpdf(datos, *params))
            k, n = len(params), len(datos)
            aic = 2 * k - 2 * logL
            bic = k * np.log(n) - 2 * logL
            D, p_ks = stats.kstest(datos, dist.name, args=params)
            filas.append({
                "Distribución": nombre,
                "AIC": aic,
                "BIC": bic,
                "KS_D": D,
                "KS_pvalue": p_ks,
                "params": params,
            })
        except Exception as e:
            print(f"{nombre}: no se pudo ajustar ({e})")

    return pd.DataFrame(filas).sort_values("AIC").reset_index(drop=True)

resultados = identificar_distribucion(datos)
print(resultados[["Distribución", "AIC", "BIC", "KS_D", "KS_pvalue"]])

10. Recomendaciones finales

Empieza por el fenómeno, no por los datos. Si mides un tiempo a falla, parte probando Weibull y exponencial, no normal.
No confíes en una sola prueba. Combina KS, AD o Shapiro con AIC/BIC y gráficos Q-Q.
Cuidado con muestras grandes. Con n > 2000, las pruebas rechazan normalidad por diferencias irrelevantes; privilegia lo visual.
Documenta los parámetros ajustados. Anótalos en el informe junto con el tamaño de muestra y la versión de SciPy — tu yo futuro lo agradecerá.
Considera transformar antes que forzar el ajuste. Datos lognormales ajustados en escala log suelen ser aproximadamente normales y permiten usar herramientas estándar.
Si nada ajusta, la distribución empírica (ECDF) también es una respuesta válida. No todo fenómeno cabe en una familia paramétrica.

Cierre

Identificar una distribución no es un ejercicio académico: condiciona la validez de todo lo que venga después — intervalos de confianza, simulaciones, cartas de control, modelos de confiabilidad. El flujo mostrado aquí (exploración → candidatos → MLE → bondad de ajuste → AIC/BIC → validación visual) es transferible a prácticamente cualquier variable continua que encuentres en planta, laboratorio o modelo de negocio.

Si trabajas con datos de proceso metalúrgico, consumo de reactivos o confiabilidad de equipos, este flujo te servirá como plantilla base.