简介
核密度估计是一种用于估计随机变量概率密度函数的统计技术。在本实验中，我们将使用 Python 的 scikit-learn 库来演示一维核密度估计的原理。
虚拟机使用提示
虚拟机启动完成后，点击左上角切换到“笔记本”标签页，以访问 Jupyter Notebook 进行练习。
有时，你可能需要等待几秒钟让 Jupyter Notebook 完成加载。由于 Jupyter Notebook 的限制，操作验证无法自动化。
如果你在学习过程中遇到问题，随时向 Labby 提问。课程结束后提供反馈，我们会及时为你解决问题。
绘制直方图和核函数
我们首先绘制直方图和核函数，以展示两者之间的差异。我们将使用高斯核密度估计来展示两者的区别。我们还将比较 scikit-learn 中可用的其他核函数。
## 导入必要的库
import numpy as np
import matplotlib.pyplot as plt
from scipy.stats import norm
from sklearn.neighbors import KernelDensity
## 生成数据
np.random.seed(1)
N = 20
X = np.concatenate(
(np.random.normal(0, 1, int(0.3 * N)), np.random.normal(5, 1, int(0.7 * N)))
)[:, np.newaxis]
X_plot = np.linspace(-5, 10, 1000)[:, np.newaxis]
bins = np.linspace(-5, 10, 10)
## 创建图形和坐标轴
fig, ax = plt.subplots(2, 2, sharex=True, sharey=True)
fig.subplots_adjust(hspace=0.05, wspace=0.05)
## 绘制直方图 1
ax[0, 0].hist(X[:, 0], bins=bins, fc="#AAAAFF", density=True)
ax[0, 0].text(-3.5, 0.31, "直方图")
## 绘制直方图 2
ax[0, 1].hist(X[:, 0], bins=bins + 0.75, fc="#AAAAFF", density=True)
ax[0, 1].text(-3.5, 0.31, "直方图，箱线图偏移")
## 绘制顶帽核密度估计
kde = KernelDensity(kernel="tophat", bandwidth=0.75).fit(X)
log_dens = kde.score_samples(X_plot)
ax[1, 0].fill(X_plot[:, 0], np.exp(log_dens), fc="#AAAAFF")
ax[1, 0].text(-3.5, 0.31, "顶帽核密度")
## 绘制高斯核密度估计
kde = KernelDensity(kernel="gaussian", bandwidth=0.75).fit(X)
log_dens = kde.score_samples(X_plot)
ax[1, 1].fill(X_plot[:, 0], np.exp(log_dens), fc="#AAAAFF")
ax[1, 1].text(-3.5, 0.31, "高斯核密度")
## 绘制数据点
for axi in ax.ravel():
axi.plot(X[:, 0], np.full(X.shape[0], -0.01), "+k")
axi.set_xlim(-4, 9)
axi.set_ylim(-0.02, 0.34)
## 设置左列的 y 轴标签
for axi in ax[:, 0]:
axi.set_ylabel("归一化密度")
## 设置底行的 x 轴标签
for axi in ax[1, :]:
axi.set_xlabel("x")
## 显示图形
plt.show()
绘制可用核函数
我们将绘制所有可用的核函数以展示它们的形状。
## 生成数据
X_plot = np.linspace(-6, 6, 1000)[:, None]
X_src = np.zeros((1, 1))
## 创建图形和坐标轴
fig, ax = plt.subplots(2, 3, sharex=True, sharey=True)
fig.subplots_adjust(left=0.05, right=0.95, hspace=0.05, wspace=0.05)
## x 轴标签的格式化函数
def format_func(x, loc):
if x == 0:
return "0"
elif x == 1:
return "h"
elif x == -1:
return "-h"
else:
return "%ih" % x
## 绘制可用核函数
for i, kernel in enumerate(
["gaussian", "tophat", "epanechnikov", "exponential", "linear", "cosine"]
):
axi = ax.ravel()[i]
log_dens = KernelDensity(kernel=kernel).fit(X_src).score_samples(X_plot)
axi.fill(X_plot[:, 0], np.exp(log_dens), "-k", fc="#AAAAFF")
axi.text(-2.6, 0.95, kernel)
axi.xaxis.set_major_formatter(plt.FuncFormatter(format_func))
axi.xaxis.set_major_locator(plt.MultipleLocator(1))
axi.yaxis.set_major_locator(plt.NullLocator())
axi.set_ylim(0, 1.05)
axi.set_xlim(-2.9, 2.9)
## 设置第二行的标题
ax[0, 1].set_title("可用核函数")
## 显示图形
plt.show()
绘制一维密度示例
我们将绘制一个一维密度示例，其中包含 100 个样本。我们将比较三种不同的核密度估计：顶帽核（tophat）、高斯核（Gaussian）和叶甫根尼科夫核（epanechnikov）。
## 生成数据
N = 100
np.random.seed(1)
X = np.concatenate(
(np.random.normal(0, 1, int(0.3 * N)), np.random.normal(5, 1, int(0.7 * N)))
)[:, np.newaxis]
X_plot = np.linspace(-5, 10, 1000)[:, np.newaxis]
true_dens = 0.3 * norm(0, 1).pdf(X_plot[:, 0]) + 0.7 * norm(5, 1).pdf(X_plot[:, 0])
## 创建图形和坐标轴
fig, ax = plt.subplots()
## 绘制输入分布
ax.fill(X_plot[:, 0], true_dens, fc="black", alpha=0.2, label="输入分布")
## 设置颜色和核函数
colors = ["海军蓝", "矢车菊蓝", "深橙色"]
kernels = ["高斯核", "顶帽核", "叶甫根尼科夫核"]
lw = 2
## 绘制核密度估计
for color, kernel in zip(colors, kernels):
kde = KernelDensity(kernel=kernel, bandwidth=0.5).fit(X)
log_dens = kde.score_samples(X_plot)
ax.plot(
X_plot[:, 0],
np.exp(log_dens),
color=color,
lw=lw,
linestyle="-",
label="核函数 = '{0}'".format(kernel),
)
ax.text(6, 0.38, "N={0} 个点".format(N))
## 设置图例并绘制数据点
ax.legend(loc="左上角")
ax.plot(X[:, 0], -0.005 - 0.01 * np.random.random(X.shape[0]), "+k")
## 设置 x 和 y 轴范围
ax.set_xlim(-4, 9)
ax.set_ylim(-0.02, 0.4)
## 显示图形
plt.show()
总结
在本实验中，我们学习了如何使用 Python 的 scikit-learn 库来演示一维核密度估计的原理。我们绘制了直方图和核函数、可用的核函数以及一个一维密度示例。我们还比较了不同的核密度估计。