https://zhuanlan.zhihu.com/p/336716858

你可能聽說過核密度估計(KDE：kernel density estimation)或非引數迴歸（non-parametric regression）。你甚至可能在不知不覺的情況下使用它。比如在Python中使用seaborn或plotly時，distplot就是這樣，在預設情況下都會使用核密度估計器。但是這些大概是什麼意思呢?

也許你處理了一個迴歸問題，卻發現線性迴歸不能很好地工作，因為特性和標籤之間的依賴似乎是非線性的。在這裡，核迴歸（kernel regression）可能是一種解決方案。

在這篇文章中，我們通過示例，並試圖對核心估計背後的理論有一個直觀的理解。此外，我們還看到了這些概念在Python中的實現。

核迴歸

圖1:全球谷歌搜尋“chocolate”;x軸:時間，y軸:搜尋百分比

讓我們從一個例子開始。假設你是一個數據科學家，在一家糖果工廠的巧克力部門工作。

你可能想要預測巧克力的需求基於它的歷史需求，作為第一步，想要分析趨勢。2004-2020年的巧克力需求可能類似於圖1中的資料。

顯然，這是有季節性的，冬天的需求會增加，但是由於你對趨勢感興趣，你決定擺脫這些波動。為此，你可以計算視窗為b個月的移動平均線，也就是說，對於每一個時刻t，你計算從t-b到t+b的時間段內需求的平均值。

更正式地說，如果我們有一段時間內觀察到的資料X(1)，…，X(n)，即一個時間序列，視窗為b的移動平均值可以定義為

從下圖（圖2）中可以看出，移動平均值是原始資料的平滑版本，平滑程度取決於頻寬。頻寬越大，函式越平滑。

圖2:視窗頻寬為6、24和42的移動平均;x軸:時間，y軸:搜尋百分比

頻寬的選擇至關重要，但不清楚如何選擇頻寬。 如果頻寬太小，我們可能無法擺脫季節性波動。 如果頻寬太大，我們可能無法捕捉到趨勢。 例如，如果我們選擇頻寬b = 0，則具有原始資料及其季節性。 相反，如果b = n，我們僅獲得所有觀測值的平均值，而看不到任何趨勢。

在此示例中，b = 6個月是“平滑”季節性因素的合理選擇，因為我們計算的是整個年度（13個月）的平均值。 但是，b = 12或b = 18是同等有效的選擇。 根據b的選擇，我們將更多的權重賦予與時刻t（b = 12）相同的季節或相反季節（b = 6或b = 18）的觀測值。 減輕此問題的可能解決方案是為觀察值賦予不同的權重，從而計算加權平均值而不是簡單平均值。

理論上講，接近時間t的觀測比更遠的觀測更重要，並且權重更大。 這使我們可以選擇更大的頻寬b，同時將更多的權重分配給更重要的觀察。 例如，對於時間t的觀測，我們可以賦予權重1；對於距離i為i = 1，…，11的觀測，我們可以賦予1-i / 12權重；對於i = 1，...， b-1和頻寬b（請參見圖3）。

圖3：頻寬為6、24和42的加權移動平均線； x軸：時間，y軸：搜尋百分比

這是核估計背後的基本思想：對不同距離的觀測值賦予不同的權重。

權重(1-i/b) 的上述選擇相當隨意，其他權重也可以理解。 我們可以將加權移動平均值歸納如下：令b為頻寬，K為核，即具有某些其他性質的函式，例如非負[K（x）≥0]對稱性[K（x）= K （-x）]和歸一化[K上的積分為1]。 許多核心都支援間隔[-1，1]，這意味著它們在間隔之外為0。

核心的重要示例是| x |的Epanechnikov 核心K(x) = 3/4 (1-x²) 當 |x| ≤ 1 和高斯核K（x）= 1 / sqrt（2π）exp（-x²/ 2）。

基於核心K和頻寬b，我們將NWE（Nadaraya-Watson Estimator）定義為：

在圖4中，我們看到具有高斯核且頻寬b = 12的NWE。 核心和頻寬的選擇仍然很重要，但是具有頻繁使用核心的估計器（例如Epanechnikov，Quartic或Gaussian）在頻寬選擇方面比移動平均估計器更健壯。

注意，簡單移動平均估計量本身是具有統一核K（x）= 1/2的核估計量。 加權移動平均估計器也是一個核心估計器，wiki裡面有一個和函式的整理列表，有興趣的可以仔細檢視https://en.wikipedia.org/wiki/Kernel(statistics)#Kernelfunctionsincommon_use

圖4：具有高斯核和頻寬12的NEW； x軸：時間，y軸：百分比搜尋

進一步說明：首先，通常基於重新定標的時間（即i / n而不是i）來定義NEW，並且公式也會相應變化。 其次，對於距離變化的觀測，可以更普遍地定義估算器。

核密度估計

讓我們考慮另一個例子。 由於某種原因，你可能會對德國的汽油價格感興趣。 因此，你可以上網搜尋所有14,000個加油站的當前價格。 圖5中是該資料的常見表示形式：直方圖。 直方圖顯示汽油價格的分佈。 塊的頻寬是關鍵引數，不同的選擇會導致不同的直方圖（類似於移動平均估計器的不同頻寬）。

圖5：直方圖顯示德國（05/12/2020）分別有10個和50個垃圾箱的天然氣價格頻率； x軸：以EUR為單位的汽油價格； y軸：頻率；

如果我們假設天然氣價格的分佈是連續的，我們可能更喜歡估計和視覺化基礎分佈的密度函式。 這個想法的原理很簡單：對於天然氣價格範圍內的任何一點，我們都會計算價格接近該點的頻率，並給予較接近的價格更多的權重，而向較遠的價格賦予較小的權重。

如果“距離決定權重”是確定正確的， 那麼我們將重點關注這個調節，這就是是核心迴歸背後的想法。

資料X（1），…，X（n）的核密度估計器的定義與NWE非常相似。 給定一個核心K且頻寬h> 0，定義

通常使用與核迴歸情況相同的核函式（例如，高斯，Epanechnikov或Quartic）。 核密度估計可以解釋為提供關於底層資料生成過程的分佈的平滑的直方圖。 核心和頻寬的選擇同樣至關重要（有關不同的估算器，請參見圖6）。

圖6:不同核心(上:Epanechnikov，下:高斯)和不同頻寬(左:0.05，右:0.1)下天然氣價格密度的KDE;x軸:天然氣價格(歐元);軸:頻率

在Python中實現

為了展示核心迴歸，我們使用Google搜尋一詞中的Chocolate一詞的數量，該詞可在Google Trends中下載。 如果儲存為kde_chocolate.csv，則以下Python指令碼將計算NWE並繪製圖4：

import numpy as np 
import pandas as pd 
import matplotlib.pyplot as plt 
from statsmodels.nonparametric.kernel_regression import KernelReg 
 
df = pd.read_csv('kde_chocolate.csv') 
 
n = df.shape[0] 
kde = KernelReg(endog=df['chocolate'], exog=np.arange(n), var_type='c', bw=[12]) 
 
estimator = kde.fit(np.arange(n)) 
estimator = np.reshape(estimator[0],n) 
 
fig, ax = plt.subplots(figsize=(10,5)) 
ax.plot(df['month'], df['chocolate'], '-',alpha=0.5) 
ax.plot(df['month'], estimator, '-', color='tab:blue', alpha=0.8) 
ax.set_ylim(bottom=30, top=105) 
ax.set_xticks(['2005-01','2010-01','2015-01','2020-01']) 
 
plt.show()

考慮到德國不同加油站的汽油價格，其中“ stationuuid”和“ e5”列儲存在kdegasdata.csv中，可通過以下指令碼獲得類似於圖6的圖。 可以通過Azure Repo獲取資料（https://dev.azure.com/tankerkoenig/git/tankerkoenig-data）。

import pandas as pd 
import matplotlib.pyplot as plt 
from sklearn.neighbors import KernelDensity 
import numpy as np 
 
# Load data and remove outliers 
df = pd.read_csv('kde_gas_data.csv')[['station_uuid', 'e5']] 
df = df.groupby(['station_uuid']).mean() 
df = df[df['e5'] > 0] 
df = df.sort_values('e5')[round(df.shape[0]*0.001):round(df.shape[0]*0.999)] 
 
# Calculate and plot KDEs 
fig, ax = plt.subplots(2, 2, sharex=True, sharey=True, figsize=(10,6)) 
X = df['e5'][:,np.newaxis] 
X_plot = np.linspace(1,1.6,1000)[:,np.newaxis] 
zeros = np.zeros(1000) 
 
kde = KernelDensity(kernel='epanechnikov', bandwidth=0.05).fit(X) 
log_dens = kde.score_samples(X_plot) 
ax[0, 0].fill_between(X_plot[:, 0], zeros, np.exp(log_dens), fc='#AAAAFF') 
 
kde = KernelDensity(kernel='epanechnikov', bandwidth=0.1).fit(X) 
log_dens = kde.score_samples(X_plot) 
ax[0, 1].fill_between(X_plot[:, 0], zeros, np.exp(log_dens), fc='#AAAAFF') 
 
 
kde = KernelDensity(kernel='gaussian', bandwidth=0.05).fit(X) 
log_dens = kde.score_samples(X_plot) 
ax[1, 0].fill_between(X_plot[:, 0], zeros, np.exp(log_dens), fc='#AAAAFF') 
 
kde = KernelDensity(kernel='gaussian', bandwidth=0.1).fit(X) 
log_dens = kde.score_samples(X_plot) 
ax[1, 1].fill_between(X_plot[:, 0], zeros, np.exp(log_dens), fc='#AAAAFF') 
 
plt.show()

總結

在兩種情況下（核心迴歸和核心密度估計），我們都必須選擇核心和視窗的頻寬。 由於常用的核心具有相似的形狀（請參見圖7），因此頻寬的選擇更為關鍵。 關於頻寬選擇，有大量文獻。 在應用程式中，兩種方法佔主導地位：拇指法則頻寬估計器（例如Silverman的經驗法則）和交叉驗證（尤其是用於核心迴歸）。

圖7：常用的核心函式：高斯（左上），Epanechnikov（右上），四次（左下），均勻（右下）

由於大多數核心函式是平滑的（或至少是連續可微的），因此得出的估計量也是平滑的。 因此，呈現形式的核迴歸僅對於足夠平滑的迴歸/密度函式才有意義。

NWE是一類更廣泛的非引數估計器的特例，即區域性多項式估計器。

基於核心的估計器非常有用，並適用於許多情況。 甚至k最近鄰演算法也可以被視為具有統一核心和可變頻寬的基於核心的估計器。 正如我們可以採用其他核心一樣，這種範例為我們提供了新的可能性，從而為不同的點賦予了不同的權重。

我相信基於核的估計的概念對資料科學家來說很重要，希望你能建立一些直觀的理解。

作者：Florian Heinrichs

deephub翻譯組