Estymacja gęstości jądrowej: od teorii do praktyki

Próbowałeś kiedyś zwizualizować rozkład, a histogram zmieniał kształt za każdym razem, gdy zmieniałeś szerokość przedziałów?

Zwykle wygląda to tak. Ustawiasz 10 przedziałów i widzisz gładką krzywą. Potem przełączasz na 30 i pojawia się kilka pików. Dane się nie zmieniły, ale różna liczba przedziałów daje różne interpretacje. To największy problem z histogramami: nie pokazują rozkładu, tylko jedną z jego wersji. A ta wersja jest zależna od parametru, który ustawiasz dość arbitralnie.

KDE podchodzi do sprawy inaczej. Zamiast dzielić dane na przedziały, nakłada na każdy punkt danych małą, gładką krzywą i sumuje je wszystkie. Dzięki temu otrzymujesz jedną, ciągłą estymację ukrytego rozkładu.

W tym artykule poznasz intuicję stojącą za KDE, przejdziesz przez wzór, zrozumiesz, jak szerokość pasma steruje gładkością, oraz zobaczysz praktyczne przykłady w Pythonie i R.

Dopiero zaczynasz z histogramami? Oto kompleksowy przewodnik po histogramach częstości, który pomoże ci zacząć.

Czym jest estymacja gęstości jądrowej?

Estymacja gęstości jądrowej to nieparametryczna metoda szacowania funkcji gęstości prawdopodobieństwa dla zbioru danych.

Właśnie ta nieparametryczność ją wyróżnia.

W metodach parametrycznych zakładasz, że dane podążają za konkretnym rozkładem — normalnym, wykładniczym — a następnie dopasowujesz parametry. Jeśli to założenie jest błędne, model też jest błędny. KDE nie robi takich założeń. Pozwala danym „mówić same za siebie” i buduje estymację rozkładu bezpośrednio z obserwacji.

Wynik to gładka krzywa pokazująca, gdzie wartości najprawdopodobniej się znajdą — i jak bardzo. Wysokie punkty krzywej oznaczają gęste obszary. Niskie — rzadkie.

Dlaczego korzystać z estymacji gęstości jądrowej?

Histogramy to domyślne narzędzie do wizualizacji rozkładów, ale mają problem: kształt, który widzisz, zależy od liczby przedziałów, które wybierzesz. A ta liczba to twój wybór. Dwie osoby mogą spojrzeć na ten sam zbiór danych i dojść do zupełnie różnych wniosków tylko dlatego, że wybrały inną liczbę przedziałów.

W KDE, zamiast wciskać dane w przedziały, dostajesz gładką, ciągłą krzywą, która nie zmienia się w oparciu o arbitralny parametr ustawiony z góry.

To przydaje się w kilku sytuacjach:

• Wizualizacja ogólnego kształtu rozkładu bez stronniczości od szerokości przedziałów
• Wykrywanie wzorców, takich jak wiele pików, skośność czy „grube ogony”
• Porównywanie rozkładów między grupami na tym samym wykresie
• Uzyskanie czystszego obrazu danych przed modelowaniem

Intuicja stojąca za KDE

Dla każdego punktu danych nakładasz na niego małą, gładką krzywą. Ta krzywa to jądro. Potem po prostu sumujesz wszystkie te pojedyncze krzywe w jedną.

W efekcie dostajesz jedną gładką krzywą, która pokazuje gęstość twoich danych. Tam, gdzie punkty się grupują, wiele jąder nakłada się i sumuje, więc krzywa rośnie. Tam, gdzie danych jest mało, jądra ledwo się nakładają i krzywa pozostaje niska. Każdy punkt wnosi równy wkład do końcowej estymacji.

Wyobraź sobie, że zapisałeś wyniki końcowego egzaminu w klasie. Zamiast wrzucać je do histogramu, KDE nakłada małą gładką krzywą na każdy wynik. Tam, gdzie wyniki się skupiają — powiedzmy, wokół 70–75 — krzywe się sumują i estymata rośnie. Pojedynczy uczeń z wynikiem 95 dodaje tylko niewielkie wybrzuszenie w ogonie.

Poniższa wizualizacja pokazuje dokładnie to. Większość uczniów zdobyła wynik bliski średniej, a jeden znacznie wyższy:

Zwizualizowane KDE

Wzór na estymację gęstości jądrowej

Wzór KDE wygląda groźniej, niż jest w rzeczywistości.

Wzór KDE

Oto co oznaczają poszczególne elementy:

• n — liczba punktów danych

• x_i — pojedyncze punkty danych w twoim zbiorze

• K — funkcja jądra, czyli gładka krzywa nakładana na każdy punkt

• h — szerokość pasma (bandwidth), czyli to, jak szerokie jest każde jądro

• x — punkt, w którym oceniasz gęstość

Po ludzku: dla dowolnego punktu x sprawdź, jak blisko są od niego wszystkie punkty danych x_i, zważ tę bliskość funkcją jądra K i uśrednij wynik po wszystkich n punktach. Zrób to dla każdego x w całym zakresie, a otrzymasz pełną krzywą gęstości.

Szerokość pasma h znajduje się w mianowniku ułamka wewnątrz K. Mniejsze h zwęża jądro, więc na estymatę wpływają tylko bardzo bliskie punkty. Większe h rozszerza ten wpływ. Więcej na ten temat później w artykule.

Czym jest funkcja jądra?

Jądro to gładka krzywa, którą nakładasz na każdy punkt danych. Określa, jak wpływ tego punktu rozchodzi się na sąsiednie wartości.

Każde jądro jest wyśrodkowane w punkcie danych i przypisuje wagi w zależności od odległości. Punkty blisko środka dostają wysoką wagę. Dalekie — niską albo żadną. Dokładny kształt ważenia zależy od wybranego jądra.

Najczęściej wybiera się jedno z trzech:

• Gaussowskie: krzywa dzwonowa. Nadaje wagę każdemu punktowi w zbiorze, a wpływ maleje wraz z odległością. To domyślna opcja w większości bibliotek i w praktyce sprawdza się dobrze
• Jednostajne: płaski prostokąt. Każdy punkt w ustalonej odległości dostaje równą wagę, a wszystko poza nią — żadną. Proste, ale daje mniej gładkie estymaty
• Epanechnikova: kształt paraboliczny, matematycznie optymalny pod względem średniego błędu kwadratowego. Jest wydajniejsze niż gaussowskie, ale rzadko daje widoczną różnicę w praktyce.

W większości przypadków wybór jądra niewiele zmienia. Dwa różne jądra zastosowane do tych samych danych i z tą samą szerokością pasma dadzą niemal identyczne krzywe. Znacznie bardziej liczy się szerokość pasma — i temu przyjrzymy się dalej.

Rola szerokości pasma w KDE

Szerokość pasma to pojedynczy parametr, który ma największy wpływ na wynik KDE — większy niż wybór jądra.

Steruje tym, jak szerokie jest każde jądro. Wąskie jądro „zbiera” wpływ tylko z pobliskich punktów. Szerokie rozciąga ten wpływ na dużo większy zakres. W efekcie dostajesz krzywą, która albo ściśle podąża za danymi, albo je wygładza.

Mała szerokość pasma

Mała szerokość pasma sprawia, że każde jądro jest ciasne i wąskie. Estymata gwałtownie reaguje na każdy punkt danych, co oznacza, że wychwytuje prawdziwą strukturę, ale też szum.

W praktyce wygląda to jak „kolczasta” krzywa z wieloma małymi pikami. Część z nich pokazuje realne skupiska w danych. Inne to artefakty zbyt słabego wygładzania. Trudno odróżnić jedno od drugiego — i w tym problem.

KDE z małą szerokością pasma

Duża szerokość pasma

Duża szerokość pasma szeroko rozciąga każde jądro. Sąsiednie jądra się nakładają, a końcowa krzywa jest gładka.

Jeśli wygładzanie jest zbyt silne, tracisz prawdziwą strukturę. Dwa odrębne skupiska mogą zlać się w jedną krzywą. Rozkład z grubym ogonem może wyglądać symetrycznie. Wizualizacja może coś przed tobą ukrywać.

KDE z dużą szerokością pasma

Jak znaleźć właściwą szerokość pasma

Nie ma jednej, zawsze poprawnej szerokości pasma. Celem jest znalezienie wartości wystarczająco gładkiej, by odfiltrować szum, ale nie tak gładkiej, by zatarła prawdziwe wzorce.

Większość bibliotek dobiera ją automatycznie. Reguła kciuka Silvermana jest najczęstsza. Wybiera szerokość pasma na podstawie liczebności próby i odchylenia standardowego danych. Dobrze działa dla w przybliżeniu normalnych rozkładów, ale może nadmiernie wygładzać rozkłady wielomodalne.

Jeśli nie masz pewności, spróbuj kilku wartości i porównaj krzywe. Różnice wiele powiedzą ci o danych.

KDE vs. histogram

Zarówno histogram, jak i KDE pokazują rozkład danych — ale robią to w zupełnie inny sposób.

Histogram

Histogram dzieli dane na dyskretne przedziały i liczy, ile punktów wpada do każdego. Jest szybki, intuicyjny i łatwy do wyjaśnienia nietechnicznej publiczności.

Problemem jest wrażliwość na liczbę przedziałów. Zmieniając ich liczbę, zmieniasz kształt. Nie ma obiektywnie poprawnej liczby, więc dwie osoby mogą wyciągnąć różne wnioski z tych samych danych tylko przez ten wybór.

Histogramy dają też schodkowy, nieciągły kształt. To w porządku na szybki rzut oka, ale może zaciemniać prawdziwy rozkład.

KDE

KDE daje gładką, ciągłą krzywą bez udziału przedziałów. Lepiej pokazuje faktyczny kształt rozkładu — takie rzeczy jak skośność, wiele pików czy grube ogony, które histogram może pominąć lub zniekształcić w zależności od wyboru przedziałów.

Minusem jest to, że KDE wprowadza własny parametr — szerokość pasma — i wymaga więcej obliczeń. Jest też mniej intuicyjne do wyjaśnienia, bo oś y pokazuje gęstość prawdopodobieństwa, a nie zliczenia, co może mylić osoby niezaznajomione z pojęciem.

Kiedy używać którego

Użyj histogramu, gdy potrzebujesz szybkiego, łatwego do interpretacji podsumowania danych lub gdy odbiorcy nie znają estymacji gęstości. Użyj KDE, gdy liczy się kształt rozkładu — na przykład porównując grupy lub próbując wykryć wielomodalność w danych.

Histogram w porównaniu z KDE

W praktyce często używa się ich razem: histogram dla zliczeń, a na nim nałożona krzywa KDE dla kształtu.

Estymacja gęstości jądrowej w Pythonie

Python oferuje kilka sposobów obliczenia i narysowania KDE — w zależności od tego, czy potrzebujesz szybkiego wykresu, czy większej kontroli nad samą estymatą.

Wizualizacja KDE z seaborn

Najszybszy sposób na wykres KDE to seaborn.kdeplot(). Wystarczy tyle:

import seaborn as sns
import numpy as np

np.random.seed(42)
scores = np.concatenate([
    np.random.normal(65, 4, 60),
    np.random.normal(80, 3, 40)
])

sns.kdeplot(scores, bw_adjust=1)

KDE w seaborn

Parametr bw_adjust skaluje automatycznie dobraną szerokość pasma. Wartości poniżej 1 zacieśniają krzywą, powyżej 1 — wygładzają. To mnożnik na tym, co seaborn dobiera wewnętrznie, więc nie musisz sam podawać surowej szerokości pasma.

Oś y pokazuje gęstość prawdopodobieństwa, a nie zliczenia. Krzywa mówi, jak bardzo dany obszar jest prawdopodobny względem reszty rozkładu. Wyżej oznacza większą koncentrację danych.

Obliczanie KDE ze scipy

Jeśli potrzebujesz rzeczywistych wartości gęstości, nie tylko wykresu, użyj scipy.stats.gaussian_kde. Zwraca obiekt wywoływalny, który możesz oceniać w dowolnym punkcie.

from scipy.stats import gaussian_kde
import numpy as np
import matplotlib.pyplot as plt

kde = gaussian_kde(scores, bw_method="scott")
x = np.linspace(45, 100, 500)
density = kde(x)

plt.plot(x, density)

bw_method="scott" używa reguły Scotta do automatycznego doboru szerokości pasma. To dobry domyślny wybór w większości przypadków. Możesz też przekazać skalar, by ustawić szerokość ręcznie.

KDE w scipy i matplotlib

Estymacja gęstości jądrowej w R

W R KDE jest wbudowane w sam język. Nie potrzebujesz dodatkowych pakietów.

Obliczanie KDE funkcją density()

Funkcja density() przyjmuje wektor liczbowy i zwraca obiekt KDE.

set.seed(42)
scores <- c(rnorm(60, mean = 65, sd = 4),
            rnorm(40, mean = 80, sd = 3))

kde <- density(scores, bw = "SJ")

Argument bw steruje wyborem szerokości pasma. "SJ" używa metody Sheathera-Jonesa, która lepiej radzi sobie z rozkładami wielomodalnymi niż domyślna. Możesz też przekazać wartość liczbową, by ustawić szerokość ręcznie.

Wynik to obiekt listy z dwoma kluczowymi elementami:

• kde$x: sekwencja punktów, w których oceniono gęstość
• kde$y: odpowiadające im wartości gęstości

Rysowanie krzywej

Po prostu przekaż wynik bezpośrednio do plot().

plot(kde,
     main = "Exam Score Distribution",
     xlab = "Score",
     ylab = "Density")

KDE narysowane w R

Aby nałożyć KDE na histogram, najpierw użyj hist() z freq = FALSE, a następnie dodaj krzywą przez lines().

hist(scores, freq = FALSE, main = "Histogram + KDE", xlab = "Score")
lines(kde, col = "blue", lwd = 2)

Histogram z KDE w R

freq = FALSE skaluje histogram do gęstości, aby słupki i krzywa miały tę samą oś y.

Zalety i ograniczenia KDE

KDE to naprawdę przydatna wizualizacja, ale jak wszystko — ma kompromisy, o których warto wiedzieć, zanim zastąpisz nim histogramy.

Zalety

Największą zaletą jest to, że KDE nie zakłada rozkładu twoich danych. Nie musisz z góry decydować, czy dane są normalne, wykładnicze czy jakiekolwiek inne. Kształt wynika bezpośrednio z danych, dzięki czemu KDE jest na tyle elastyczne, by poradzić sobie z rozkładami wielomodalnymi i innymi, które nie pasują do standardowych form parametrycznych.

Wynik to gładka, ciągła krzywa, a nie schodkowe przybliżenie. Dzięki temu łatwiej dostrzec wzorce — takie jak wiele pików czy długie ogony — które histogram może ukryć w zależności od wyboru przedziałów.

A ponieważ KDE działa na surowych danych i nie wymaga najpierw dopasowania modelu, to dobry pierwszy krok w każdej analizie eksploracyjnej.

Ograniczenia

Wybór szerokości pasma to główna słabość. Jeśli się pomylisz, estymata albo goni szum, albo wygładza prawdziwe wzorce w danych. Automatyczne metody, jak reguła Silvermana, dobrze działają przy w przybliżeniu normalnych danych, ale przy złożonych rozkładach mogą wprowadzać w błąd. Często warto ręcznie sprawdzić kilka wartości, zanim zaufasz wynikowi.

Wydajność może stać się problemem w dużej skali. KDE ocenia funkcję jądra dla każdego punktu danych w każdym miejscu oceny, więc liczba obliczeń szybko rośnie wraz z wielkością zbioru. Do większości analiz eksploracyjnych to nie problem, ale przy setkach tysięcy punktów może być wolne.

Efekty brzegowe to subtelniejsza kwestia. Standardowe KDE zakłada, że dane mogą rozciągać się nieskończenie w obie strony. Gdy zmienna ma twardą granicę — np. wartości nie mogą być poniżej zera — estymata „wycieka” z masą prawdopodobieństwa poza dopuszczalny zakres, przez co krzywa sztucznie opada przy krawędziach. Istnieją wersje KDE z korektą brzegów, ale rzadziej są zaimplementowane w standardowych bibliotekach.

Wnioski

KDE daje czystszy sposób spojrzenia na rozkład twoich danych niż histogramy. Nie ma wyboru przedziałów ani założeń parametrycznych — jest tylko gładka krzywa pokazująca to, co faktycznie jest w zbiorze.

Szerokość pasma to jedyny parametr, który naprawdę się liczy. Wypróbuj kilka wartości, porównaj krzywe, skorzystaj z opcji automatycznych i upewnij się, że estymata zgadza się z twoją wiedzą o danych, zanim wyciągniesz wnioski.

Najlepszy sposób na zbudowanie intuicji do KDE to użycie go na prawdziwych danych. Wybierz zbiór, który już znasz, zastosuj KDE i porównaj z histogramem, by zobaczyć, co ci umykało.

Interesuje cię wizualizacja danych? Sprawdź nasz kurs Data Visualization with Seaborn, jeśli korzystasz z Pythona, lub Data Visualization with ggplot2, jeśli używasz R.