68

mxn行列の場合、列のすべてのペア ( nxn )の相互情報を計算する最適な (最速の) 方法は何ですか?

相互情報とは、次のことを意味します。

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

ここで、H(X)はXのシャノン エントロピーを表します。

現在、私は and を使用np.histogram2dしnp.histogramて、関節(X、Y)と個々の(X または Y)の数を計算しています。特定の行列A(たとえば、250000 X 1000 の float 行列) に対して、ネストされたforループを実行しています。

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

確かにこれを行うためのより良い/より速い方法があるに違いありませんか?

余談ですが、配列の列に対するマッピング関数 (列方向または行方向の操作) も探しましたが、まだ適切な一般的な答えは見つかりませんでした。

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)

ネストされたループを使用した私の作業バージョンは妥当な速度で実行されますが、 (ペアごとの相互情報を計算するために) のすべての列にfor適用するより最適な方法があるかどうかを知りたいですか?calc_MIA

私も知りたいです：

  1. 関数をマップして列 (または行) を操作する効率的な方法があるかどうか (np.arraysデコレータnp.vectorizeのように見える のように)?

  2. この特定の計算に最適な実装が他にあるかどうか (相互情報)?

4

1 に答える 1

67

n*(n-1)/2 ベクトルに対する外側のループの高速な計算を提案することはできませんが、 scipy バージョン 0.13 またはscikit-learncalc_MI(x, y, bins)を使用できる場合、実装を簡素化できます。

scipy 0.13 では、lambda_引数が追加されましたscipy.stats.chi2_contingency 。この引数は、関数によって計算される統計を制御します。lambda_="log-likelihood"(または)を使用するlambda_=0と、対数尤度比が返されます。これは、G または G 2統計とも呼ばれます。係数 2*n (n は分割表のサンプルの総数) を除いて、これは相互情報量です。したがって、次のように実装できますcalc_MI 。

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

これとあなたの実装の唯一の違いは、この実装では底 2 の対数ではなく自然対数を使用することです (したがって、「ビット」ではなく「nats」で情報を表現しています)。本当にビットを好む場合は、milog(2) で割ります。

持っている (またはインストールできる) sklearn(つまり、scikit-learn) 場合は、 を使用 して、次のようsklearn.metrics.mutual_info_scoreに実装できますcalc_MI。

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
于 2013-12-10T21:20:10.690 に答える