Python中计算Kullback-Leibler散度的高效方法

Question

Python中计算Kullback-Leibler散度的高效方法

5

我需要计算成千上万个离散概率向量之间的Kullback-Leibler Divergence（KLD）。目前我正在使用以下代码，但对于我的目的来说太慢了。我想知道是否有更快的方法来计算KL散度？

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, :])

- Amir

我该如何在Python中使用Scipy来获取具有最小KL散度的概率分布生成器？ - yishairasowsky

1个回答

网页内容由stack overflow 提供, 点击上面的

可以查看英文原文，
原文链接

- Divakar · Accepted Answer

Scipy的 stats.entropy 在其默认情况下接受作为输入的一维数组并给出一个标量，这在所列出的问题中已经完成。该函数在内部还允许broadcasting，我们可以在此处滥用它以获得矢量化解决方案。

从文档中-

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)，使用broadcasting为我们提供一个形状为(M,M) 的输出数组。

因此，解决我们的情况的矢量化和更有效的方法是-

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