在Python中计算Kullback-Leibler分歧的有效方法

Question

我必须计算数千个离散概率向量之间的Kullback-Leibler发散(KLD).目前我正在使用以下代码,但它对我的目的来说太慢了.我想知道是否有更快的方法来计算KL Divergence？

import numpy as np
import scipy.stats as sc

    #n is the number of data points
    kld = np.zeros(n, n)
        for i in range(0, n):
            for j in range(0, n):
                if(i != j):
                    kld[i, j] = sc.entropy(distributions[i, :], distributions[j, :])

Answer 1

stats.entropy默认意义上的Scipy将输入作为一维数组给我们一个标量,这是在列出的问题中完成的.在内部,此功能也允许broadcasting,我们可以在这里滥用矢量化解决方案.

来自docs-

scipy.stats.entropy(pk,qk = None,base = None)

如果仅给出概率pk,则将熵计算为S = -sum(pk*log(pk),axis = 0).

如果qk不是None,则计算Kullback-Leibler散度S = sum(pk*log(pk/qk),axis = 0).

在我们的例子中,我们针对所有行对每一行进行这些熵计算,执行求和减少以在每次迭代时使用这两个嵌套循环具有标量.因此,输出数组将具有形状(M,M),其中M是输入数组中的行数.

现在,这里的结果是stats.entropy()总和axis=0,所以我们将提供两个版本,两个版本distributions都将axis=0沿着它减少行和尺寸以及其他两个轴交错 - (M,1)并(1,M)给我们一个(M,M)整形输出数组使用broadcasting.

因此,解决我们案件的矢量化和更有效的方法是 -

from scipy import stats
kld = stats.entropy(distributions.T[:,:,None], distributions.T[:,None,:])

运行时测试和验证 -

In [15]: def entropy_loopy(distrib):
    ...:     n = distrib.shape[0] #n is the number of data points
    ...:     kld = np.zeros((n, n))
    ...:     for i in range(0, n):
    ...:         for j in range(0, n):
    ...:             if(i != j):
    ...:                 kld[i, j] = stats.entropy(distrib[i, :], distrib[j, :])
    ...:     return kld
    ...: 

In [16]: distrib = np.random.randint(0,9,(100,100)) # Setup input

In [17]: out = stats.entropy(distrib.T[:,:,None], distrib.T[:,None,:])

In [18]: np.allclose(entropy_loopy(distrib),out) # Verify
Out[18]: True

In [19]: %timeit entropy_loopy(distrib)
1 loops, best of 3: 800 ms per loop

In [20]: %timeit stats.entropy(distrib.T[:,:,None], distrib.T[:,None,:])
10 loops, best of 3: 104 ms per loop