stats - mutual information python

Forma óptima de calcular información mutua por pares usando numpy

Question

Para una matriz mxn , ¿cuál es la forma óptima (más rápida) de calcular la información mutua para todos los pares de columnas ( nxn )?

Por información mutua , quiero decir:

I (X, Y) = H (X) + H (Y) - H (X, Y)

donde H (X) se refiere a la entropía de Shannon de X.

Actualmente estoy usando np.histogram2d y np.histogram para calcular los np.histogram de articulaciones (X, Y) e individuales (X o Y) . Para una matriz A dada (por ejemplo, una matriz de flotadores de 250000 X 1000), estoy haciendo un ciclo anidado for

 n = A.shape[1]
for ix = arange(n)
for jx = arange(ix+1,n):
matMI[ix,jx]= calc_MI(A[:,ix],A[:,jx])

¿Seguramente debe haber formas mejores / más rápidas de hacer esto?

Como comentario adicional, también he buscado funciones de mapeo en columnas (operaciones en columnas o filas) en arreglos, pero aún no he encontrado una buena respuesta general.

Aquí está mi implementación completa, siguiendo las convenciones en la página Wiki :

import numpy as np
def calc_MI(X,Y,bins):
c_XY = np.histogram2d(X,Y,bins)[0]
c_X = np.histogram(X,bins)[0]
c_Y = np.histogram(Y,bins)[0]
H_X = shan_entropy(c_X)
H_Y = shan_entropy(c_Y)
H_XY = shan_entropy(c_XY)
MI = H_X + H_Y - H_XY
return MI
def shan_entropy(c):
c_normalized = c / float(np.sum(c))
c_normalized = c_normalized[np.nonzero(c_normalized)]
H = -sum(c_normalized* np.log2(c_normalized))
return H
A = np.array([[ 2.0, 140.0, 128.23, -150.5, -5.4 ],
[ 2.4, 153.11, 130.34, -130.1, -9.5 ],
[ 1.2, 156.9, 120.11, -110.45,-1.12 ]])
bins = 5 # ?
n = A.shape[1]
matMI = np.zeros((n, n))
for ix in np.arange(n):
for jx in np.arange(ix+1,n):
matMI[ix,jx] = calc_MI(A[:,ix], A[:,jx], bins)

Aunque mi versión de trabajo con anidados for bucles lo hace a una velocidad razonable, me gustaría saber si hay una forma más óptima de aplicar calc_MI en todas las columnas de A (para calcular su información mutua por pares).

También me gustaría saber:

1. ¿Hay formas eficientes de asignar funciones para operar en columnas (o filas) de np.arrays (tal vez como np.vectorize , que se parece más a un decorador)?

2. ¿Hay otras implementaciones óptimas para este cálculo específico (información mutua)?

Answer 1

No puedo sugerir un cálculo más rápido para el bucle externo sobre los vectores n * (n-1) / 2, pero su implementación de calc_MI(x, y, bins) puede simplificarse si puede usar scipy versión 0.13 o scikit-learn

En scipy 0.13, el argumento lambda_ se agregó a scipy.stats.chi2_contingency Este argumento controla la estadística que calcula la función. Si usa lambda_="log-likelihood" (o lambda_=0 ), se devuelve la relación log-likelihood. Esto también se llama a menudo la estadística G o G ² . Aparte de un factor de 2 * n (donde n es el número total de muestras en la tabla de contingencia), esta es la información mutua. Entonces podrías implementar calc_MI como:

from scipy.stats import chi2_contingency
def calc_MI(x, y, bins):
c_xy = np.histogram2d(x, y, bins)[0]
g, p, dof, expected = chi2_contingency(c_xy, lambda_="log-likelihood")
mi = 0.5 * g / c_xy.sum()
return mi

La única diferencia entre esto y su implementación es que esta implementación utiliza el logaritmo natural en lugar del logaritmo base-2 (por lo que expresa la información en "nats" en lugar de "bits"). Si realmente prefiere bits, simplemente divida mi por log (2).

Si tiene (o puede instalar) sklearn (es decir, scikit-learn), puede usar sklearn.metrics.mutual_info_score e implementar calc_MI como:

from sklearn.metrics import mutual_info_score
def calc_MI(x, y, bins):
c_xy = np.histogram2d(x, y, bins)[0]
mi = mutual_info_score(None, None, contingency=c_xy)
return mi