核密度估计(Kernel Density Estimation)
密度评估器是一种利用D维数据集生成D维概率分布估计的算法,GMM使用不同的高斯分布的加权汇总来表示概率分布估计。
核密度估计算法将高斯混合理念扩展到了逻辑极限,它通过对每一个点生成高斯分布的混合成分,获得本实质上是无参数的密度评估器。
核密度估计的自由参数是核类型和核带宽,前者指定每个点核密度分布的形状,后者指定每个点核的大小。
一维数据的密度估计——直方图,是一个简单的密度评估器,直方图将数据分成若干区间,统计落入每个区间内的点的数量,然后用直观的方式将结果可视化。
核密度估计示例
核密度估计的自由参数是核类型(kernel)参数,他可以指定每个点核密度分布的形状。
核带宽(kernel bandwidth)参数控制每个点的核的大小
核密度估计算法在sklearn.neighbors.KernelDensity评估器中实现,借助六个核中的任意一个核、二三十个距离量度就可以处理具有多个维度的KDE。
由于KDE计算量非常大,因此Scikit-Learn评估器底层使用了一种基于树的算法,可以利用atol(绝对容错)和rtol(相对容错)参数来平衡计算时间与准确性,可以用Scikit-Learn的标准交叉检验工具来确定自由参数核带宽。
以下代码生成一组高斯分布的随机数据,分别绘制其直方图与KDE图,查看效果;
import numpy as np
from matplotlib import pyplot as plt
from sklearn.model_selection import GridSearchCV, LeaveOneOut
from sklearn.neighbors import KernelDensity
x_train = np.hstack((np.random.normal(2, 1.66, 200), np.random.normal(8, 2.33, 200))) # 两个高斯分布组合到一起
bandwidth = 1.0
model = KernelDensity(bandwidth=bandwidth, kernel='gaussian')
model.fit(x_train[:, np.newaxis])
x_range = np.linspace(x_train.min() - 1, x_train.max() + 1, 500)
x_log_prob = model.score_samples(x_range[:, np.newaxis]) # 这个东西返回概率的对数
x_prob = np.exp(x_log_prob)
plt.figure(figsize=(10, 10))
r = plt.hist(
x=x_train,
bins=50,
density=True,
histtype='stepfilled',
color='red',
alpha=0.5,
label='histogram',
)
plt.fill_between(
x=x_range,
y1=x_prob,
y2=0,
color='green',
alpha=0.5,
label='KDE',
)
plt.plot(x_range, x_prob, color='gray')
plt.vlines(x=2, ymin=0, ymax=r[0].max() + 0.01, color='k', linestyle='--', alpha=0.7) # 分布中心
plt.vlines(x=8, ymin=0, ymax=r[0].max() + 0.01, color='k', linestyle='--', alpha=0.7) # 分布中心
plt.ylim(0, r[0].max() + 0.011)
plt.legend(loc='upper right')
plt.title('histogram and kde')
plt.show()
结果:
使用交叉检验确定带宽
kernel参数指定需要用到的核函数,来适应不同的需求;
在使用核密度估计时,如果带宽设置过小,会出现过拟合的现象,如果带宽设置过大,会出现欠拟合的现象,因此需要确定好最佳的带宽;
sklearn中的KernelDensity类支持使用GridSearchCV来寻找最佳参数,以下是一个示例:
def get_suit_bandwidth(x_train):
grid = GridSearchCV(
estimator=KernelDensity(kernel='gaussian'),
param_grid={'bandwidth': 10 ** np.linspace(-1, 1, 100)},
cv=LeaveOneOut(),
)
grid.fit(x_train[:, np.newaxis])
print(f'最佳带宽:{grid.best_params_["bandwidth"]}')
return grid.best_params_["bandwidth"]
bandwidth = get_suit_bandwidth(x_train)
bandwidth 大约0.7左右,再次运行上述代码:
可以看出宽度变窄后你和效果更完美一些。
用二维数据也试验了一下,也不错: