# How to decide whether a Kernel Density Estimate is good?

Consider this samples set, feel free to look at it by using your own favourite tools (e.g. R, Python, whatever other stats tools).

The problem is that I lack experience in deciding whether any of these PDF estimations are reasonable, or whether there is even a possibly better estimation.

Thus the question in the title:

• Q1: how to decide whether a KDE is good?

From this question raises another sub question:

• Q2: is it a problem that my tails aren't touching 0 on the Y axis?

Using the data above, here are various KDEs that I made using 3 different kernels and bandwidths as follows:

And the code that generated those plots (replace PATH/TO/SAMPLES/FILE by the correct path to the samples set):

import numpy as np
from sklearn.neighbors.kde import KernelDensity
import matplotlib.pyplot as plt
import math

fle = 'PATH/TO/SAMPLES/FILE'
with open(fle, 'r') as f:
X = [float(s) for s in f.read().splitlines()]
X.sort()

X_np = np.array(X)
X_np = X_np.reshape(-1,1)

# kernel density estimations
kernel='epanechnikov'
bw=0.001
kde = KernelDensity(kernel=kernel, bandwidth=bw).fit(X_np)
plt.plot(X_np, np.exp(kde.score_samples(X_np)), label='%s, bw=%s' % (kernel, bw))
plt.legend(loc=0)

kernel='gaussian'
bw=0.001
kde = KernelDensity(kernel=kernel, bandwidth=bw).fit(X_np)
plt.plot(X_np, np.exp(kde.score_samples(X_np)), label='%s, bw=%s' % (kernel, bw))
plt.legend(loc=0)

kernel='tophat'
bw=0.001
kde = KernelDensity(kernel=kernel, bandwidth=bw).fit(X_np)
plt.plot(X_np, np.exp(kde.score_samples(X_np)), label='%s, bw=%s' % (kernel, bw))
plt.legend(loc=0)

plt.show()

kernel='epanechnikov'
bw=0.01
kde = KernelDensity(kernel=kernel, bandwidth=bw).fit(X_np)
plt.plot(X_np, np.exp(kde.score_samples(X_np)), label='%s, bw=%s' % (kernel, bw))
plt.legend(loc=0)

kernel='gaussian'
bw=0.01
kde = KernelDensity(kernel=kernel, bandwidth=bw).fit(X_np)
plt.plot(X_np, np.exp(kde.score_samples(X_np)), label='%s, bw=%s' % (kernel, bw))
plt.legend(loc=0)

kernel='tophat'
bw=0.01
kde = KernelDensity(kernel=kernel, bandwidth=bw).fit(X_np)
plt.plot(X_np, np.exp(kde.score_samples(X_np)), label='%s, bw=%s' % (kernel, bw))
plt.legend(loc=0)

plt.show()

kernel='epanechnikov'
bw=0.05
kde = KernelDensity(kernel=kernel, bandwidth=bw).fit(X_np)
plt.plot(X_np, np.exp(kde.score_samples(X_np)), label='%s, bw=%s' % (kernel, bw))
plt.legend(loc=0)

kernel='gaussian'
bw=0.05
kde = KernelDensity(kernel=kernel, bandwidth=bw).fit(X_np)
plt.plot(X_np, np.exp(kde.score_samples(X_np)), label='%s, bw=%s' % (kernel, bw))
plt.legend(loc=0)

kernel='tophat'
bw=0.05
kde = KernelDensity(kernel=kernel, bandwidth=bw).fit(X_np)
plt.plot(X_np, np.exp(kde.score_samples(X_np)), label='%s, bw=%s' % (kernel, bw))
plt.legend(loc=0)

plt.show()

kernel='epanechnikov'
bw=0.09
kde = KernelDensity(kernel=kernel, bandwidth=bw).fit(X_np)
plt.plot(X_np, np.exp(kde.score_samples(X_np)), label='%s, bw=%s' % (kernel, bw))
plt.legend(loc=0)

kernel='gaussian'
bw=0.09
kde = KernelDensity(kernel=kernel, bandwidth=bw).fit(X_np)
plt.plot(X_np, np.exp(kde.score_samples(X_np)), label='%s, bw=%s' % (kernel, bw))
plt.legend(loc=0)

kernel='tophat'
bw=0.09
kde = KernelDensity(kernel=kernel, bandwidth=bw).fit(X_np)
plt.plot(X_np, np.exp(kde.score_samples(X_np)), label='%s, bw=%s' % (kernel, bw))
plt.legend(loc=0)

plt.show()


• While not very related (but I find it interesting) is your comment on the set of delta functions: do you mean the Dirac delta functions? If so, I assume the PDF that maximizes the likelihood would look like a comb (many tiny spikes) such that the length of each comb is $1/500$ (since I have $500$ data points) -- is this right? – caveman Jun 12 '16 at 10:23