Tegyük fel, hogy rendelkezésre állnak nem egyenletes eloszlású adatminták. Bár az eloszlást nem ismerjük, de annyit tudunk, hogy folytonos. Tudni szeretnénk, hogy mi az eloszlás módusza.
Folytonos eloszlások móduszának, vagyis a sűrűségfüggvény maximumának ismerete a hétköznapi vagy üzleti életben fontos lehet. Például egy városban a nappali hőmérséklet eloszlásának maximuma nyáron a légkondicionálás vagy szabadtéri rendezvények időzítése szempontjából lehet hasznos. Egy szupermarketben a napi vásárlási összegek eloszlásának csúcsa a kosárérték-optimalizálás szempontjából érdekes, például promóciós vagy árazási stratégiák kialakításához.
Bármilyen célból is vagyunk kíváncsiak a móduszra, az adatmintákból közvetlenül nem tudjuk azt meghatározni, mert az eloszlás nem diszkrét, hanem folytonos. Az is kérdés lehet, hogy van-e egyáltalán módusza az eloszlásnak.
A válasz: igen. Minden folytonos, nem egyenletes eloszlásnak van legalább egy módusza. Vagyis a sűrűségfüggvénynek létezik legalább egy globális maximuma, sőt akár több lokális maximuma is lehet.
Ez a megállapítás szemlélet alapján is belátható. Ugyanis a sűrűségfüggvény egy olyan folytonos, nemnegatív értékű függvény, amelynek görbe alatti területe (integrálja) véges (konkrétan 1). Ha az eloszlás nem egyenletes, akkor a sűrűségfüggvény nem lehet konstans, így:
- vagy csak egy korlátos és zárt intervallumban vesz fel pozitív értékeket, másutt 0,
- vagy a függvényértéke a pozitív/negatív végtelenben nullához tart, különben az integrálja nem lenne véges.
Mindkét esetben szükségszerű, hogy a sűrűségfüggvénynek legyen maximuma.
Természetesen a fenti állítás egzakt matematikai úton is bizonyítható, de ehhez felsőbb szintű matematikai ismeretek kellenek (az ezt ismerőknek a Weierstrass tételre kell emlékezni).
Ez a kis matematikai előzetes azért szükséges, hogy esetleg ne hiába keressünk móduszt ott, ahol nincs. De most már tudjuk, hogy kell lenni.
Most, hogy ezt tisztáztuk, már ténylegesen megkereshetjük az eloszlás móduszát egy algoritmikus eljárással.
Az egyszerűség kedvéért tételezzük még fel, hogy a sűrűségfüggvénynek csak egyetlen maximuma van, vagyis az eloszlás egymóduszú (unimodális). Ez nem túl szűkítő korlátozás, mert a nevezetes eloszlások többsége (pl. háromszög, béta, gamma, Weibull, normális, lognormális) unimodális. De az alábbiakban bemutatott módszer kis módosítással kiterjeszthető többmóduszú (multimodális) folytonos eloszlásokra is.
A becslési eljárás két fő lépésből áll:
- A minták alapján meghatározunk egy közelítő sűrűségfüggvényt
- Ennek megkeressük a maximumhelyét, ami az eloszlás becsült módusza lesz.
A sűrűségfüggvény becslésére a magfüggvényes, vagy más szóhasználattal kernelalapú sűrűségbecslési módszert (kernel density estimation, KDE) alkalmazzuk. Ennek elvéről és példákkal illusztrált működéséről a „Magfüggvényes valószínűségi sűrűségbecslés” című korábbi cikkben már volt szó, ezért e módszert itt nem ismertetjük.
A maximumkeresésre az intervallumharmadolós keresési eljárást (ternary search) alkalmazzuk. Ez hasonlít a klasszikus intervallumfelezéses módszerhez, de ennél az eljárásnál – ahogy a neve is utal rá – a keresés során az intervallumot mindig három részre osztjuk. Részletes leírása megtalálható például a „Ternary search” Wikipédia szócikkben, ahol Python nyelvű megvalósítás is szerepel. Ez azonban hibás eredményt adhat, ha a harmadolópontokban a függvényértékek megegyeznek. Ezért mi egy javított változatot használunk, amely kiküszöböli ezt a problémát. A módosítás lényege, hogy egyenlőség esetén a vizsgált intervallumot a harmadolóppontok által határoltra szűkítjük és innen folytatjuk a keresést.
Ami a konkrét megvalósítást illeti, definiálunk egy olyan KDEEstimator nevű osztályt, amelynek példánya fogadja a mintákat, és ezek alapján nem csak a becsült sűrűségfüggvényt tudja szolgáltani egy metódushívással, hanem ennek alapján egy másik metódussal a becsült móduszt is. A magfüggvényes sűrűségbecsléshez alkalmazandó paramétereket a példányosításkor a konstruktorban lehet megadni, de azokat később meg lehet változtatni. Ez lehetővé teszi, hogy adott mintasorozathoz több fajta kernellel is végezzünk számítást. Az osztály definíciója alább látható.
|
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 |
# Python 3.13+ from random import triangular, betavariate, gammavariate, weibullvariate, normalvariate, lognormvariate from statistics import kde, stdev from typing import Callable, Sequence from matplotlib import pyplot as plt from math import isclose type RealNum = int | float class KDEEstimator: def __init__(self, kernel: str, bandwidth_adjustment: RealNum = 1): self._samples = () self.kernel = kernel self.bandwidth_adjustment = bandwidth_adjustment @property def samples(self) -> Sequence[RealNum]: return self._samples @samples.setter def samples(self, data: Sequence[RealNum]) -> None: self._samples = tuple(data) def set_kde_params(self, kernel: str, bandwidth_adjustment: RealNum = 1) -> None: """Beállítja a kernel súrúségbecslés paramétereit. Paraméterek: kernel: A használandó kernel típusa. Megadható értékek: "normal", "gauss", "logistic", "sigmoid", "rectangular", "uniform", "triangular", "parabolic", "epanechnikov", "quartic", "biweight", "triweight", "cosine" bandwidth_adjustment: A számított sávszélességhez alkalmazott korrekció szorzótényező. """ self.kernel = kernel self.bandwidth_adjustment = bandwidth_adjustment def estimated_pdf(self) -> Callable: """Valószínűségsűrűség-függvény (PDF) becslése kernel sűrűségbecsléssel. Az alkalmazott kernel függvény a self.kernel értéke lesz. A sávszélesség automatikusan számolódik a Scott-féle szabály alapján, ami a self.bandwidth_adjustment értékével módosul. Visszatérési érték a becsült sűrűségfüggvény. """ # A sávszélesség automatikus kiszámítása a Scott-féle "ökölszabály" szerinti képlettel. bandwidth = 1.06 * stdev(self.samples) * len(self.samples) ** (-1 / 5) # A sávszélesség kiigazítása. bandwidth = bandwidth * self.bandwidth_adjustment return kde(self.samples, bandwidth, self.kernel) @staticmethod def find_max(func: Callable[[RealNum], RealNum], left: RealNum, right: RealNum) -> float: """Visszaadja a [left, right] intervallumban egyetlen maximummal rendelkező func függvény maximumának helyét.""" # Maximumkeresés intervallum harmadolós módszerrel (ternary search). while not isclose(right, left): left_third = left + (right - left) / 3 right_third = right - (right - left) / 3 if func(left_third) < func(right_third): left = left_third elif func(left_third) > func(right_third): right = right_third else: left = left_third right = right_third return (left + right) / 2 def find_mode(self) -> float: """Visszaadja a becsült súrúségfüggvény móduszát.""" sample_min, sample_max = min(self.samples), max(self.samples) pdf_est = self.estimated_pdf() return self.find_max(pdf_est, sample_min, sample_max) |
A becslőobjektum működésének teszteléséhez készítünk egy függvényt random_samples() néven, amely a vizsgálandó mintasorozatot állítja elő. Ennek első argumentuma egy olyan hívható objektum, amely minden egyes híváskor egy bizonyos folytonos eloszlásból egy véletlen értéket szolgáltat. Ennek paramétereit lehet felsorolni a random_samples() első argumentuma után. Utolsó argumentumként lehet megadni, hogy hány minta legyen generálva. Ennek alapértelmezett értékét igény szerint állíthatjuk be. A függvény definíciója:
|
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 |
def random_samples(sampler_func: Callable, *parameters, n: int = 200) -> list[float]: """Véletlenszerű minták generálása egy adott mintavételi függvénnyel. Paraméterek: sampler_func: függvény, amely egyetlen mintát (valós számot) generál egy folytonos valószínűségi eloszlásból a megadott paraméterek alapján (pl. a szabványos könyvtár random moduljának valamely véletlenszám generáló függvénye) *parameters: A sampler_func függvénynek átadandó paraméterek. n: A generálandó minták száma. Visszatérési érték: A generált minták listája. """ return [sampler_func(*parameters) for _ in range(n)] |
A becsült sűrűségfüggvényt egy grafikonon vizuálisan is meg kívánjuk jeleníteni, amin bejelöljük a módusz helyét és kiírjuk a becsült értékét. Ezt valósítja meg a következő függvény:
|
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 |
def plot_pdf_with_mode(kde_mode_estimator: KDEEstimator, low_end=-1, high_end=1, title=''): """Grafikonon ábrázolja az első argumentum becslője által szolgáltatott sűrűségfüggvényt a low_end és high_end értékek között és kiírja a becsült módusz értékét. Az utolsó opcionális argumentummal a diagram címszövege adható meg. """ pdf_est = kde_mode_estimator.estimated_pdf() mode = kde_mode_estimator.find_mode() num_points = 10000 delta = (high_end - low_end) / (num_points - 1) xs = [low_end + n * delta for n in range(0, num_points)] pdf_curve = [pdf_est(x) for x in xs] plt.plot(xs, pdf_curve) # Függőleges vonal rajzolása a módusz értékénél. plt.axvline(x=mode, color='red', linestyle='--', label=f'Módusz: {mode:.2f}') plt.title(title) plt.xlabel('X') plt.ylabel('Sűrűség') plt.legend() plt.grid() plt.show() |
A tesztelést többfajta ismert, nevezetes eloszlás mintáival végezzük. Tesszük ezt azért, mert ekkor előre tudhatjuk az elméleti módusz értékét, és így meg tudjuk ítélni a móduszbecslés pontosságát különböző esetekben.
Ehhez első lépésben létrehozzuk a becslőobjektumot mint a KDEEstimator példányát adott kernellel, és sávszélességkorrekciós tényezővel. A nevezetes folytonos eloszlások közül a háromszög, béta, gamma, Weibull, normális és lognormális eloszlásokat fogjuk használni a vizsgálatokhoz. Ezek mintagenerátor függvényeit a szabványos könyvtár random moduljából importáljuk. Az egyes függvényneveket az eloszlást meghatározó paraméterekkel, valamint a megjelenítési intervallum alsó és felső szélével egy-egy tuple-ban fogjuk össze. Ezeket pedig egy listában gyűjtjük össze.
Mindezen előkészítés után egymás után vesszük a tuple konténereket, és kicsomagoljuk a tartalmukat. Az aktuális mintagenerátor függvénnyel előállítjuk a vizsgálandó mintasorozatot, amit átadunk a becslőobjektumnak. Végül, ezt használva megjelenítjük a becsült sűrűségfüggvényt és a becsült móduszt. A tesztkód és az eredmények alább láthatók. Érdemes több futtatást is végezni, hogy lássuk a becslések megfelelőségét.
|
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 |
# TESZT distr_names = {triangular: 'Háromszög eloszlás', betavariate: 'Béta-eloszlás', gammavariate: 'Gamma-eloszlás', weibullvariate: 'Weibull-eloszlás', normalvariate: 'Normális eloszlás', lognormvariate: 'Log-normális eloszlás'} estimator = KDEEstimator('cosine', 2) distributons_data = [(triangular, 0, 50, 35, 0, 50), # elméleti módusz = 35 (betavariate, 2, 5, -0.2, 1.2), # elméleti módusz = 0.2 (gammavariate, 4, 1, -2, 17), # elméleti módusz = 3 (weibullvariate, 1, 1.5, 0, 3), # elméleti módusz = 0.4807 (normalvariate, 0, 1, -4, 4), # elméleti módusz = 0 (lognormvariate, 0, 1, -2, 10) # elméleti módusz = 0.367 ] for distr_data in distributons_data: sampler, *params, plot_xmin, plot_xmax = distr_data estimator.samples = random_samples(sampler, *params, n=300) plot_pdf_with_mode(estimator, plot_xmin, plot_xmax, distr_names[sampler]) |
Megállapíthatjuk, hogy a becsült módusz elég jól közelíti az elméleti értéket. Egyetlen kivétel a lognormális eloszlás, ahol jelentősebb az eltérés.
Ebben a bejegyzésben a nem-paraméteres móduszbecslés egy lehetséges módszerét mutattuk be. A megvalósításhoz a szabványos könyvtár random és statistics moduljainak eszközkészletét használtuk. Ezekről a Python tudásépítés lépésről lépésre című e-könyv „Készétel fogyasztás – a szabványos könyvtár moduljainak használata” fejezet „Matematikai számítások támogatása” alfejezetén belül a „A véletlen használatba vétele” és „Statisztikai eszköztár” alcímek alatt van részletesen szó, példákkal segítve a megértést. A bemutatott tesztkódot egyszerűsítette a tuple konténerek elemeinek rugalmas kicsomagolási lehetősége. Erről az „Iterálható objektumok értékeinek ki- és becsomagolása” című részben lehet bővebben olvasni.