Comment puis-je utiliser SciPy pour trouver la distribution théorique la mieux adaptée à un ensemble de données empiriques et calculer les probabilités dépassant un seuil donné ?

Ajustement des distributions empiriques aux distributions théoriques avec Scipy

Aperçu du problème

Considérons un ensemble de données de valeurs entières échantillonnées à partir d'une inconnue distribution continue. Nous cherchons à déterminer la probabilité (valeur p) de rencontrer des valeurs supérieures à un seuil donné. Pour estimer avec précision ces probabilités, il est essentiel d'adapter notre distribution empirique à une distribution théorique appropriée. Cet article explique comment effectuer un tel ajustement à l'aide de Scipy en Python.

Ajustement de distribution

Pour évaluer la qualité de l'ajustement, nous pouvons utiliser une somme d'erreurs au carré (SSE) métrique pour comparer les histogrammes des données empiriques et les distributions ajustées. La distribution avec le SSE le plus bas est considérée comme la meilleure solution.

Mise en œuvre de Scipy

Le module statistique de Scipy fournit un large éventail de classes de distribution continue. Nous pouvons parcourir chaque distribution, estimer ses paramètres, calculer le SSE et stocker les résultats.

Exemple : El Niño Dataset

Illustrons le processus en utilisant les données de température de surface de la mer (SST) de l'ensemble de données El Niño.

1. Chargez les données et tracer son histogramme.
2. Effectuez un ajustement de distribution à l'aide de la métrique SSE.
3. Identifiez la distribution la mieux adaptée en fonction du SSE le plus bas.
4. Tracez la fonction de densité de probabilité (PDF) de la distribution la mieux adaptée ainsi que l'histogramme empirique.

Le code ci-dessous présente ceci implémentation :

import numpy as np
import pandas as pd
import scipy.stats as st
import matplotlib.pyplot as plt
from scipy.stats._continuous_distns import _distn_names
import warnings

# El Niño SST data
data = pd.Series(sm.datasets.elnino.load_pandas().data.set_index('YEAR').values.ravel())

# Function to fit distributions based on SSE
def best_fit_distribution(data):
    return sorted(
        [
            (getattr(st, distribution), distribution.fit(data), np.sum(np.power(data.hist(bins=50).values - distribution.pdf(data.index), 2.0))) 
            for distribution in _distn_names 
            if not distribution in ['levy_stable', 'studentized_range']
        ], 
    key=lambda x:x[2]
)

# Find best fit
best_dist = best_fit_distribution(data)[0]

# Plot distribution
fig, ax = plt.subplots(figsize=(12,8))
ax.plot(data.hist(bins=50, density=True, alpha=0.5, color='gray'))
param_names = best_dist[0].shapes + ', loc, scale' if best_dist[0].shapes else ['loc', 'scale']
param_str = ', '.join(['{}={:0.2f}'.format(k, v) for k, v in zip(param_names, best_dist[1])])
dist_str = '{}({})'.format(best_dist[0].name, param_str)

ax.plot(best_dist[0].pdf(data.index, **best_dist[1]), lw=2, label=dist_str)
ax.set_title('Fitted Distribution: ' + dist_str)
ax.set_xlabel('SST (°C)')
ax.set_ylabel('Frequency')
ax.legend()

La sortie montre la distribution la mieux adaptée comme la distribution de Weibull avec les paramètres :

scale=0.64, loc=15.59

Ce qui précède est le contenu détaillé de. pour plus d'informations, suivez d'autres articles connexes sur le site Web de PHP en chinois!