Ядерное сглаживание для оценки плотности

Идея aka Краткое содержание

Проблема: нужен непараметрический метод для оценки плотности.

Решение: метод будет основан на локальной оценке плотности в окрестности интересующей точки по известной выборке. Локальная оценка опирается на само определение плотности распределения: $p(x) = \lim_{h \to 0} \frac{1}{2h}P[x - h, x + h]$, где $P[a, b]$ — вероятностная мера отрезка $[a, b]$.

Так и родился на свет один из непараметрических способов оценки плотности распределения — ядерное сглаживание (KDE или Kernel Density Estimation). В отличие от метода гистограмм, блоки (окна), по которым оценивается распределение, не фиксированы, а центрируются по точке-представителю.

Общая формула KDE (для одномерного и многомерного случая) представлена ниже.

Два важных параметра метода: ядро и ширина окна. Выбор ядра в основном влияет на гладкость итогового распределения, но на точность аппроксимации намного большее влияние оказывает второй параметр, поэтому подбор ширины окна является важной и не всегда тривиальной задачей (прибегают к кросс-валидации, различным эвристикам или динамическому выбору ширины окна: см далее), но основное правило приблизительно таково: чем плотнее выборочное распределение, тем уже должно быть окно.

Итак...

Условные обозначения

$N$ — количество объектов в выборке.

$D$ — размер признакового пространства, $D \geq 1$.

$C$ — количество классов.

$\mathbb{X}=\{(x_i, y_i)\}_{i=1}^{N}$ — выборка, $x_i \in \mathbb{R}^D$, $y_i \in \{1, \dots, C\}$. $x \in \mathbb{R}^D$.

$h$ — ширина окна (bandwidth), $h \geq 0$.

$\hat{p}(x)$ — оценка плотности распределения $p(x)$.

$\mathbb{I}[$условие$]$ — равняется 1, если условие выполнено, иначе равняется 0.

$\hat{y}(x)$ — оценка зависимости $y(x)$.

Гистограммы

Недостаток: необходимо фиксировать отрезки, на которые разбивается интервал. Проблема: выбор количества корзинок и ширины корзинок.

Две гистограммы для одной выборки

Ядерное сглаживание

Идея: каждый $x_i$ выборки будет центром блока.

Блок может иметь следующий вид: $\tfrac{1}{2}\cdot\mathbb{I}[|\tfrac{x - x_i}{h}| \leq 1], D = 1$.

Одномерный случай $D=1$

Kernel Density Estimation (KDE, локальная непараметрическая оценка Парзена-Розенблатта) — $\hat{p}(x)=\tfrac{1}{N\cdot h}\cdot\sum_{i=1}^{N}K(\tfrac{x-x_i}{h})$, $K(\cdot)$ — ядро, чётная и нормированная функция: $\int_{-\infty}^{+\infty}K(u)du=1$. Следствие: $\hat{p}(x)$ обладает той же степенью гладкости, что и ядро $K(\cdot)$.

Виды ядер

• Прямоугольное ядро (tophat kernel): $K(u)=\tfrac{1}{2}\cdot\mathbb{I}[|u| \leq 1]$ соответствует эмпирической оценке плотности (доля точек выборки, лежащих внутри отрезка $[x-h, x+h]$). Одно из простейших ядер, но не учитывает расстояние между объектами, а также итоговое распределение не будет являться непрерывным.
• Гауссово ядро: $K(u)=\tfrac{1}{\sqrt{2\cdot\pi}}\cdot\exp{\tfrac{-u^2}{2}}$
• Ядро Епанечникова: $K(u)=\tfrac{3}{4}\cdot\max\{1 - u^2, 0\}$
• Треугольное ядро: $K(u)=\max\{1 - |u|, 0\}$
• Косинусное ядро: $K(u)=\tfrac{\pi}{4}\cos(\tfrac{\pi\cdot u}{2})\cdot\mathbb{I}[|u| \leq 1]$
• Экспоненциальное ядро: $K(u)=\tfrac{1}{2}\cdot\exp(-|u|)$
• Квартическое ядро: $K(u)=\tfrac{15}{16}\cdot(1 - u^2)^2\cdot\mathbb{I}[|u| \leq 1]$

Состоятельность оценки $\hat{p}(x)$

Оценка $\hat{p}(x)$ состоятельна, если $\forall x~ \mathbb{E}[(\hat{p}(x) - p(x))^2] \xrightarrow[]{N \rightarrow \infty} 0$.

Достаточные условия состоятельности оценки $\hat{p}(x)$:

$\lim_{N \to \infty} h(N) = 0$, $\lim_{N \to \infty} N \cdot h(N) = \infty$

$\int_{-\infty}^{+\infty}|K(u)|du < \infty$, $\int_{-\infty}^{+\infty}K(u)du=1$, $\sup_{u \in R} K(u) < \infty$, $\lim_{u \to \infty} |u\cdot K(u)| = 0$

Многомерный случай $D\geq2$

$\hat{p}(x)=\tfrac{1}{N\cdot h^D}\cdot\sum_{i=1}^{N}K(\tfrac{x-x_i}{h}), x \in \mathbb{R}^D$

Виды ядер

• Гауссово ядро: $K(u)=\tfrac{1}{(2\cdot\pi)^{\tfrac{D}{2}}}\cdot\exp(\tfrac{-u^T\cdot u}{2})$
• Ядро Епанечникова: $K(u)={\tfrac {(D+2)!!}{2^{\left\lceil {\tfrac {D+2}{2}}\right\rceil }\cdot \pi ^{\left\lfloor {\tfrac {D}{2}}\right\rfloor }}}\cdot \max\{1-u^{T}\cdot u,0\}$
• Произведение одномерных ядер: $K(\tfrac{x - x_i}{h})=\prod_{d=1}^{D} K_d(\tfrac{x^d - x_i^d}{h}), x=(x^1, \dots, x^D), x_i=(x_i^1, \dots, x_i^D)$

Зависящие от метрики $\rho(\cdot, \cdot)$ ядра

Примечание: в случае метрики, отличной от евклидовой, коэффициенты перед ядрами могут быть другими.

$\hat{p}(x)=\tfrac{1}{N\cdot h^D}\cdot\sum_{i=1}^{N}K(\tfrac{\rho(x, x_i)}{h})$

• Гауссово ядро: $K(\tfrac{\rho(x, x_i)}{h})=\tfrac{1}{(2\cdot\pi)^{\tfrac{D}{2}}}\cdot\exp(\tfrac{-\rho^2(x, x_i)}{2\cdot h^2})$
• Ядро Епанечникова: $K({\tfrac {\rho (x,x_{i})}{h}})={\tfrac {(D+2)!!}{2^{\left\lceil {\tfrac {D+2}{2}}\right\rceil }\cdot \pi ^{\left\lfloor {\tfrac {D}{2}}\right\rfloor }}}\cdot \max\{1-{\tfrac {\rho ^{2}(x,x_{i})}{h^{2}}},0\}$

Выбор ширины окна (bandwidth)

При $h \to 0$ плотность концентрируется вблизи точек выборки, $\hat{p}(x)$ претерпевает резкие скачки. При $h \to \infty$ более гладкая плотность, происходит вырождение в константу. При построении KDE ширина окна $h$ важнее, чем функция ядра $K(\cdot)$, так как тип ядра влияет на гладкость, а не на точность аппроксимации.

Стратегия выбора: чем более плотное распределение объектов выборки, тем меньше должно быть $h$

• Постоянное значение h, примеры стратегий:
  • $h=\tfrac{1}{N}\cdot\sum_{i=1}^{N}d_{iK}, d_{iK}$ — расстояние от $x_i$ до $K$-го ближайшего соседа ($K$ можно вычислять по скользящему контролю).
  • $h$ вычисляется по скользящему контролю (Leave-one-out, например), можно найти по максимальному правдоподобию на отложенной выборке (поиск максимального значения правдоподобия производится по заданному списку значений $h$).
• Переменное значение $h(x)$, например: $h(x)$ — расстояние от $x$ до $K$-го ближайшего соседа ($K$ можно найти по скользящему контролю).

Метод Парзеновского окна

Метод Парзеновского окна — метод байесовской классификации, основанный на непараметрическом восстановлении плотности по имеющейся выборке.

Оценка условной плотности $p(x|y)$ через KDE ($y \in \{1, \dots, C\}, D \geq 1$):

$p(x|y)=\tfrac{1}{N_y\cdot h^D}\cdot\sum_{i:y_i=y}K(\tfrac{\rho(x, x_i)}{h}), N_y$ — число объектов класса $y, \rho(\cdot, \cdot)$ — метрика.

Байесовское решающее правило даёт следующий классификатор: $\hat{y}(x)=\arg\max_{y}p(y|x)=\arg\max_{y}\tfrac{p(x|y)\cdot p(y)}{p(x)}=\arg\max_{y}p(x|y)\cdot p(y)$

Оценим $p(x|y)$ с помощью KDE, $p(y)$ как $\tfrac{N_y}{N}$:

$\hat{y}(x)=\arg\max_{y}\tfrac{1}{N_y\cdot h^D}\cdot\sum_{i:y_i=y}K(\tfrac{\rho(x, x_i)}{h})\cdot\tfrac{N_y}{N}=\arg\max_{y}\sum_{i:y_i=y}K(\tfrac{\rho(x, x_i)}{h})$

Преобразование метода Парзеновского окна в метод ближайших соседей

Обозначим $h(x)=\rho(x, x_{i(K)}), i(K)$ — индекс $K$-го ближайшего соседа для $x$, $K(u)=\mathbb{I}[|u| \leq 1]$. Тогда:

$\hat{y}(x)=\arg\max_{y}\sum_{j:y_j=y}\mathbb{I}[\rho(x, x_j) \leq \rho(x, x_{i(K)})]=\arg\max_{y}\sum_{j:\rho(x, x_j) \leq \rho(x, x_{i(K)})}\mathbb{I}[y_j=y]=$$=\arg\max_{y}\sum_{j=1}^{K}\mathbb{I}[y_{i(j)}=y]$

Ссылки

• лекция Виктора Китова по ядерной оценке плотности
• К. В. Воронцов Математические методы обучения по прецедентам