Python：用于多个分布计算的 Kullback Leibler 散度的快速/高效实现答案

【问题标题】：Python:Fast/efficient implementation of the Kullback Leibler divergence for multiple distributions calculationPython：用于多个分布计算的 Kullback Leibler 散度的快速/高效实现
【发布时间】：2014-09-14 19:27:55
【问题描述】：

我有以下问题：我有一个矩阵，比如说 20K 离散分布（直方图），我需要计算每对之间的 KL 散度 (KLD)。最简单的方法是使用两个 for 循环并通过标准 KLD 计算计算每两个分布之间的 KLD。这需要时间。很多时间。我想知道上面是否有某种基于矩阵/数组的计算。当然，我不是第一个遇到这个问题的人。不幸的是，我是 python 新手。任何帮助将不胜感激。谢谢！ A.

【问题讨论】：

到目前为止你做了什么？你有一些代码要看吗？
单个分布是什么样子的 - 即 20 个浮点值总和为 1.0？

标签： python performance algorithm statistics histogram

【解决方案1】：

这应该可行。

def kullback_leibler_divergence(X):

    """
    Finds the pairwise Kullback-Leibler divergence
    matrix between all rows in X.

    Parameters
    ----------
    X : array_like, shape (n_samples, n_features)
        Array of probability data. Each row must sum to 1.

    Returns
    -------
    D : ndarray, shape (n_samples, n_samples)
        The Kullback-Leibler divergence matrix. A pairwise matrix D such that D_{i, j}
        is the divergence between the ith and jth vectors of the given matrix X.

    Notes
    -----
    Based on code from Gordon J. Berman et al.
    (https://github.com/gordonberman/MotionMapper)

    References:
    -----------
    Berman, G. J., Choi, D. M., Bialek, W., & Shaevitz, J. W. (2014). 
    Mapping the stereotyped behaviour of freely moving fruit flies. 
    Journal of The Royal Society Interface, 11(99), 20140672.
    """

    X_log = np.log(X)
    X_log[np.isinf(X_log) | np.isnan(X_log)] = 0

    entropies = -np.sum(X * X_log, axis=1)

    D = np.matmul(-X, X_log.T)
    D = D - entropies
    D = D / np.log(2)
    D *= (1 - np.eye(D.shape[0]))

    return D

【讨论】：

【解决方案2】：

itertools 是一个强大的工具，用于处理不同实体之间的各种排列等。除其他外，它通过 product 函数执行双重“for”循环：

from itertools import product
KL = [kl(hist_mat[:,i],hist_mat[:,j]) for i,j in product( range(0,K), range(0,K) )]

K 是以成对方式比较的直方图的数量，
hist_mat 包含每个列中的直方图。
kl = kullback leibler 散度的选择实现

唯一的缺点是，对于大型数据集，它仍然运行缓慢。

干杯， A.

【讨论】：

【解决方案3】：

我发现 Computation of Kullback-Leibler (KL) distance between text-documents using numpy 指出 SciPy 有一个 KLD 的实现

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

可以在http://docs.scipy.org/doc/scipy-dev/reference/generated/scipy.stats.entropy.html 找到文档；这应该会使计算本身更快。

或者，an implementation in numpy:

import numpy as np

def kl(p, q):
"""Kullback-Leibler divergence D(P || Q) for discrete distributions

Parameters
----------
p, q : array-like, dtype=float, shape=n
Discrete probability distributions.
"""
p = np.asarray(p, dtype=np.float)
q = np.asarray(q, dtype=np.float)

return np.sum(np.where(p != 0, p * np.log(p / q), 0))

请注意，将数据与 numpy 数组相互转换相对较慢；我不知道是保留一个 20k 1D numpy 数组的列表（可能非常占用内存），还是保留一个 2D numpy 数组并在切片上操作是否会更好。

【讨论】：