Keverék eloszlások

Az EM algoritmus különösen jól használható keverék eloszlások paramétereinek maximum likelihood becslésére, ezt a [5] cikk is elsők között említi. Ilyen esetben a számolások jelentősen leegyszrűsödnek, mivel a feladat kezelhető hiányos adat problémának. Most a kevert komponensek keverési aránya becslésének azt az esetét vizsgáljuk, amikor a komponensek sűrűsége ismert. Általánosabb, és nehezebb eset, amikor a komponensek sűrűségében is szerepelnek paraméterek, és ezeket együtt kell becsülni a keverési arányokkal. Ez az egyszerűbb eset sem valószerűtlen, mert gyakorlati példákban előfordul, hogy rendelkezésre állnak a komponensekből külön-külön vett minták is, ami által még a kevert eloszlás paramétereinek becslése előtt lehetőség nyílik a komponensek paramétereinek tetszőlegesen pontos becslésére.

Tehát tételezzük fel, hogy egy $\boldsymbol Y$ valószínűségi változó sűrűség függvényének $M$ komponensű kevert alakja van:

$$g_{\boldsymbol\theta }(y)=\sum_{i=1}^M \boldsymbol\theta _i g_i(y),$$

ahol $\boldsymbol\theta =(\boldsymbol\theta _1,\ldots,\boldsymbol\theta _{M-1})$ a paramétervektor. Ha nem akarunk feltételes szélsőértéket keresni, akkor célszerű egyel kevesebb $\boldsymbol\theta _i$ paramétert venni, a legutolsó paraméterre pedig teljesüljön

$$\boldsymbol\theta _M=1-\sum_{i=1}^{M-1}\boldsymbol\theta _i.$$

A $g_i$ komponens sűrűségek pedig teljesen meghatározottak.

Ez a modell például olyan gyakorlati eseteket foglal magába, ahol például egy populáció $G_1,\ldots,G_M$ $M$ különböző csoportra osztható bizonyos ismeretlen $\boldsymbol\theta _1,\ldots,\boldsymbol\theta _{M}$ arányban. és ahol az $i$ -edik csoportba való tartozás feltétele mellett az $\boldsymbol Y$ feltételes eloszlása $g_i.$ Például egy Do és McLachlan által 1984-ben tekintett gyakorlati problémában, a tekintett populáció $M$ különböző faj $G_1,\ldots,G_M$ csoportjába, részpopulációjába tartozó patkányok, melyek közül a baglyok minden csoportból valamilyen ismeretlen $\boldsymbol\theta _1,\ldots,\boldsymbol\theta _{M}$ valószínűséggel fogyasztanak el egy egyedet. Ez teljesülhet például, ha a $G_i$ csoport részaránya a teljes populációban $\boldsymbol\theta _i$ , és a baglyok minden egyedet egyforma valószínűséggel találnak meg, és fogyasztanak el, vagy lehet a $G_i$ részaránya tetszőleges, de a baglyok egy ide eső egyedet olyan valószínűséggel lelnek fel, hogy összességében $\boldsymbol\theta _i$ valószínűséggel jussanak hozzá egy egyedhez ebből a csoportból. A feladat a $\boldsymbol\theta _1,\ldots,\boldsymbol\theta _{M}$ paraméterek becslése bagolyköpetekből vett $n$ patkánykoponyán végzett mérések alapján, ahol egy ilyen koponyán mért adatokat tartalmazza az $\boldsymbol Y$ valószínűségi változó. A patkány hozzátartozik a baglyok étrendjéhez, és az emészthetetlen részeket bagolyköpet formájában öklendezik vissza.

Tehát csupán a koponyákon végzett mérések alapján nem tudjuk meghatározni a patkányok faját, viszont tudjuk, hogy az egyes fajokon belül a mért értékek milyen eloszlást követnek. Az eloszlások ismeretében csak annak feltételes valószínűségét tudjuk meghatározni, hogy a mért értékeket feltéve mennyi a valószínűsége annak, hogy a patkány egy adott fajba tartozott, ami a csoportba tartozás indikátor változójának a várható értéke. Ha tudnánk, hogy melyik patkány milyen fajba tartozik, egyszerűen vennénk a relatív gyakoriságokat, és ezek adnák a paraméterek maximum likelihood becslését. Most azonban csak a csoportba tartozás feltételes várható értékét tudjuk, ezért a gyakoriságok számolásánál ezzel a valószínűséggel súlyozzuk az adott csoportba tartozások számát. Egy egyed esetében tehát nem vagyunk biztosak benne, hogy hová tartozik és ezért nem egy egyest adunk ahhoz a csoporthoz tartozó mintából vett gyakorisághoz, amelyik csoporthoz tartozik, hanem ezt az egyest szétosztjuk az összes csoport között annak arányában, hogy mennyi az egyed odatartozásának valószínűsége feltéve a méréseink eredményét. Így a gyakoriság már nem egész lesz, hanem tetszőleges valós szám. Persze a feltételes valószínűség számolásához ismernünk kellene a paramétereket. Ezen az önmagába visszakanyarodó problámán úgy tudunk úrrá lenni, hogy veszünk egy kiindulási paramétert, amivel elvégezzük a számolást, majd az így kapott új paraméterekkel megismételjük, amíg a paraméterértékek állandósulnak. Ez az EM algoritmus. A konvergenciát az EM algoritmus tétel biztosítja.

A keverék eloszlásból vett független $n$ elemszámú megfigyelhető mintát jelölje

$${\bf y}=(y_1,\ldots,y_n),$$

ami a koponyákon mért értékek sorozata. ${\bf y}$ log likelihoodja

$$\sum_{j=1}^n \log g_{\boldsymbol\theta }(y_j)=\sum_{j=1}^n \log \sum_{i=1}^M \boldsymbol\theta _i g_i(y_j)$$

Ezt differenciálva a $\boldsymbol\theta _1,\ldots,\boldsymbol\theta _{M-1}$ paraméterek szerint, és egynlővé téve nullával kapjuk a

$$\sum_{j=1}^n\left({g_i(y_j)\over g_{\boldsymbol\theta }(y_j)}-{g_M(y_j)\over g_{\boldsymbol\theta }(y_j)}\right)=0$$

likelihood egyenletet kapjuk, ami általános esetben nehezen oldható meg a paraméterekre.

Ahhoz, hogy hiányos adat problémaként kezelhessük a feladatot vezessük be a hiányzó, illetve rejtett

$${\bf z}=(z_1,\ldots,z_n)$$

adatrendszert, ahol $z_i$ egy $M$ elemű indikátorváltozó és a $z_{ij}$ koordinátája 1, ha az $i$ -edik minta a $g_j$ sűrűségfüggvény szerinti eloszlásból származik. Példánkban ez azt jelenti, hogy az $i$ -edik koponya egy $j$ -edik fajba tartozó patkányé volt. De néha a $z_i$ változók mögött nem állnak valóságos csoportok, hanem a sűrűség ilyen felbontása csak egy elméleti segédeszköz ahhoz, hogy a problémát az EM algoritmus keretei közé helyezhessük.

Ha a $z_{ij}$ változók megfigyelhetőek lennének, akkor $\boldsymbol\theta _i$ becslése egyszerűen

$$\sum_{j=1}^n z_{ij}\over n$$

volna, ami a sűrűség $i$ -edik komponenséből származó minta részaránya. Amint az előző példában is, az EM algoritmus a megfigyelhető minta szerinti feltételes várható értéken keresztül veszi figyelembe a hiányzó adatokat. Legyen a teljes adatrendszer ennek megfelelően a hiányos adatrendszer kiegészítve a hiányzó adatrendszerrel:

$${\bf x}=({\bf y},{\bf z}).$$

A teljes adatrendszer log likelihoodja polinomiális alakú:

$$\sum_{i=1}^M\sum_{j=1}^n z_{ij}\log\boldsymbol\theta _i+C,$$

ahol

$$C=\sum_{i=1}^M\sum_{j=1}^n z_{ij}\log g_i(y_j),$$

ami nem függ a paraméterektől. Mivel ez a kifejezés lineáris $z_{ij}$ adatokban, ezért az EM algoritmusban szereplő feltételes várható értékhez egyszerűen csak kell venni $Z_ij$ feltételes várható értékét a megfigyelt adatok mellett. Ahol $z_{ij}$ a $Z_ij$ valószínűségi változó eloszlásából származik. Ki kell tehát számolni a feltételes várható értéket, amlely

$$\boldsymbol E_{\boldsymbol\theta }(Z_{ij}\vert{\bf y})=\boldsymbol P_{\boldsymbol\theta }(Z_{ij}=1\vert{\bf y}).$$

Láttuk korábban, hogy minden esetben, amikor a két változónak van sűrűségfüggvénye, akkor a feltételes várható értéket az együttes sűrűségfüggvény és a feltételben álló változó marginális sűrűségfüggvényének hányadosa szerint kell venni. Az együttes sűrűségfüggvény a $Z_{ij}=1, {\bf y}$ pontban $\boldsymbol\theta _i g_i(y_j)$ , az ${\bf y}$ pontban pedig a marginális sűrűség $g_{\boldsymbol\theta }(y_j).$ Ezért a feltételes várható érték

$$\boldsymbol E_{\boldsymbol\theta }(Z_{ij}\vert{\bf y})={\boldsymbol\theta _i g_i(y_j)\over g_{\boldsymbol\theta }(y_j)}$$

Ez tulajdonképpen a Bayes tétel sűrűségfüggvényekre vonatkozó alakja, mert $g_{\boldsymbol\theta }(y_j)=\sum_{i=1}^M \theta_i g_i(y_j)$ és ezért

$$\boldsymbol P_{\boldsymbol\theta }(Z_{ij}=1\vert{\bf y})={\boldsy... ..._j)\over{\boldsymbol\theta _1 g_1(y_j)+\cdots +\boldsymbol\theta _M g_M(y_j)}},$$

vagy feltételes sűrűségfüggvényekkel kifejezve

$$g(Z_{ij}=1\vert y_j)={g(y_j\vert Z_{ij}=1)\over{g(y_j\vert Z_{1j}=1)+\cdots +g(y_j\vert Z_{Mj}=1)}},$$

mert az előző egyenletben a bal oldal csak az $y_j$ értékétől függ, ugyanis ${\bf y}$ elemei egymástól függetlenek. A maximalizációs lépés a már felírt relatív gyakoriságok számolása lesz, csak éppen az adott csoporthoz tartozás feltételes várható értékeinek kell venni a relatív gyakoriságát. Összefoglalva, az előző paraméterek ismeretében a következő paramétereket a

$$\boldsymbol\theta _i^{k+1}={1\over n}\sum_{j=1}^n {\boldsymbol\th... ...over{\boldsymbol\theta _1^k g_1(y_j)+\cdots +\boldsymbol\theta _M^k g_M(y_j)}}$$

egyenlet definiálja.

Az EM algoritmus tétel biztosítja a konvergenciát a hiányos adatrendszer likelihoodjának egy stacionárius pontjához, de ebben az esetben be is bizonyítható, hogy az EM algoritmus tényleg a likelihood egyenlet megoldásához konvergál. Ha $\boldsymbol\theta ^*$ a likelihood egyenlet megoldása, akkor a paraméter helyébe a $\boldsymbol\theta ^*$ paramétert helyettesítve és a likelihood egyenlet mindkét oldalát beszorozva $\boldsymbol\theta _i^*$ paraméterrel, ahol az $i$ index $1$ és $M-1$ között változik, akkor a

$$\sum_{j=1}^n\left(\boldsymbol\theta _i^*{g_i(y_j)\over g_{\boldsy... ...}- \boldsymbol\theta _i^*{g_M(y_j)\over g_{\boldsymbol\theta ^*}(y_j)}\right)=0$$

egyenletet kapjuk. Az egyenlet $i=M$ esetén is fennál, ezért összegezhetjük minden $i$ -re, ami a

$$\sum_{j=1}^n\left({\sum_{i=1}^M\boldsymbol\theta _i^*g_i(y_j)\ove... ...}^M\boldsymbol\theta _i^*{g_M(y_j)\over g_{\boldsymbol\theta ^*}(y_j)}\right)=0$$

egyenletre vezet. Azaz

$$n=\sum_{j=1}^n{g_M(y_j)\over g_{\boldsymbol\theta ^*}(y_j)}$$

Ezt visszahelyettesítve az első egyenletbe, azt kapjuk, hogy

$$\sum_{j=1}^n\boldsymbol\theta _i^*{g_i(y_j)\over g_{\boldsymbol\theta ^*}(y_j)}= \boldsymbol\theta _i^*n,$$

ami azt jelenti, hogy $\boldsymbol\theta ^*$ megoldása az EM algoritmus által adott iterációnak, vagyis az EM algoritmus $\boldsymbol\theta ^*$ optimális paraméterhez konvergál.

További részletek és numerikus példa található a [4] könyvben.