python - scipy gaussian_kde 根据使用的方法(权重与无权重)产生不同的结果

Question

我有一系列要应用 KDE 的坐标，并且一直在使用 scipy.stats.gaussian_kde这样做。这里的问题是这个函数需要一组离散的坐标，然后它会执行密度估计。
当我希望记录我的数据时，这会导致问题(对于坐标特别空闲的集合，并且使用未触及的数据提供的信息很少)。可以想象，如果您必须处理离散数量的点，如果 2 个点出现 18 次，其他 24 次，取 18 和 24 的对数将使它们相同，因为它们必须按顺序四舍五入到最接近的整数保持离散。
为了解决这个问题，我一直在使用 weights scipy.stats.gaussian_kde 中的参数功能。不是有一个数组，其中每个点出现的次数等于其密度，而是每个点出现一次，而是由其密度加权。所以现在，使用之前的例子，密度为 18 和 24 的 2 个点将不相同，因为这些密度可以是连续的。
这有效并产生了似乎是一个很好的估计，但是使用这两种不同的方法，它们都会产生具有细微差别的图形。如果我只使用一种方法，我会保持无知，但现在我已经使用了两种方法，我无法确定估计是否合理。
这两种方法产生不同的结果是否有原因？
请参阅下面一些重现该问题的示例代码:

from scipy.stats import gaussian_kde
import numpy as np
import matplotlib.pyplot as plt

np.random.seed(1)

discrete_points = np.random.randint(0,10,size=(2,400))

continuous_points = [[],[]]
continuous_weights = []

recorded_points = []
for i in range(discrete_points.shape[1]):
    p = discrete_points[:,i]
    if tuple(p) in recorded_points:
        continuous_weights[recorded_points.index(tuple(p))] += 1
    else:
        continuous_points[0].append(p[0])
        continuous_points[1].append(p[1])
        continuous_weights.append(1)
        recorded_points.append(tuple(p))

resolution = 1

kde = gaussian_kde(discrete_points)
x, y = discrete_points
# https://www.oreilly.com/library/view/python-data-science/9781491912126/ch04.html
x_step = int((max(x)-min(x))/resolution)
y_step = int((max(y)-min(y))/resolution)
xgrid = np.linspace(min(x), max(x), x_step+1)
ygrid = np.linspace(min(y), max(y), y_step+1)
Xgrid, Ygrid = np.meshgrid(xgrid, ygrid)
Z = kde.evaluate(np.vstack([Xgrid.ravel(), Ygrid.ravel()]))
Zgrid = Z.reshape(Xgrid.shape)

ext = [min(x), max(x), min(y), max(y)]
earth = plt.cm.gist_earth_r

plt.imshow(Zgrid,
    origin='lower', aspect='auto',
    extent=ext,
    alpha=0.8,
    cmap=earth)

plt.title("Discrete method (no weights)")
plt.savefig("noweights.png")

kde = gaussian_kde(continuous_points, weights=continuous_weights)
x, y = continuous_points
# https://www.oreilly.com/library/view/python-data-science/9781491912126/ch04.html
x_step = int((max(x)-min(x))/resolution)
y_step = int((max(y)-min(y))/resolution)
xgrid = np.linspace(min(x), max(x), x_step+1)
ygrid = np.linspace(min(y), max(y), y_step+1)
Xgrid, Ygrid = np.meshgrid(xgrid, ygrid)
Z = kde.evaluate(np.vstack([Xgrid.ravel(), Ygrid.ravel()]))
Zgrid = Z.reshape(Xgrid.shape)

ext = [min(x), max(x), min(y), max(y)]
earth = plt.cm.gist_earth_r

plt.imshow(Zgrid,
    origin='lower', aspect='auto',
    extent=ext,
    alpha=0.8,
    cmap=earth)

plt.title("Continuous method (weights)")
plt.savefig("weights.png")

产生以下图:

和

Answer 1

kde 的一个进口方面是使用的带宽。 Scipy's gaussian_kde 用途 "Scott's factor"作为对带宽的猜测。
特别是gaussian_kde用途 n**(-1./(d+4))哪里d是维度(在本例中为 2)，而 n

非加权版本的数据点数

加权版本的“有效数据点数”；计算为 neff = sum(weights)^2 / sum(weights^2)

在帖子的例子中 n = 400和 neff = sum(continuous_weights)**2 / sum([w**2 for w in continuous_weights]) = 84.0336 .
为了获得相同的结果，在两种情况下都应该使用相同的带宽。它可以明确设置为 gaussian_kde(..., bw_method=bandwidth .

bandwidth = discrete_points.shape[1]**(-1./(2+4))

# kde without weights
kde = gaussian_kde(discrete_points, bw_method=bandwidth)

# kde for the weighted points
kde = gaussian_kde(continuous_points, weights=continuous_weights, bw_method=bandwidth)

如果您计划创建多个图，您可能希望对所有图使用相同的带宽，而不受点数或权重的影响。您可能想要试验分辨率和带宽。更高的带宽可以在更大的距离上平滑所有内容，更小的带宽更忠实于给定的数据。