python - 为什么相同数据的 scipy.stats.gaussian_kde() 比 seaborn.kde_plot() 慢？

Question

在 python 3.7 中，我有一个 shape=(2, 34900) 的 numpy 数组。该数组是一个坐标列表，其中索引 0 表示 X 轴，索引 1 表示 y 轴。

当我使用 seaborn.kde_plot() 可视化此数据的分布时，在第 7 代 i5 上运行时，我能够在大约 5-15 秒内得到结果。

但是当我尝试运行以下代码时:

#Find the kernel for 
k = scipy.stats.kde.gaussian_kde(data, bw_method=.3)
#Define the grid
xi, yi = np.mgrid[0:1:2000*1j, 0:1:2000*1j]
#apply the function
zi = k(np.vstack([xi.flatten(), yi.flatten()]))

为该数据找到高斯核并将其应用于我定义的网格，这需要更多时间。我无法运行完整阵列，但在大小为 140 的切片上运行时，大约需要 40 秒才能完成。

140 大小的切片确实产生了一个有趣的结果，我可以使用 plt.pcolormesh() 将其可视化。

我的问题是我在这里缺少什么。如果我理解正确发生了什么，我正在使用 scipy.stats.kde.gaussian_kde() 来创建由数据定义的函数的估计。然后我将该函数应用于 2D 空间并得到它的 Z 分量作为结果。然后我正在绘制 Z 分量。但是这个过程与seaborn.kde_plot() 这使得代码花费的时间更长。

Scipy 的实现只是遍历每个点:

for i in range(self.n):
    diff = self.dataset[:, i, newaxis] - points
    tdiff = dot(self.inv_cov, diff)
    energy = sum(diff*tdiff,axis=0) / 2.0
    result = result + exp(-energy)

Answer 1

Seaborn 通常有两种方法来计算双变量 kde。如果可用，它会使用 statsmodels，如果不可用，它会回退到 scipy。

scipy 代码类似于问题中显示的内容。它使用 scipy.stats.gaussian_kde . statsmodels 代码使用 statsmodels.nonparametric.api.KDEMultivariate .

但是，为了公平比较，我们需要为两种方法采用相同的网格大小。 seaborn 的标准 gridsize 是 100 点。

import numpy as np; np.random.seed(42)
import seaborn.distributions as sd

N = 34900
x = np.random.randn(N)
y = np.random.randn(N)
bw="scott"
gridsize=100
cut=3
clip = [(-np.inf, np.inf), (-np.inf, np.inf)]

f = lambda x,y : sd._statsmodels_bivariate_kde(x, y, bw, gridsize, cut, clip)
g = lambda x,y : sd._scipy_bivariate_kde(x, y, bw, gridsize, cut, clip)

如果我们对这两个函数计时，

# statsmodels
%timeit f(x,y)  # 1 loop, best of 3: 16.4 s per loop
# scipy
%timeit g(x,y)  # 1 loop, best of 3: 8.67 s per loop

因此，Scipy 的速度是 statsmodels(seaborn 默认值)的两倍。问题中的代码花费这么长时间的原因是使用了大小为 2000 的网格而不是大小为 100 的网格。

看到这些结果，人们实际上很想使用 scipy 而不是 statsmodels。不幸的是，它不允许选择使用哪一个。因此需要手动设置相应的标志。

import seaborn.distributions as sd
sd._has_statsmodels = False
# plot kdeplot with scipy.stats.kde.gaussian_kde
sns.kdeplot(x,y)