Estimation de la densité de noyau simple en 1D

Introduction

L'estimation de la densité de noyau est une technique statistique utilisée pour estimer la fonction de densité de probabilité d'une variable aléatoire. Dans ce laboratoire, nous utiliserons la bibliothèque Python scikit-learn pour démontrer les principes de l'estimation de la densité de noyau en une dimension.

Conseils sur la machine virtuelle

Une fois le démarrage de la machine virtuelle terminé, cliquez dans le coin supérieur gauche pour basculer vers l'onglet Carnet de notes pour accéder au carnet Jupyter pour pratiquer.

Parfois, vous devrez peut-être attendre quelques secondes pour que le carnet Jupyter ait fini de charger. La validation des opérations ne peut pas être automatisée en raison des limitations du carnet Jupyter.

Si vous rencontrez des problèmes pendant l'apprentissage, n'hésitez pas à demander à Labby. Donnez des commentaires après la session, et nous réglerons rapidement le problème pour vous.

Skills Graph

Tracer des histogrammes et des noyaux

Nous traçons d'abord des histogrammes et des noyaux pour démontrer la différence entre les deux. Nous utiliserons une estimation de la densité de noyau gaussienne pour montrer la différence entre les deux. Nous comparerons également d'autres noyaux disponibles dans scikit-learn.

## Import necessary libraries
import numpy as np
import matplotlib.pyplot as plt
from scipy.stats import norm
from sklearn.neighbors import KernelDensity

## Generate data
np.random.seed(1)
N = 20
X = np.concatenate(
    (np.random.normal(0, 1, int(0.3 * N)), np.random.normal(5, 1, int(0.7 * N)))
)[:, np.newaxis]
X_plot = np.linspace(-5, 10, 1000)[:, np.newaxis]
bins = np.linspace(-5, 10, 10)

## Create figure and axes
fig, ax = plt.subplots(2, 2, sharex=True, sharey=True)
fig.subplots_adjust(hspace=0.05, wspace=0.05)

## Plot histogram 1
ax[0, 0].hist(X[:, 0], bins=bins, fc="#AAAAFF", density=True)
ax[0, 0].text(-3.5, 0.31, "Histogramme")

## Plot histogram 2
ax[0, 1].hist(X[:, 0], bins=bins + 0.75, fc="#AAAAFF", density=True)
ax[0, 1].text(-3.5, 0.31, "Histogramme, bins décalés")

## Plot tophat KDE
kde = KernelDensity(kernel="tophat", bandwidth=0.75).fit(X)
log_dens = kde.score_samples(X_plot)
ax[1, 0].fill(X_plot[:, 0], np.exp(log_dens), fc="#AAAAFF")
ax[1, 0].text(-3.5, 0.31, "Densité de noyau tophat")

## Plot Gaussian KDE
kde = KernelDensity(kernel="gaussien", bandwidth=0.75).fit(X)
log_dens = kde.score_samples(X_plot)
ax[1, 1].fill(X_plot[:, 0], np.exp(log_dens), fc="#AAAAFF")
ax[1, 1].text(-3.5, 0.31, "Densité de noyau gaussienne")

## Plot data points
for axi in ax.ravel():
    axi.plot(X[:, 0], np.full(X.shape[0], -0.01), "+k")
    axi.set_xlim(-4, 9)
    axi.set_ylim(-0.02, 0.34)

## Set y-axis label for left column
for axi in ax[:, 0]:
    axi.set_ylabel("Densité normalisée")

## Set x-axis label for bottom row
for axi in ax[1, :]:
    axi.set_xlabel("x")

## Show plot
plt.show()

Tracer les noyaux disponibles

Nous allons tracer tous les noyaux disponibles pour montrer leurs formes.

## Generate data
X_plot = np.linspace(-6, 6, 1000)[:, None]
X_src = np.zeros((1, 1))

## Create figure and axes
fig, ax = plt.subplots(2, 3, sharex=True, sharey=True)
fig.subplots_adjust(left=0.05, right=0.95, hspace=0.05, wspace=0.05)

## Format function for x-axis labels
def format_func(x, loc):
    if x == 0:
        return "0"
    elif x == 1:
        return "h"
    elif x == -1:
        return "-h"
    else:
        return "%ih" % x

## Plot available kernels
for i, kernel in enumerate(
    ["gaussien", "tophat", "épanechnikov", "exponentielle", "linéaire", "cosinus"]
):
    axi = ax.ravel()[i]
    log_dens = KernelDensity(kernel=kernel).fit(X_src).score_samples(X_plot)
    axi.fill(X_plot[:, 0], np.exp(log_dens), "-k", fc="#AAAAFF")
    axi.text(-2.6, 0.95, kernel)

    axi.xaxis.set_major_formatter(plt.FuncFormatter(format_func))
    axi.xaxis.set_major_locator(plt.MultipleLocator(1))
    axi.yaxis.set_major_locator(plt.NullLocator())

    axi.set_ylim(0, 1.05)
    axi.set_xlim(-2.9, 2.9)

## Set title for second row
ax[0, 1].set_title("Noyaux disponibles")

## Show plot
plt.show()

Tracer un exemple de densité en 1D

Nous allons tracer un exemple de densité en 1D avec 100 échantillons en une dimension. Nous comparerons trois estimations de densité de noyau différentes : tophat, gaussienne et épanechnikov.

## Generate data
N = 100
np.random.seed(1)
X = np.concatenate(
    (np.random.normal(0, 1, int(0.3 * N)), np.random.normal(5, 1, int(0.7 * N)))
)[:, np.newaxis]

X_plot = np.linspace(-5, 10, 1000)[:, np.newaxis]

true_dens = 0.3 * norm(0, 1).pdf(X_plot[:, 0]) + 0.7 * norm(5, 1).pdf(X_plot[:, 0])

## Create figure and axes
fig, ax = plt.subplots()

## Plot input distribution
ax.fill(X_plot[:, 0], true_dens, fc="black", alpha=0.2, label="distribution d'entrée")

## Set colors and kernels
colors = ["navy", "cornflowerblue", "darkorange"]
kernels = ["gaussien", "tophat", "épanechnikov"]
lw = 2

## Plot kernel density estimations
for color, kernel in zip(colors, kernels):
    kde = KernelDensity(kernel=kernel, bandwidth=0.5).fit(X)
    log_dens = kde.score_samples(X_plot)
    ax.plot(
        X_plot[:, 0],
        np.exp(log_dens),
        color=color,
        lw=lw,
        linestyle="-",
        label="kernel = '{0}'".format(kernel),
    )

ax.text(6, 0.38, "N={0} points".format(N))

## Set legend and plot data points
ax.legend(loc="upper left")
ax.plot(X[:, 0], -0.005 - 0.01 * np.random.random(X.shape[0]), "+k")

## Set x and y limits
ax.set_xlim(-4, 9)
ax.set_ylim(-0.02, 0.4)

## Show plot
plt.show()

Sommaire

Dans ce laboratoire, nous avons appris à utiliser la bibliothèque Python scikit-learn pour démontrer les principes de l'estimation de la densité de noyau en une dimension. Nous avons tracé des histogrammes et des noyaux, les noyaux disponibles et un exemple de densité en 1D. Nous avons également comparé différentes estimations de densité de noyau.