Calcul de la quantité d'informations mutuelles (valeur continue) avec numpy

motivation

Je veux calculer le montant d'informations mutuelles $ I (X; Y) $ des variables de probabilité continues $ X $ et $ Y $ en Python. $ I(X;Y) = \int_Y \int_X p(x, y) \log \frac{p(x,y)}{p(x)p(y)} dx dy $

code

import numpy

def mutual_information(X, Y, bins=10):
    #Distribution de probabilité simultanée p(x,y)Calculs de
    p_xy, xedges, yedges = np.histogram2d(X, Y, bins=bins, density=True)
    
    # p(x)p(y)Calculs de
    p_x, _ = np.histogram(X, bins=xedges, density=True)
    p_y, _ = np.histogram(Y, bins=yedges, density=True)
    p_x_y = p_x[:, np.newaxis] * p_y
    
    #dx et dy
    dx = xedges[1] - xedges[0]
    dy = yedges[1] - yedges[0]
    
    #Éléments d'intégration
    elem = p_xy * np.ma.log(p_xy / p_x_y)
    #Montant d'informations mutuelles et p(x, y), p(x)p(y)Production
    return np.sum(elem * dx * dy), p_xy, p_x_y

point

Si vous souhaitez calculer la quantité d'informations mutuelles pour le moment, vous pouvez utiliser la fonction ci-dessus. Soit dit en passant, je laisserai quelques points importants pour la mise en œuvre.

«Densité» de «np.histogram2d»

J'étais un peu impatient car np.sum (p_xy) ne devenait pas 1 quand je pensais vaguement que la probabilité serait renvoyée si densité = True était défini. Le point à noter est que «p_xy» est la ** densité de probabilité **, pas la probabilité.

Puisque $ X $ et $ Y $ sont des variables continues, l'approximation dans l'histogramme est la densité de probabilité. Par conséquent, si vous les additionnez en tenant compte de la largeur du bac, ce sera 1.

«np.histogram» et «np.histogram2d» renvoient la densité de probabilité et les intervalles (arêtes dans le code). Il est nécessaire de calculer «dx» et «dy» à partir de ce bac.

import numpy as np

N = 1000
X = np.random.normal(loc=0, scale=1, size=N)

p_x, edges = np.histogram(X, bins=10, density=True)

#Si vous prenez la somme des densités de probabilité sans réfléchir, ce ne sera évidemment pas 1.
print(np.sum(p_x))  #Exemple de sortie: 1.580769264599771

#Si vous prenez la somme en tenant compte de la largeur du bac, elle devient 1.
dx = edges[1] - edges[0]
print(np.sum(p_x * dx))  #Exemple de sortie: 1.0000000000000002

Calcul de p_x_y

P_x_y dans le code essaie de calculer $ p (x) p (y) $. En fait, j'ai d'abord calculé avec le code suivant et cela n'a pas fonctionné.

p_x_y = p_x * p_y

Correctement

p_x_y = p_x[:, np.newaxis] * p_y

est. Dans le premier, p_x_y est le tableau principal, et dans le second, p_x_y est le tableau secondaire.

Exemple d'exécution

Exemple d'exécution 1 (deux ondes sin)

Comme ils ne sont pas indépendants, il y a une différence entre $ p (x, y) $ et $ p (x) p (y) $, et la quantité d'informations mutuelles augmente.

import matplotlib.pyplot as plt

#sin wave et cos wave
t = np.linspace(-5, 5, num=1000)
X = np.sin(2 * np.pi * t)
Y = np.cos(3 * np.pi * t)

#Calcul du montant d'informations mutuelles
mi, p_xy, p_x_y = mutual_information(X, Y, bins=30)

#Sortie de résultat
plt.figure(dpi=100)
ax1 = plt.subplot(121)
ax2 = plt.subplot(122)
ax1.set_title(r'$P_{XY}(x, y)$')
ax1.imshow(p_xy)
ax2.set_title(r'$P_{X}(x) P_{Y}(y)$')
ax2.imshow(p_x_y)
plt.suptitle('MI = {}'.format(mi))
plt.show()

image.png

Exemple d'exécution 2 (distribution normale indépendante)

Lorsque les deux variables sont indépendantes, $ p (x, y) $ et $ p (x) p (y) $ correspondent, et la quantité d'informations mutuelles devient petite.

import matplotlib.pyplot as plt
#Deux distributions normales indépendantes
N = 10000
X = np.random.normal(size=N)
Y = np.random.normal(size=N)

#Calcul du montant d'informations mutuelles
mi, p_xy, p_x_y = mutual_information(X, Y, bins=30)
Exemple d'exécution
#Sortie de résultat
plt.figure(dpi=100)
ax1 = plt.subplot(121)
ax2 = plt.subplot(122)
ax1.set_title(r'$P_{XY}(x, y)$')
ax1.imshow(p_xy)
ax2.set_title(r'$P_{X}(x) P_{Y}(y)$')
ax2.imshow(p_x_y)
plt.suptitle('MI = {}'.format(mi))
plt.show()

image.png

Recommended Posts

Calcul de la quantité d'informations mutuelles (valeur continue) avec numpy
1. Statistiques apprises avec Python 1-2. Calcul de diverses statistiques (Numpy)
Calcul séquentiel de la valeur moyenne avec l'algorithme en ligne
Calcul sans erreur avec le big.Float de Golang
L'histoire du calcul numérique des équations différentielles avec TensorFlow 2.0
1. Statistiques apprises avec Python 1-3. Calcul de diverses statistiques (statistiques)
Calcul en temps réel de la valeur moyenne avec corroutine
Calcul de la probabilité de concentration d'étoile FGO / valeur attendue
Convertissez les données avec la forme (nombre de données, 1) en (nombre de données,) avec numpy.
Partez numpy? !! Différencier partiellement la matrice avec Sympy
Vitesse de calcul de l'indexation pour un tableau quadratique numpy
Effectue le calcul à grande vitesse de descripteurs spécifiques uniquement avec mordred
Ajoutez des informations au bas de la figure avec Matplotlib
Prenez la valeur du thermo-hygromètre SwitchBot avec Raspberry Pi
Changer les valeurs du thermo-hygromètre Bot avec Raspberry Pi