Markov-modulált Poisson folyamatok statisztikai vizsgálata

Diplomamunka

Írta: Nagy Tamás

Matematikus szak

Témavezető:

Arató Miklós, egyetemi docens Valószínűségelméleti és Statisztika Tanszék

Eötvös Loránd Tudományegyetem, Természettudományi Kar

Eötvös Loránd Tudományegyetem Természettudományi Kar

Tartalomjegyzék

1. Bevezető 1 1.1. Bevezető . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1 1.2. A használt jelölések és rövidítések . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3

2. A Markov-modulált Poisson folyamat 4 2.1. A pontfolyamatok elméletének rövid áttekintése . . . . . . . . . . . . . . . . . . . . 4 2.2. Duplán sztochasztikus (Cox) folyamatok . . . . . . . . . . . . . . . . . . . . . . . . 14 2.3. Néhány szó mátrix értékű függvények differenciálásáról és a mátrix exponenciális

függvényről . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 21 2.4. Véges állapotterű, folytonos idejű Markov folyamatok . . . . . . . . . . . . . . . . 26 2.5. Kvázi-Markov folyamatok és Markov felújítási sorozatok . . . . . . . . . . . . . . . 28 2.6. A Markov-modulált Poisson folyamat, mint speciális Cox folyamat . . . . . . . . . 32

3. Az MMPP paramétereinek becslése az EM algoritmus segítségével 37 3.1. Az EM algoritmus . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 37 3.2. Egy EM algoritmus Markov-modulált Poisson folyamatokra . . . . . . . . . . . . . 39 3.3. Az algoritmus implementálása . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 43

4. Markov-modulált Poisson folyamat a gyakorlatban 50 4.1. Minta generálása adott paraméterű MMPP-ből . . . . . . . . . . . . . . . . . . . . 50 4.2. Egy biztosításmatematikai feladat . . . . . . . . . . . . . . . . . . . . . . . . . . . 52 4.3. Numerikus eredmények . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 58

5. Záró megjegyzések 62 5.1. Összegzés . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 62 5.2. További nyitott kérdések . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 62 5.3. Köszönetnyilvánítás . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 63

Függelék 63

A. A futási eredmények összefoglaló táblázatai 64 A.1. A valós paraméterek és az EM algoritmus által kapott becslés . . . . . . . . . . . . 64 A.2. A várható kárszám becslése . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 67 A.3. Kvantilisbecslések . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 69

Irodalomjegyzék 72

Ábrák jegyzéke

2.1. A t ↦ N_{t} trajektóriák egy-egy realizációja μ = 1 és μ = 1.5 esetén . . . . . . . . . 13 2.2. Egy SMP lehetséges realizációja . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 30 2.3. Egy MMPP lehetséges realizációja . . . . . . . . . . . . . . . . . . . . . . . . . . . 33

4.1. Napi kárszámok egy MMPP egyik realizációja esetén . . . . . . . . . . . . . . . . . 53

Táblázatok jegyzéke

1.1. Főbb jelölések és rövidítések . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3

3.1. A kiinduló becslés pontossága . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 46 3.2. A kiinduló becslés pontossága (folytatás) . . . . . . . . . . . . . . . . . . . . . . . . 46

4.1. A Rydén által vizsgált 2A és 2B eset . . . . . . . . . . . . . . . . . . . . . . . . . . 58 4.2. Szimulált kárszámok a 2A és 2B esetben . . . . . . . . . . . . . . . . . . . . . . . . 58 4.3. Rydén becslése és a mostani becslés . . . . . . . . . . . . . . . . . . . . . . . . . . 59 4.4. Az 𝔼(ζ_{k+1}) előrejelzések a három különböző modell esetén . . . . . . . . . . . . . . 59 4.5. A 0.95-os megbízhatósági szinthez tartozó kvantilisbecslések . . . . . . . . . . . . . 59 4.6. A 0.99-os megbízhatósági szinthez tartozó kvantilisbecslések . . . . . . . . . . . . . 60 4.7. Az MMPP becslésből számított kvantilisek eltéréseinek a gyakorisága . . . . . . . . 60 4.8. A kevert Poisson becslésből számított kvantilisek eltéréseinek a gyakorisága …. 60 4.9. MMPP eltérések gyakorisága 0.99-es mb. szinten . . . . . . . . . . . . . . . . . . . 61 4.10. Kevert Poisson gyakoriságok 0.99-es mb. szinten . . . . . . . . . . . . . . . . . . . 61

A.1. A becsülendő paraméterek és a kiinduló becslések . . . . . . . . . . . . . . . . . . . 64 A.2. A becsülendő paraméterek és a kiinduló becslések (folytatás) . . . . . . . . . . . . 65 A.3. A becsülendő paraméterek és az EM algoritmus által kapott becslés . . . . . . . . 66 A.4. A becsülendő paraméterek és az EM algoritmus által kapott becslés (folytatás) . . 67 A.5. Az 𝔼(ζ_{k+1}) előrejelzések valós és becsült értékei . . . . . . . . . . . . . . . . . . . . 67 A.6. Az 𝔼(ζ_{k+1}) előrejelzések valós és becsült értékei (folytatás) . . . . . . . . . . . . . 68 A.7. A 95%-os megbízhatósági szinthez tartozó kvantilisbecslések . . . . . . . . . . . . . 69 A.8. A 95%-os megbízhatósági szinthez tartozó kvantilisbecslések (folytatás) . . . . . . 70 A.9. A 99%-os megbízhatósági szinthez tartozó kvantilisbecslések . . . . . . . . . . . . . 70 A.10.A 99%-os megbízhatósági szinthez tartozó kvantilisbecslések (folytatás) . . . . . . 71

1. fejezet

Bevezető

1.1. Bevezető

Az olyan pontfolyamatoknak, melyek intenzitása az időben véletlenszerűen változik, számos gyakorlati alkalmazása van, a sorbanállási és kiszolgálási rendszerektől kezdve az orvostudományon és pénzügyi modelleken át a számítógépes hálózatokig bezárólag. A Markov-modulált Poisson folyamat (MMPP) - szemléletesen kifejezve - egy olyan általánosított Poisson folyamat, melynek intenzitása egy adott időpillanatban egy folytonos időparaméterű és véges állapotterű Markov folyamat aktuális állapota által meghatározott, s mint ilyen, a véletlentől függő érték.

Az MMPP első használata a kiszolgálási rendszerek elméletében Naor és Yechiali [15] cikkéhez fűződik, melyet aztán a Neuts által írt [16] publikáció követett. Azóta több kiváló értekezés is született a témakörben. Az egyik legátfogóbb - és egyben talán a leggyakrabban hivatkozott - mű a W. Fisher és K. Meier-Hellstern által 1993-ban közölt cikk [6]. Szintén alapvetőnek számítanak a nem sokkal később napvilágot látó [22] és [23] publikációk.

Az MMPP gyakorlati modellezés során való alkalmazhatósága miatt természetes módon merül fel a paraméterek becslésének kérdése. Ez egyrészt a mögöttes Markov folyamat kiinduló valószínűségi vektorának és infinitezimális generátormátrixának, másrészt a különböző intenzitások becsléséből áll. A jelenleg ismert módszerek azonban nem koncentrálnak a kiinduló valószínűségi vektor becslésére, hanem azzal a feltételezéssel élnek, hogy az nem egy független paraméter, hanem jól definiált módon kapható a többi paraméterből, pl. ez az eset áll fenn akkor, ha az említett kiinduló valószínűségi vektor éppen a mögöttes Markov folyamat stacionárius eloszlása. Az MMPP rendjét szintén ismertnek feltételezzük. A becslés során alkalmazható technikákról jó áttekintés található [22] -ben, melyben T. Rydén egyebek mellett bemutatott egy általánosított EM (GEM) algoritmust a maximum-likelihood becslés kiszámítására a 2-állapotú MMPP estében. Ennek az algoritmusnak a hátránya az, hogy a maximalizációs - ún. M - lépés nem explicit. 2 évre rá azonban Rydén módosította korábbi algoritmusát (lásd [23]) úgy, hogy az már egy explicit M-lépést és egy valamivel bonyolultabb, ám még kezelhető E-lépést eredményezett. További fejlesztésként megszüntette a csupán 2-állapotú esetre való szorítkozást, az új algoritmus már tetszőleges rendű MMPP-re is alkalmazhatóvá vált.

A [23]-ben közölt algoritmusnak azonban továbbra is maradt két gyenge pontja. Egyrészt, az E-lépés kiszámításában szerepet játszó ún. előre tekintő, illetve hátráló rekurziók kalkulációja során olyan értékek merültek fel, melyek meghaladják egy „átlagos” személyi számítógép számítási kapacitását. Rydén az említett cikkben felhívta erre a figyelmet, megjegyezve, hogy a probléma a rekurziók alkalmas átskálázásával elkerülhető, azonban nem közölte ennek konkrét kivitelezését. Másrészt - szintén az E-lépés kapcsán - mátrix exponenciális függvények bizonyos integráljainak kiértékelése merült fel, melyhez Rydén akkor mátrix-diagonalizációt használt, ami jelentősen lelasítja az algoritmust. Ezen problémás részek kezelésére csupán a közelmúltban született megoldás, W.J.J. Roberts, Y. Ephraim és E.Dieguez 2006-ban megjelent publikációjukban mindkét kérdés részletes tárgyalásra került, a szerzők bemutattak egy megfelelő skálázási eljárást az előre tekintő, illetve hátráló rekurziók normalizálására, illetve hatékony és szép megoldást adtak a másik problémára is, azaz a mátrix exponenciális függvény integráljainak kiszámítására C. Van Loan egy korábbi [29] eredménye alapján.

Jelen dolgozat egyik tárgya T. Rydén és W.J.J. Roberts et. al. eredményeinek áttekintése és ismertetése, illetve a Rydén féle algoritmus újbóli implementálásának részletes leírása a [21]-ben közölt eredmények ismeretében. Tekintünk továbbá egy olyan biztosításmatematikai problémát, ahol egy biztosítótársaság valamely adott, homogén kockázatú üzletágában felmerülő kárigények számának időbeli fejlődésének leírására MMPP modellt használ és a rendelkezésére álló minta alapján készített becslés ismeretében kíván előrejelzéseket adni egy soron következő időperiódusban bekövetkező károk várható számára, illetve az időperiódusra eső károk számának az eloszlásának a különböző elvárt megbízhatósági szintekhez tartozó kvantilis értékekre.

A dolgozat további felépítése a következők szerint alakul. A 2. fejezet a Markov-modulált Poisson folyamatok kellően megalapozott bevezetésére törekszik. A fő definíciók a 6. alfejezetben vannak, az ezt megelőző alfejezetek ezek előkészítését szolgálják. A 3. fejezetben kerül bemutatásra az EM algorimus, mint általános iteratív paraméterbecslési eljárás nem megfigyelhető adattal rendelkező feladatok esetén, majd speciálisan a már korábban említett Rydén-féle EM algoritmus MMPP-k becslésére a W.J.J. Roberts et. al. által eszközölt javításokkal. A 4. fejezetben egy biztosításmatematikában felmerülő, MMPP modellt használó feladatot vizsgálunk a 3. fejezetben megismert algoritmus segítségével, s ennek kapcsán konkrét numerikus eredményeket is közlünk. Az 5. fejezetben összefoglaljuk a dolgozat eredményét, illetve megemlítjük az esetleges további kutatási lehetőségeket és nyitva maradt kérdéseket.

1.2. A használt jelölések és rövidítések

Jelölés vagy rövidítés Magyarázat
F (∞)  lim_{x↦∞} F (x)
f.a.e.  független, azonos eloszlású (valószínűségi változók)
SMP  kvázi-Markov folyamat
MRS  Markov felújítási sorozat
MMPP  Markov-modulált Poisson folyamat
m.b.  majdnem biztosan
M_{ij}  az M mátrix i-edik sorának j-edik eleme
(m_{ij})  az a mátrix, melynek i. sorának j. oszlopában m_{ij} áll
M_{i}  Az M mátrix i-edik sora
M^{(j)}  Az M mátrix j-edik oszlopa
M_{i·}  Az M mátrix i-edik sorában lévő elemek összege
I  identitás mátrix
0  csupa nulla mátrix
F (t+)  lim_{s↦0+} F (t + s)
1  csupa egy oszlopvektor
a ∨ b  max(a, b)
M(S, B)  (S, B)-pontmértékek tere
\overline{ℕ}_{0}  ℕ_{0} ∪ \{∞\}
ν|A  az A σ-algebrán értelmezett mérték, vagy a ν mérték A-ra vett

megszorítása

\frac{dν}{dΛ}  a ν mérték Radon-Nykodim deriváltja a Λ mértékre nézve
M_{R}(S, B)  véletlen mértékek tere
=^{d}  eloszlásbeli egyenlőség
Q_{X}, F_{X}, f_{X}  az X valószínűségi változó eloszlása, eloszlásfüggvénye, illetve

sűrűségfüggvénye

k.v.v.  kiindulási valószínűségi vektor
diag(x), ahol x vektor  Az x vektorból képzett diagonális mátrix
diag(M), ahol M mátrix Az M mátrix diagonálisában álló elemekből képzett vektor

1.1. táblázat. Főbb jelölések és rövidítések

A Markov-modulált Poisson folyamat

Az 1. alfejezet egy általános bevezetést nyújt a pontfolyamatok elméletéről. Itt az alapvető definíciók R.-D. Reiss [20] munkájában leírt felépítésből származnak, az alfejezet túlnyomó része azonban saját munka, különös tekintettel a pontfolyamatok és számlálófolyamatuk kapcsolatának megvilágítására. A 2. alfejezet az ún. duplán sztochasztikus folyamatok általános elméletét ismerteti, azon okból kifolyólag, hogy egy Markov-modulált Poisson folyamat egy speciális duplán sztochasztikus folyamat. Itt kétféle definíciót ismertetünk, a 2.2.2 definíció saját, míg a 2.2.4 definíció Reiss-től származik. Az alfejezet további része - mely megvizsgálja például a két definíció kapcsolatát - mind saját eredmény. Mátrixok exponenciálisával kapcsolatos összefüggések folyton előkerülnek majd az egész dolgozat során, ezért az ehhez kapcsolódó legfontosabb információkat a 3. alfejezetben gyűjtöttük össze. Itt a 2.3.8 lemma saját eredmény, mely C. Van Loan egy korábbi lemmájának (lásd [29] , 2. fejezet, 1. tétel) analógiájára készült. Ez a lemma csak egy későbbi fejezetben kerül felhasználásra. A 4. alfejezet a folytonos idejű, véges állapotú Markov folyamatok elméletéből ismertet néhány szükséges fogalmat. Az 5. alfejezet a kvázi-Markov folyamatok és Markov felújítási sorozatok alapvető definícióit ismerteti R. Pyke [18] és [19] felépítése alapján. Az alfejezet szükségességét az indokolja, hogy - amint azt később látni fogjuk - az MMPP-k tárgyalása beilleszthető az így kapott általános elméletbe. A 6. alfejezet definálja a Markov-modulált Poisson folyamat fogalmát. Az alfejezetben saját eredményként levezetünk egy összefüggést a ℙ(N_{T +t} − N_{T} = n) mennyiségekre (N_{t} definícióját lásd az alfejezetben) a W. Fisher és K. S. Meier-Hellstern által felírt Chapman-Kolmogorov egyenletek segítségével.

2.1. A pontfolyamatok elméletének rövid áttekintése

Jelen alfejezet a pontfolyamatok elméletének és a hozzájuk kapcsolódó fogalmak rövid áttekintését tűzi ki célul. Először ismertetünk egy általános, mértékelméleti felépítést, mely R.-D. Reiss [20] munkáját követi, majd ezeket ültetjük át a gyakorlati alkalmazásoknak megfelelő valós esetre. Nem törekszünk a pontfolyamatok szép, ám igencsak szerteágazó elméletéről teljes képet adni, csak azon elemi összetevők ismertetésére szorítkozunk, amelyek a további fejezetek megértéséhez szükségesek. Tegyük fel, hogy adott egy S alaphalmaz, melyet a továbbiakban állapottérnek nevezünk.

2.1.1. Definíció. Egy C ⊂ S legfeljebb megszámlálható számosságú halmazt konfigurációnak hívunk.

2.1.2. Definíció. Legyen C = {x_{i} : i ∈ I} ⊂ S egy adott konfiguráció az állapottéren. Ekkor a μ_{C}(B) = |C ∩ B| képlettel meghatározott diszkrét mértéket a B σ-algebrán értelmezett, és a C konfigurációhoz tartozó pontmértéknek nevezzük.

Könnyen láthatóan μ_{C} = ∑_{i∈I} δ_{x_{i}}, ahol δ_{x} a Dirac-delta.

2.1.3. Definíció. Legyen (S, B) egy mérhető tér, és legyen μ egy B-n értelmezett halmazfügvény. Ha létezik olyan C ⊂ S konfiguráció, hogy μ előáll μ(B) = |C ∩ B| alakban, akkor azt mondjuk, hogy a μ halmazfüggvény (S, B)-pontmérték.

2.1.1. Megjegyzés. Az elnevezés jogos, hiszen μ(∅) = |C ∩ ∅| = 0, μ(B) ≥ 0 és páronként diszjunkt B_{1,}B_{2,} … ∈ B esetén

μ ( ⋃_{i=1}^{∞}B_{i}) = |C ∩ ( ⋃_{i=1}^{∞}B_{i})| = | ⋃_{i=1}^{∞}(C ∩ B_{i})| = ∑_{i=1}^{∞} |C ∩ B_{i}| = ∑_{i=1}^{∞} μ(B_{i}),

tehát μ valóban mérték.

Az S-en értelmezett konfigurációk terét jelölje C(S) = {C ⊂ S : C konfiguráció}, az (S, B)- pontmértékek terét pedig M ≡ M(S, B). Legyen κ az a C(S) ↦ M(S, B) leképezés, ami egy adott C konfigurációhoz hozzárendeli a μ_{C} = |C ∩ ·| pontmértéket. Ekkor definíció szerint μ szürjektív, viszont nem feltétlen injektív. Ha azonban a B σ-algebra kellően gazdag struktúrájú - például, ha B tartalmazza az egypontú halmazokat -, akkor az injektivitás igaz, tehát ezekben az esetekben κ egy bijekció.

A következőkben az lesz a célunk, hogy alkalmas módon definiáljunk egy σ-algebrát az M téren. Ehhez tekintsük minden B ∈ B -re a

π_{B} : M ↦ \overline{ℕ}_{0}
μ ↦ μ(B)

leképezést - ahol \overline{ℕ}_{0} = ℕ_{0} ∪ \{∞\} -, majd ezek után legyen M a legszűkebb olyan σ-algebra, melyre nézve az összes π_{B} leképezés (B ∈ B) - mint (M, M) ↦ (\overline{ℕ}_{0}, P (\overline{ℕ}_{0})) leképezés mérhető, azaz

M ≡ M(S, B) = σ\{π_{B} : B ∈ B\} = σ\{π_{B}^{−1}(C) : B ∈ B, C ⊂ \overline{ℕ}_{0}\}

A továbbiakban - ha külön nem jelezzük - \overline{ℕ}_{0}-ba képező függvények mérhetőségét mindig P (\overline{ℕ}_{0})-ra nézve értjük.

Most már rátérhetünk a pontfolyamat fogalmának definíciójára. Ennek érdekében rögzítsünk egy (Ω, A, ℙ) valószínűségi mezőt.

2.1.4. Definíció. Egy N : Ω ↦ M(S, B) leképezést (S, B) -pontfolyamatnak nevezünk, ha (A, M)- mérhető.

N tehát egy véletlen mérték B-n. Mielőtt megvizsgálnánk egy konkrét példát, bebizonyítunk egy karakterizációs tételt, mely segítségével könyebben ellenőrizhető lesz majd a definícióban szereplő mérhetőségi tulajdonság. Először is, minden B ∈ B -re definiálható az

\begin{matrix}N(B) : Ω & ↦ & \overline{ℕ}_{0} \\ ω & ↦ & N(ω)(B) \end{matrix}

leképezés. N(B) tehát nem más, mint a B-be eső pontok véletlen száma. N(S) pedig az N pontfolyamat mintaelemszáma. Könnyen láthatóan igaz a következő előállítás:

N(B) = π_{B} ∘ N.  (2.1)

Mivel N minden ω ∈ Ω-ra egy mérték, ezért fennállnak az alábbi tulajdonságok:

(1) N ≥ 0

(2) Ha A_{1}, A_{2}, … ∈ B páronként diszjunktak, akkor N( ⋃_{k=1}^{∞}A_{k}) = ∑_{k=1}^{∞} N(A_{k})

(3) Ha B ⊂ C , N(B) < ∞ , akkor N(C\B) = N(C) − N(B)

(4) Ha B ⊂ C , akkor N(B) ≤ N(C)

A továbbiakban N(B)-re, mint a N leképezés B-re vonatkozó egydimenziós marginálisára fogunk hivatkozni.

2.1.2. Megjegyzés. Természetesen a marginálisok definiálásához nem szükséges feltenni az (A, M)- mérhetőséget, illetve az sem szükséges, hogy N minden ω-ra véges értékű mértéket vegyen fel.

Ezek után az alábbi egyszerű, ám fundamentális jellegű kritérium fogalmazható meg:

2.1.1. Állítás. Legyen adott az N : Ω ↦ M(S, B) leképezés. Ekkor a következő két állítás ekvivalens:

(i) N pontfolyamat.

(ii) N(B) : Ω ↦ \overline{ℕ}_{0} (A, P (\overline{ℕ}_{0}))-mérhető minden B ∈ B -re.

Bizonyítás. (i)⇒(ii) M konstrukciója miatt triviális, hiszen az (2.1) előállítás értelmében N(B) két mérhető leképezés kompozíciója.

(ii)⇒(i) A \{π_{B}^{−1}(C) : B ∈ B, C ⊂ \overline{ℕ}_{0}\} halmazrendszer M egy generátorrendszere, ezért N mérhetőségéhez elég belátni, hogy minden B ∈ B és C ⊂ \overline{ℕ}_{0}-ra N^{−1}(π_{B}^{−1}(C)) = (π_{B} ∘ N)^{−1}(C) = N(B)^{−1}(C) ∈ A, ami viszont N(B) mérhetőségi tulajdonsága miatt nyilvánvaló.

2.1.1. Példa. Legyenek adottak az X_{1,}X_{2}, …, X_{k} : Ω ↦ (S, B) valószínűségi változók, és N legyen egyenlő ∑_{i=1}^{k} δ_{X_{i}}-vel. Így pontfolyamatot kapunk, hiszen ha B∈ B, akkor tekintve, hogy N a δ_{X_{i}}, i = 1, 2, … függvények összege, elég belátni a δ_{X_{i}} függvények mérhetőségét, ami viszont nyilvánvaló abból a tényből kifolyólag, hogy mérhető halmaz indikátora is mérhető:

δ_{X_{i}}(B) = χ(X_{i} ∈ B).

A következőkben azt tűzzük ki célul, hogy a 2.1.1 állítás második pontjának feltevését gyengítsük oly módon, hogy azt a B halmazrendszer elemei helyett elegendő legyen egy kevésbé gazdag struktúrájú halmazrendszer elemeire leellenőrizni. Ehhez először igazolunk egy lemmát, amit a későbbiekben is többször felhasználunk.

2.1.1. Lemma. Legyen adott egy F ⊂ P (X) halmazrendszer, mely rendelkezik az alábbi három tulajdonsággal:

1. X előállítható F-beli halmazok megszámlálható uniójaként,

2. F különbségzárt és

3. zárt a megszámlálható diszjunkt unióra.

Ekkor F σ-algebra.

Bizonyítás. Elég belátni, hogy F σ-uniózárt és komplementerzárt. Az előbbihez tekintsük a tetszőleges A_{i} ∈ F, i = 1, 2, … halmazokat. Ezeket diszjunktizáljuk, azaz megadunk egy olyan páronként diszjunkt halmazokból álló B_{i} ∈ F, i = 1, 2, … halmazrendszert, amire ∪_{i=1}^{∞}B_{i} = ∪_{i=1}^{∞}A_{i}. Ezt a következő módon tehetjük meg: B_{1} := A_{1} ∈ F, majd i ≥ 2-re legyen B_{i} := A_{i}\ ∪_{k=1}^{i−1} B_{k} ∈ F, ez nyilván teljesíti a kívánt tulajdonságot. Tehát az A_{i} halmazok unióját előállítottuk a páronként diszjunkt B_{i} halmazok uniójaként, ami viszont F-beli, hiszen F zárt a megszámlálható diszjunkt unióra. Tehát F valóban σ-uniózárt. A lemma 1. feltétele és a most bizonyított uniózártság miatt X ∈ F. Ezután, ha A ∈ F, akkor a különbségzártság miatt A^{C} = X\A ∈ F, tehát F valóban különbségzárt. Ezt kellett bizonyítanunk.

Ezek után már belátható az alábbi állítás.

2.1.2. Állítás. Legyen adott az N : Ω ↦ M(S, B) leképezés. Legyen a G ⊂ P(S) halmazrendszer a B σ-algebra egy generátorrendszere és tegyük fel, hogy S előállítható G-beli halmazok megszámlálható uniójaként. Ekkor az alábbi három állítás ekvivalens.

(i) N pontfolyamat.

(ii) N(B) : Ω ↦ \overline{ℕ}_{0} (A, P (\overline{ℕ}_{0}))-mérhető minden B ∈ B -re.

(iii) N(B) : Ω ↦ \overline{ℕ}_{0} (A, P (\overline{ℕ}_{0}))-mérhető minden B ∈ G-re.

Bizonyítás. (i)⇒(ii)⇒(iii) nyilvánvaló, be kell látnunk még a (iii)⇒(i) irányt. Ehhez legyen

E = σ{π_{B} : B ∈ G}.

{π_{B} : B ∈ G} ⊂ {π_{B} : B ∈ B}-ből következik, hogy E ⊂ σ{π_{B} : B ∈ B}. Legyen F = {B ∈ B : π_{B} E-mérhető}. Ekkor nyilvánvalóan G ⊂ F ⊂ B. Másrészt fenállnak az alábbiak.

1. S előállítható megszámlálható sok G-beli - s így egyúttal F-beli - halmaz uniójaként.

2. F különbségzárt, hiszen A, B ∈ F, A ⊂ B esetén

π_{B\A}(ν) = ν(B\A) = (ν(B) − ν(A))χ(ν(A) < ∞) = (π_{B}(ν) − π_{A}(ν))χ(π_{A}(ν) < ∞),

ami viszont mérhető, hiszen ν ↦ π_{B}(ν), ν ↦ π_{A}(ν), illetve ν ↦ χ(π_{A}(ν) < ∞) egyaránt mérhető.

3. F zárt a megszámlálható diszjunkt unióra, hiszen ha B_{i} ∈ F, i = 1, 2, … páronként diszjunkt

halmazok, akkor

π_{∪B_{i}}(ν) = ν(∪_{i}B_{i}) = ∑_{i=1}^{∞} ν(B_{i}) = ∑_{i=1}^{∞} π_{B}(ν),

azaz π_{∪B_{i}} = ∑_{i=1}^{∞} π_{B}, megszámlálható sok mérhető függvény összege pedig megint csak mérhető.

Így a 2.1.2 lemma segítségével azt kapjuk, hogy F σ-algebra, ekkor azonban csak az F ≡ B eset lehetséges. Ez éppen azt jelenti, hogy π_{B} E-mérhető minden B ∈ B-re, vagyis σ{π_{B} : B ∈ B} ⊂ E. Ekkor tehát

σ{π_{B} : B ∈ G} = σ{π_{B} : B ∈ B},

vagy másképpen írva

σ\{π_{B}^{−1}(C) : B ∈ G, C ⊂ \overline{ℕ}_{0}\} = σ\{π_{B}^{−1}(C) : B ∈ B, C ⊂ \overline{ℕ}_{0}\},

azaz a \{π_{B}^{−1}(C) : B ∈ G, C ⊂ \overline{ℕ}_{0}\} halmazrendszer is M egy generátorrendszere. Így elég belátnunk, hogy minden B ∈ G és C ⊂ \overline{ℕ}_{0}-ra N^{−1}(π_{B}^{−1}(C)) = (π_{B} ∘ N)^{−1}(C) = N(B)^{−1}(C) ∈ A, ami viszont N(B) mérhetőségi tulajdonsága miatt nyilvánvaló.

A következő definíciók szintén alapvetőek a pontfolyamatok elméletében:

2.1.5. Definíció. Az N (S, B)-pontfolyamat intenzitásmértéke alatt a ν(B) = 𝔼N(B), (B ∈ B) összefüggéssel definiált mértéket értjük.

ν(B) tehát a B halmazba eső pontok várható számát adja meg. Könnyű ellenőrízni, hogy ν valóban mérték B-n.

2.1.2. Példa. A 2.1.1 példában szereplő folyamat intenzitásmértékére

ν(B) = 𝔼N(B) = ∑_{i=1}^{k} 𝔼(δ_{X_{i}}(B)) = ∑_{i=1}^{k} ℙ(X_{i} ∈ B) = ∑_{i=1}^{k} ℙ(X_{i} ∈ B) = ∑_{i=1}^{k} Q_{X_{i}}(B)

adódik, azaz ν = ∑_{i=1}^{k} Q_{X_{i}}.

2.1.3. Példa. Legyen ν|B egy egészen σ-véges mérték. Azt mondjuk, hogy az N : Ω ↦ M(S, B) pontfolyamat Poisson folyamat, ha az alábbi két tulajdonság teljesül:

1. N(B) Poisson eloszlású ν(B) paraméterrel ∀B ∈ B, ν(B) < ∞ -re, és

2. N(B_{1}), N(B_{2}), …, N(B_{k}) függetlenek minden k ∈ ℕ_{0} és diszjunkt B_{1,}B_{2,} …, B_{k} ∈ B , ν(B_{i}) < ∞ esetén.

Ekkor az intenzitásmérték éppen 𝔼N(B) = ν(B), amennyiben ν(B) < ∞. ν egészen σ-véges tulajdonságát kihasználva egyszerűen megmutatható, hogy ekkor a Poisson folyamat előállítható megszámlálható sok páronkét független, véges intenzitású Poisson folyamat összegeként.

Koncentráljunk most arra az esetre, ami a tipikusan felmerülő szituáció a gyakorlati alkalmazásokban, vagyis amikor az S állapottér a pozitív valós félegyenes - azaz S = ℝ^{+} = (0, ∞) - és B az S Borel-halmazainak a σ-algebrája: B = B(ℝ^{+}) ≡ B_{0}. Mivel ekkor B_{0} tartalmazza az egypontú halmazokat, ezért korábbi észrevételünk alapján ekkor a κ : C(ℝ^{+}) ↦ M(ℝ^{+}, B_{0}) leképezés bijektív, azaz egy pontmérték és az őt meghatározó konfiguráció kölcsönösen egyértelműen meghatározzák egymást. Tekintsük most a következő definíciókat.

2.1.6. Definíció. Azt mondjuk, hogy a μ ∈ M(ℝ^{+}, B_{0}) pontmérték torlódó, ha a hozzá tartozó konfigurációnak létezik ∞-től különböző torlódási pontja. Egyéb esetben a μ pontmértéket nem torlódónak mondjuk.

Jelölje M_{0}(ℝ^{+}, B_{0}) a nem torlódó pontmértékek terét, azaz

M_{0} := {μ ∈ M : κ^{−1}(μ)-nek nem létezik ∞-től különböző torlódási pontja}.

Tehát M_{0}(ℝ^{+}, B_{0}) azon pontmértékek terét jelöli, amelyekre igaz az, hogy a hozzájuk tartozó konfigurációnak nem létezik ∞-től különböző torlódási pontja. Nyilvánvaló, hogy μ ∈ M_{0} akkor és csak akkor, ha μ korlátos intervallumokhoz véges értéket rendel. Kevés utánagondolással belátható, hogy fennáll az alábbi karakterizáció is.

2.1.3. Állítás. Legyen adott a μ pontmérték és tekintsük ℝ^{+} felbontását megszámlálható sok korlátos intervallumra: ℝ^{+} = ⋃_{i=1}^{∞}B_{i}. Ekkor

μ ∈ M_{0} ⇔ μ(B_{i}) < ∞ ∀B_{i}, i = 1, 2, … esetén.

Hogy M_{0} is mérhető tér legyen, definiáljunk rajta egy M_{0} σ-algebrát, mégpedig az M(ℝ^{+}, B_{0}) σ-algebra M_{0}(ℝ^{+}, B_{0})-ra vett nyomaként, azaz:

M_{0} := M|M_{0} = {H ∩ M_{0} : H ∈ M}.

Persze az is világos, hogy M_{0} = σ{π_{B}|M_{0} : B ∈ B_{0}}. Ekkor tekintsük a következő definíciót:

2.1.7. Definíció. Azt mondjuk, hogy az N : Ω ↦ M(ℝ^{+}, B_{0}) pontfolyamat nem torlódó, ha N ∈ M_{0}(ℝ^{+}, B_{0}) majdnem biztosan.

A 2.1.3 állítás segítségével belátható a következő állítás.

2.1.4. Állítás. Legyen adott az ℝ^{+} állapottér ℝ^{+} = ⋃_{i=1}^{∞}B_{i} felbontása megszámlálható sok korlátos intervallumra és tekintsük az N : Ω ↦ M(ℝ^{+}, B_{0}) pontfolyamatot. Ekkor

N nem torlódó ⇔ ℙ(N(B_{i}) < ∞) = 1 ∀B_{i}, i = 1, 2, … esetén.

Bizonyítás. N nem torlódó ⇔ ℙ(N ∈ M_{0}) = 1 ⇔ ℙ (⋂_{i=1}^{∞}(N(B_{i}) < ∞)) = 1 ⇔ ℙ (⋃_{i=1}^{∞}(N(B_{i}) = ∞)) = 0 ⇔ ℙ(N(B_{i}) = ∞) = 0, i = 1, 2, … ⇔ ℙ(N(B_{i}) < ∞) = 1, i = 1, 2, … esetén.

Az alábbiakban következő valós esetre vonatkozó definíciókban, illetve állításokban ahol szükséges, mindig feltesszük, hogy a pontfolyamatunk nem torlódó. Ez nem jelent számunkra nagy megszorítást, hiszen a későbbi vizsgálódásaink tárgyát képező pontfolyamatok úgyis nem torlódóak lesznek majd.

Tekintsük most az alábbi függvényt:

\begin{aligned}μ_{(·)} : ℝ^{+} & & ↦ & \overline{ℕ}_{0} \\ & t & ↦ & μ((0, t]) \end{aligned}

Ekkor nyilvánvaló, hogy μ_{(·)} jobbról folytonos módon kiterjeszthető ℝ_{0}^{+}-ra a μ_{0} = lim_{t↦0+} μ_{t} értelmezéssel. Egyben az is világos, hogy μ_{0} csak ∞ vagy 0 lehet attól függően, hogy a μ pontmérték torlódó-e vagy sem. Vezessünk is be erre egy definíciót.

2.1.8. Definíció. A μ ∈ M(ℝ_{0}^{+}, B_{0}) pontmérték által meghatározott számlálófüggvény alatt azt a μ_{(·)} : ℝ_{0}^{+} ↦ \overline{ℕ}_{0} függvényt értjük, amelyre

μ_{t} = \{ \begin{matrix}μ((0, t]) & , ha & t > 0 \\ lim_{t↦0+} μ_{t} & , ha & t = 0 \end{matrix} .

2.1.5. Állítás. Tegyük fel, hogy a μ ∈ M_{0}(ℝ_{0}^{+}, B_{0}). Ekkor μ számlálófüggvénye rendelkezik az alábbi tulajdonságokkal.

1. μ_{t} ∈ ℕ_{0}.

2. μ_{0} = 0.

3. μ_{(·)} monoton növő, jobbról folytonos, balról határértékkel rendelkező függvény 1 nagyságú ugrásokkal.

Bizonyítás.

1. Korlátos intervallum mértéke véges.

2. Világos, hiszen μ ∈ M_{0}, így μ nem torlódik a nullában, így létezik s > 0, hogy μ_{t} = μ((0, t]) = 0 minden 0 < t < s esetén.

3. t ≥ s esetén

μ_{t} − μ_{s} = μ((0, t]) − μ((0, s]) = μ((s, t]) ≥ 0,

hiszen μ korlátos intervallumokon véges értéket vesz fel, így igaz a szubtraktív tulajdonság. Szintén ez okból kifolyólag μ monoton folytonos is, tehát lim_{s↘0} μ_{t+s} = lim_{s↘0} μ((0, t+s]) = μ((0, t]) = μ_{t}, tehát μ_{(·)} jobbról folytonos. Az x > 0 -beli baloldali határérték vizsgálatához pedig lim_{t↘0} μ_{x} − μ_{x−t} = lim_{t↘0} μ_{x} − μ_{x−t} = lim_{t↘0} μ((x − t, x]) = μ({x}) =

\{ \begin{matrix}0 & ,ha & x ∉ C \\ 1 & ,ha & x ∈ C \end{matrix} . Ezzel beláttuk a baloldali határérték létezését, illetve azt is, hogy az

ugrások csak 1 nagyságúak lehetnek.

Nyilván a 2. és 3. tulajdonság együttesen implikálja az elsőt. Ezzel az iménti lemma úgy foglalható össze, hogy egy μ ∈ M_{0} pontmérték által meghatározott számlálófüggvény egy olyan nullában eltűnő monoton növő függvény, amely jobbról folytonos, balról határértékkel rendelkezik és ugrásai 1 nagyságúak. Egy ilyen tulajdonságú f : ℝ_{0}^{+} ↦ ℕ_{0} függvényt a továbbiakban korlátos számlálófüggvénynek nevezünk. A korlátos számlálófüggvények terére bevezetjük a BCF jelölést. Megmutatjuk, hogy minden korlátos számlálófüggvény nem torlódó pontmérték által meghatározott és azt, hogy különböző nem torlódó pontmértékekhez különböző számlálófüggvény tartozik. Jelölje α azt a leképezést, ami egy μ ∈ M pontmértékhez hozzárendeli az általa meghatározott számlálófüggvényt, akkor az előzőek alapján világos, hogy az α|M_{0} függvény BCF -be képez. Az egyszerűség kedvéért a továbbiakban a megszorított leképezést is α-val jelöljük, amennyiben nem okoz félreértést. Ekkor fennáll az alábbi állítás.

2.1.6. Állítás. Az α : M_{0}(ℝ_{0}^{+}, B_{0}) ↦ BCF leképezés bijekció.

Bizonyítás. α(μ_{1}) = α(μ_{2}) ⇔ μ_{1}((0, t]) = μ_{2}((0, t]) minden t > 0-ra⇔ μ_{1}((s, t]) = μ_{2}((s, t]) minden 0 < s ≤ t-re. Tekintsük a G = {(s, t] : 0 ≤ s ≤ t} félgyűrűt. Ekkor σ(G) = B_{0} és közismert, hogy félgyűrűn értelmezett egészen σ-véges mérték egyértelműen terjed ki a generált σalgebrára (márpedig μ_{1} és μ_{2} M_{0}-beli elemek lévén egyaránt egészen σ-végesek), így az ekvivalens átalakítások fenti sora így folytatható:

⇔ μ_{1}|G =μ_{2}|G ⇔ \overline{μ_{1}|G} = \overline{μ_{2}|G} ⇔ μ_{1} = μ_{2}.

Ezzel beláttuk az injektivitást.

Legyen most f egy korlátos számlálófüggvény és legyen G az iménti intervallum-félgyűrű. Értelmezzük a μ pontmértéket az alábbiak szerint: 0 < s ≤ t esetén μ((s, t]) := f(t) − f(s). Legyen C = dis(f), ahol dis(f) jelöli az f függvény szakadási helyeinek a halmazát. Nyilván C legfeljebb megszámlálható számosságú halmaz. Ekkor μ((s, t]) = μ_{t} − μ_{s} = ∑_{x∈C∩(s,t]} 1 = ∑_{x∈C} δ_{x}((s, t]). Vagyis ezzel egy G-n értelmezett mértéket definiáltunk. C-nek nincs torlódási pontja végesben, ezért μ egészen σ-véges. Mivel félgyűrűről egyértelmű a kiterjesztés σ(G) = B_{0}-ra - és ez az egyértelmű kiterjesztés éppen a μ = ∑_{x∈C} δ_{x} képlettel meghatározott -, azért ezzel hozzá tudtunk rendelni egy μ ∈ M_{0}(ℝ_{0}^{+}, B_{0}) pontmértéket az f függvényhez. Most már csak az van hátra, hogy megmutassuk, hogy a μ által meghatározott számlálófüggvény éppen f. Ez azonban világos, hiszen t > 0 esetén μ_{t} = μ((0, t]) = μ(⋃_{n=1}^{∞}(\frac{1}{n}, t]) = lim_{n↦∞} μ((\frac{1}{n}, t]) = lim_{n↦∞}(f(t) − f(\frac{1}{n})) = f(t) − lim_{n↦∞} f(\frac{1}{n}) = f(t) − f(0) = f(t). Innen már μ_{0} = f(0) is világos, amivel bebizonyítottuk az állítást.

Ha pedig N egy M(ℝ^{+}, B_{0})-beli véletlen elem, akkor a fenti definícióval összhangban N -hez hozzárendelhető az alábbi sztochasztikus folyamat.

2.1.9. Definíció. Az N : Ω ↦ M(ℝ^{+}, B_{0}) pontfolyamat által meghatározott számlálófolyamat alatt az N_{t}(ω) := (N(ω))_{t}, (t ≥ 0) képlettel értelmezett sztochasztikus folyamatot értjük.

Világos, hogy így valóban egy sztochasztikus folyamatot értelmeztünk, hiszen 2.1.2 alapján N_{t} = N((0, t]) : Ω ↦ \overline{ℕ}_{0} mérhető leképezés. Ha most t-re, mint időparaméterre gondolunk, akkor tehát N_{t} szemlétesen a t időpontig megfigyelt érkezések véletlen számát fejezi ki, ami jól magyarázza az elnevezést. Az alábbi állítás a számlálófolyamatok tulajdonságait viszgálja.

2.1.7. Állítás. Tegyük fel, hogy az N pontfolyamat nem torlódó és legyen N_{t}, t ≥ 0 az általa meghatározott számlálófolyamat. Ekkor N_{(·)} 1 valószínűséggel az alábbi tulajdonságokkal rendelkezik:

1. N_{0} ∈ 0,

2. N_{(·)} : Ω ↦ ℕ_{0},

3. N_{(·)} trajektóriái monoton növő, jobbról folytonos, balról határértékkel rendelkező függvények 1 nagyságú ugrásokkal.

4. N_{t} : Ω ↦ ℕ_{0} valószínűségi változó minden t ≥ 0-ra.

Bizonyítás. Az 1.-3. pontok bizonyítása világos a determinisztikus esetre vonatkozó korábbi állítás alapján. Csupán az szorul igazolásra, hogy N_{t} = N((0, t]) : Ω ↦ ℕ_{0} mérhető leképezés minden t > 0-ra, ez pedig a 2.1.2 fennállása miatt világos. Ebből N_{0} mérhetősége is következik, mivel mérhető függvények limesze.

2.1.10. Definíció. Az N_{t}, t ≥ 0 sztochasztikus folyamatot korlátos számlálófolyamatnak nevezzük, ha majdnem biztosan teljesül rá, hogy N_{0} = 0 és, hogy a trajektóriái monoton növő, jobbról folytonos, balról határétékkel rendelkező függvények 1 nagyságú ugrásokkal.

A korlátos számlálófolyamatok terére bevezetjük a BCP jelölést. Az előzőek alapján világos, hogy ha az N pontfolyamat nem torlódó, akkor az általa meghatározott számlálófolyamat eleme lesz BCP -nek. Jelölje β azt a leképezést, ami egy adott M-beli N pontfolyamathoz hozzárendeli az általa meghatározott számlálófolyamatot. Ekkor β(N) = α ∘ N

. \begin{matrix}M_{0}(ℝ^{+}, B_{0}) & ⟶^{α} & BCF \\ ↑ N \\ Ω \end{matrix}

2.1.8. Állítás. β bijekció a nem torlódó pontfolyamatok és BCP között.

2.1.3. Megjegyzés. Itt abban az értelemben beszélünk bijekcióról, hogy két valószínűségi változót azonosnak tekintünk, ha azok 1 valószínűséggel megegyeznek.

Bizonyítás. Legyen N_{1} és N_{2} nem torlódó. Ekkor β(N_{1}) = β(N_{2}) m.b. ⇔ α(N_{1}) = α(N_{2}) m.b. ⇔ α(N_{1}) = α(N_{2}) és N_{1}, N_{2} ∈ M_{0} m.b.⇒ N_{1} = N_{2} és N_{1}, N_{2} ∈ M_{0} m.b ⇒ N_{1} = N_{2} m.b.

Ha pedig ξ ∈ BCP, akkor módosítsuk úgy ξ-t, hogy a fennmaradó nullmértékű halmazon legyen egyenlő az azonosan nulla függvénnyel. Jelölje a módosított folyamatot ξ′. Ekkor N(ω) = α^{−1}(ξ′(ω)) olyan nem torlódó pontfolyamatot értelmez, amely által meghatározott számlálófüggvény β(N) = α ∘ N = α ∘ α^{−1} ∘ ξ′ = ξ′, ez pedig m.b. egyenlő ξ-vel.

Újabb definíciókat ismertetünk az intenzitásmértékkel kapcsolatban, továbbra is a valós esetben maradva.

2.1.11. Definíció. Egy pontfolyamat G középértékfüggvénye alatt a G(t) = ν((0, t]) függvényt értjük.

Világos, hogy ekkor ν((a, b]) = (G(b) − G(a))χ(G(a) < ∞) és G(t) = 𝔼N((0, t]) = 𝔼N_{t}.

2.1.12. Definíció. A ν mérték Radon-Nikodym deriváltját a Lebesgue-mértékre nézve - amennyiben létezik - intenzitásfüggvénynek, vagy röviden intenzitásnak hívjuk. Ha az intenzitásfüggvény konstans és egyenlő α-val, akkor az α értéket inteniztáshányadnak nevezzük.

Megjegyezzük, hogy az általános esetben is definiálhatunk intenzitásfüggvényt a Lebesguemérték helyett egyéb alkalmas domináló mértéket választva.

Ha G(t) abszolút folytonos, akkor a μ(t) intenzitásfüggvényre μ(t) = G′(t) - azaz μ = \frac{dν}{dΛ}, vagy másképpen G(b) − G(a) = ∫_{a}^{b} μ(t)Λ(dt) - adódik, ahol Λ|B_{0} jelöli a Lebesgue-mértéket.

Tekintsünk most egy N Poisson folyamatot a pozitív valós számokon. Ha N egy olyan ν(B) mérték által meghatározott, ami a korlátos intervallumokon véges, akkor

P (N(B) < ∞) = ∑_{k=1}^{∞} \frac{ν(B)^{k}}{k!} e^{−ν(B)} = 1

minden B korlátos intervallumra. Legyen B_{i} = (i − 1, i], i = 1, 2, . . .. Ekkor a 2.1.4 állítás alapján azt kapjuk, hogy a ν által meghatározott N Poisson folyamat nem torlódó. Ekkor azonban, mint ahogy azt korábban láttuk, kölcsönönösen egyértelműen meghatározzák egymást a számlálófolyamatával. Emiatt ekvivalens módon beszélhetünk a két matematikai objektumról. Ez alapján ebben az esetben a Poisson folyamat számlálófolyamata a következőképpen jellemezhető.

1. ℙ(N_{0} = 0) = 1,

2. N_{t} − N_{s} Poisson eloszlású G(t) − G(s) = ∫_{s}^{t} μ(x)Λ(dx) paraméterrel ∀ 0 ≤ s ≤ t -re és

3. N_{t} független növekményű.

Szokás a Poisson folyamatot a számlálófolyamatával definiálni az imént látott módon. Ha most az intenzitásmérték ν(B) = λΛ(B), ahol Λ|B, akkor nyilván fenáll, hogy Λ véges a korlátos intervallumokon. Továbbá G(t) = λt, és

ℙ(N_{t} − N_{s} = k) =\frac{[λ(t − s)]^{k}}{k!} e^{−λ(t−s)}.

Ebben az esetben homogén Poisson folyamatról beszélünk. A 2.1 ábra egy homogén Poisson számlálófolyamat realizációit illetve a hozzájuk tartozó középértékfüggvényeket mutatja két különböző intenzitáshányad esetén.

Tetszőleges N valós állapotterű pontfolyamat esetén definiálhatók a

τ_{i} = inf{t : N_{t} ≥ i}, i = 0, 1, …

2.1. ábra. A t ↦ N_{t} trajektóriák egy-egy realizációja μ = 1 és μ = 1.5 esetén

valószínűségi változók, melyet a pontfolyamat i. érkezési idejének nevezünk. Könnyen látható, hogy nem torlódó pontfolyamat esetén τ_{i} nem más, mint a t ↦ N_{t} leképezés i-edik ugrásának "időpontja". Ilyenkor a számlálófolyamat trajektóriáinak jobbról folytonossága miatt inf helyett min is írható. Ekkor (τ_{i} < x) = (N_{x} ≥ i), így τ_{i} eloszlásfüggvénye

F_{τ_{i}}(x) = ℙ(τ_{i} < x) = ℙ(N_{x} ≥ i) = ∑_{j=i}^{∞} ℙ(N_{x} = j).

2.1.4. Példa. Legyen az N_{t} Poisson számlálófolyamat a véges értékű G középértékfüggvény által

meghatározott, ekkor

F_{τ_{i}}(x) = ℙ(τ_{i} < x) = ℙ(N_{x} ≥ i) = ∑_{j=i}^{∞} ℙ(N((0, x]) = j) = ∑_{j=i}^{∞} \frac{G(x)^{j}}{j!} e^{−G(x)}.

Speciálisan F_{τ_{1}}(x) = 1−ℙ(N((0, x]) = 0) = 1−e^{−G(x)}. Innen látható, hogy amennyiben N homogén Poisson pontfolyamat λ intenzitáshányaddal, úgy τ_{1} eloszlása λ paraméterű exponenciális.

2.1.13. Definíció. A ξ_{i} = τ_{i} − τ_{i−1} valószínűségi változót az i-edik köztes időnek nevezzük (i = 1, 2, …).

Megmutatható, hogy λ intenzitáshányadú homogén Poisson folyamat esetén a köztes idők páronként függetlenek és exponenciális eloszlásúak λ paraméterrel.

2.2. Duplán sztochasztikus (Cox) folyamatok

Az előző fejezetben megismerkedtünk a Poisson folyamat fogalmával, most pedig megvizsgáljuk ennek egyfajta általánosítási lehetőségét, amikor is a ν intenzitásmérték nem determinisztikus, hanem szintén a véletlentől függő mennyiség. Ehhez szükségünk lesz az úgynevezett véletlen mérték fogalmára. Ezen általánosítás eredményeként a pontfolyamatok egy bővebb osztályát kapjuk, az ún. duplán sztochasztikus, vagy más néven Cox folyamatokét. Ennek segítségével olyan modellek építhetők, melyek bizonyos esetekben pontosabban közelítik a valóságot a hagyományos Poisson modelleknél. A duplán sztochasztikus folyamat elnevezést Cox mutatta be először 1955-ös cikkében [3].

Legyen (S, B) mérhető tér és jelölje M_{R} ≡ M_{R}(S, B) a B -n értelmezett egészen σ-véges mértékek terét. Ahogy azt az előző fejezetben tettük, itt is definiálhatóak minden B ∈ B -re a

π_{B} : M_{R} ↦ [0, ∞)
μ ↦ μ(B)

projekciók. Ezután szintén az előző fejezetben látottakkal analóg módon legyen M_{R} a legszűkebb olyan σ-algebra, melyre nézve az összes π_{B}, B ∈ B leképezés (M_{R}, B([0, ∞))) -mérhető, azaz

M_{R}≡ M_{R}(S, B) = σ{π_{B} : B ∈ B}.

A továbbiakban használni fogjuk a B′ = B([0, ∞)) jelölést is ℝ_{0}^{+} Borel halmazainak σ-algebrájára.

2.2.1. Definíció. Egy ρ : Ω ↦ M_{R}(S, B) mérhető leképezést véletlen mértéknek nevezünk.

Most is ugyanúgy definiálhatóak a ρ(B) = π_{B} ∘ ρ egydimenziós marginálisok, mely segítségével a 2.1.1 állítással összhangban az alábbi karakterizáció foggalmazható meg:

2.2.1. Állítás. A következő két állítás ekvivalens:

(i) ρ : Ω ↦ M_{R} egy véletlen mérték

(ii) ρ(B) : Ω ↦ [0, ∞) (A, B([0, ∞))) -mérhető minden B ∈ B -re.

Tekintsük a 0 < λ < ∞ paraméterhez tartozó Poisson eloszlást a nemnegatív egészeken. Ez tehát i-hez a \frac{λ^{k}}{k!} e^{−λ} valószínűségeket rendeli. Ez kiterjeszthető az \overline{ℕ}_{0} halmazra úgy, hogy a ∞- hez 0 valószínűséget rendelünk. Ezután már értelmezhetjük a λ = 0 illetve a λ = ∞-hez tartozó elfajuló eseteket a lim_{λ↦0} \frac{λ^{k}}{k!} e^{−λ} és lim_{λ↦∞} \frac{λ^{k}}{k!} e^{−λ} határértékek segítségével. Ezáltal a 0-ra illetve a ∞-re koncentrált eloszlásokat kapjuk. Legyen tehát 0 < λ < ∞ esetén p_{k}(λ) = \frac{λ^{k}}{k!} e^{−λ}, ha k = 1, 2, … és p_{∞}(λ) = 0. Legyen továbbá

p_{k}(0) = \{ \begin{matrix}1 & , ha k = 0 \\ 0 & egyébként \end{matrix}

és

p_{k}(∞) = \{ \begin{matrix}1 & , ha k = ∞ \\ 0 & egyébként \end{matrix} .

Ha most ν ∈ M_{R}(S, B) tetszőleges, akkor bevezethetjük a P_{k}^{B}(ν) = p_{k}(ν(B)) jelölést. P_{k}^{B} tehát egy M_{R} ↦ ℝ függvény, ami a π_{B} leképezések mérhetőségéből kifolyólag szintén mérhető. Ezek után meghatározzuk a Cox folyamat fogalmát:

2.2.2. Definíció. Legyen ρ : Ω ↦ M_{R}(S, B) egy véletlen mérték. Ekkor a C : Ω ↦ M(S, B) pontfolyamatot duplán sztochasztikus vagy más néven Cox folyamatnak nevezzük, ha a diszjunkt B_{1}, B_{2}, …, B_{k} halmazokhoz tartozó egydimenziós marginálisai - C(B_{1}), C(B_{2}), …, C(B_{k}) - ρ-ra nézve feltételesen függetlenek és C(B) ρ melletti feltételes eloszlása ν(B) partaméterű Poisson eloszlású, azaz

ℙ(C(B) = k|ρ) = P_{k}^{B}(ρ).

Ekkor a teljes várható érték tételével felírhatjuk C(B) eloszlását:

ℙ(C(B) = k) = 𝔼(ℙ(C(B) = k|ρ)) = 𝔼(P_{k}^{B}(ρ)) = ∫_{M_{R}} P_{k}^{B}(ν)Q_{ρ}(dν),  (2.2)

Vagyis C(B) eloszlása keverék Poisson eloszlás a ρ keverő valószínűségi változóval. Innen meghatározhatjuk a C Cox folyamat γ intenzitásmértékét, hiszen amennyiben 𝔼(C(B)) < ∞, úgy igaz az alábbi összefüggés:

γ(B) = 𝔼(C(B)) = ∑_{k=0}^{∞} kℙ(C(B) = k) = ∑_{k=0}^{∞} k ∫_{M_{R}} P_{k}^{B}(ν)Q_{ρ}(dν) =

∫_{M_{R}} ∑_{k=0}^{∞} kP_{k}^{B}(ν)Q_{ρ}(dν) = ∫_{M_{R}} ν(B)Q_{ρ}(dν) = 𝔼(ρ(B))

A Cox folyamat fogalmának létezik egy másik definíciója, amely sokszor célszerűbb tárgyalást tesz lehetővé, ám kevésbé szemléletes, ezért most csak vázlatosan ismertetjük, rámutatva a két definíció közötti kapcsolatra.

2.2.3. Definíció. Legyen adott egy ζ : Ω ↦ Θ valószínűségi változó és egy N^{ϑ} pontfolyamat (S, B)- n minden ϑ ∈ Θ -ra úgy, hogy minden M ∈ M esetén a ϑ ↦ Q_{N^{ϑ}}(M) leképezés mérhető legyen. Azt mondjuk, hogy a

Q_{N} = ∫_{Θ} Q_{N^{ϑ}}(·)Q_{ζ}(dϑ)

eloszlással értelmezett N : Ω ↦ M(S, B) pontfolyamat keverék pontfolyamat Q_{ζ} keverő eloszlással (vagy ζ keverő valószínűségi változóval).

2.2.4. Definíció. Legyen adott egy ρ : Ω ↦ M_{R} véletlen mérték és jelölje N^{ν} minden ν ∈ M_{R} esetén azt a Poisson folyamatot, aminek intenzitásmértéke éppen ν. Ekkor az N^{ν} folyamatokból a ρ véletlen mérték, mint keverő valószínűségi változó segítségével képzett C : Ω ↦ M(S, B) keverék Poisson-folyamatot sztenderd paraméterezésű Cox folyamatnak nevezzük.

Ekkor tehát C eloszlása a

Q_{C} = ∫_{M_{R}} Q_{N^{ν}} (·)Q_{ρ}(dν)

összefüggéssel definiált.

2.2.1. Megjegyzés. Még általánosabb tárgyaláshoz jutunk, ha nem követeljük meg, hogy a keverő valószínűségi változó véletlen mérték legyen, hanem tetszőleges ζ : Ω ↦ Θ valószínűségi változó szerint keverjük az N^{ϑ}, ϑ ∈ Θ Poisson folyamatokat, feltéve persze, hogy ϑ ↦ Q_{N^{ϑ}}(M) mérhető minden M ∈ M esetén.

Most megvizsgáljuk, hogy milyen viszonyban állnak egymással a 2.2.4 és a 2.2.2 definíciók. Ehhez definiáljuk a ρ^{∗} : Ω ↦ M_{R} véletlen mértéket úgy, hogy a (C, ρ^{∗}) együttes eloszlása az alábbi legyen:

Q_{(C,ρ^{∗})}(M × B) = ∫_{B} Q_{N^{ν}} (M)Q_{ρ}(dν)  (2.3)

minden M ∈ M, B ∈ M_{R} esetén. Ekkor világos, hogy Q_{ρ^{∗}}(B) = Q_{(C,ρ^{∗})}(M ×B) = Q_{ρ}(B), vagyis ρ^{∗} =^{d} ρ. Bebizonyítjuk a következő lemmát:

2.2.1. Lemma. A fenti jelölések mellett C feltételes eloszlása a ρ^{∗} valószínűségi változóra nézve éppen Q_{N_{ρ^{∗}}} (.), azaz:

ℙ(C ∈ M|ρ^{∗} = ν) = Q_{N^{ν}} (M), M ∈ M.

Bizonyítás. A feltételes eloszlás definíciója szerint azt kell belátnunk, hogy

𝔼(Q_{N^{ρ^{∗}}} (M)χ(ρ^{∗} ∈ B)) = ℙ(C ∈ M, ρ^{∗} ∈ B),

∀M ∈ M, B ∈ M_{R} esetén. Ez azonban világos, hiszen a baloldal nem más, mint

∫_{B} Q_{N^{ν}} (M)Q_{ρ^{∗}}(dν),

míg a jobboldal 2.3 alapján ∫_{B} Q_{N^{ν}} (M)Q_{ρ}(dν), de ez a két integrál egyenlő, hiszen Q_{ρ^{∗}} = Q_{ρ}.

Ha most az N^{ν} paraméterezést megtartva a ρ^{∗} véletlen mérték segítségével definiálok egy C^{∗} Cox folyamatot, akkor nyilván C^{∗} =^{d} C, ezért

ℙ(C^{∗} ∈ M|ρ^{∗} = ν) = Q_{N^{ν}} (M)  (2.4)

is igaz lesz. Egy ilyen tulajdonságú (C^{∗}, ρ^{∗}) vektort szokás Cox-párnak nevezni. Határozzuk most meg a C^{∗}(B) valószínűségi változó ρ^{∗} = ν feltétel melletti eloszlását:

ℙ(C^{∗}(B) = k|ρ^{∗} = ν) = ℙ(π_{B} ∘ C^{∗} = k|ρ^{∗} = ν) = ℙ(C^{∗} ∈ π_{B}^{−1}({k})|ρ^{∗} = ν) =
Q_{N^{ν}} (π_{B}^{−1}({k})) = ℙ(N^{ν} ∈ π_{B}^{−1}({k})) = ℙ(π_{B} ∘ N^{ν} = k) = ℙ(N^{ν}(B) = k) = P_{k}^{B}(ν),

vagyis éppen a Cox folyamat 2.2.2 meghatározásában szerepelt definiáló tulajdonságot kaptuk vissza, mely rámutat a 2.2.4 definíció megközelítésének eredetére. A továbbiakban kimondunk majd egy tételt, mellyel a 2.2.2. definícióban szereplő F^{ρ} = σ{ρ} σ-algebrát jellemezhetjük. Ehhez először szükségünk lesz néhány előkészítő lemmára, amelyeket az alábbiakban bizonyítunk.

2.2.5. Definíció. Ha S nemüres halmaz és H ⊂ P(S), akkor legyen

H^{δ}= { ∩_{i=1}^{∞} A_{i} : A_{i} ∈ H (i = 1, 2, …)} és H^{σ}= { ∪_{i=1}^{∞} A_{i} : A_{i} ∈ H (i = 1, 2, …)}.

H_{δ} és H_{σ} tehát rendre a H halmazrendszert tartalmazó legszűkebb σ-metszetzárt illetve σuniózárt halmazrendszer. Legyen P_{0} = H, Q_{0} = H^{C}. Ha pedig α > 0 megszámlálható rendszám és P_{β}, illetve Q_{β} már definiált minden β < α esetén, akkor legyen

P_{α} = ( ⋃_{β<α}Q_{β})^{δ} és Q_{α} = ( ⋃_{β<α}P_{β})^{σ}.

Ekkor transzfinit indukcióval belátható, hogy P_{α} ⊂ B(H), Q_{α} ⊂ B(H) és P_{α}^{C} = Q_{α} minden releváns α rendszámra. A definíció közvetlen következménye még, hogy β < α < ω_{1} esetén P_{β} ⊂ Q_{α} és Q_{β} ⊂ P_{α}, amiből az is nyilvánvalóan látható, hogy ⋃_{α<ω_{1}}P_{α} = ⋃_{α<ω_{1}}Q_{α}. Jelölje ezt a közös uniót F. A fent említettek miatt H ⊂ F ⊂ B(H). Ekkor F komplementerzárt, hiszen

H ∈ F ⇔ ∃α < ω_{1} : H ∈ P_{α} ⇔ ∃α < ω_{1} : H^{C} ∈ P_{α}^{C} = Q_{α} ⇔ H^{C} ∈ F.

Másrészt F zárt a σ-unió-ra is, hiszen ha H_{i} ∈ F (i = 1, 2, …), akkor minden i-hez létezik α_{i} rendszám úgy, hogy H_{i} ∈ P_{α_{i}}. Ha most a γ megszámlálható rendszámot úgy választjuk meg, hogy nagyobb legyen minden α_{i}-nél (megtehetjük, hiszen például γ = α_{1} + α_{2} + … + 1 megfelelő), akkor H_{i} ∈ Q_{γ} ⊂ F. Tehát F olyan σ-algebra, amely tartalmazza H-t és része a legszűkebb H-t tartalmazó σ-algebrának, ezért csak azonos lehet ez utóbbival, azaz F ≡ B(H).

2.2.2. Lemma. Ha f : Ω ↦.X és G ⊂ P (X) egy nemüres halmazrendszer, akkor

1. {f^{−1}(H) : H ∈ G}^{δ} = {f^{−1}(H) : H ∈ G^{δ}} és {f^{−1}(H) : H ∈ G}^{σ} = {f^{−1}(H) : H ∈ G^{σ}}.

2. ⋃_{i∈I}{f^{−1}(H) : H ∈ G_{i}} = {f^{−1}(H) : H ∈ ⋃_{i∈I}G_{i}} tetszőleges I indexhalmazra.

Bizonyítás.

1. Csak az állítás első felét mutatjuk meg, a második ezzel teljesen analóg módon kivitelezhető. {f^{−1}(H) : H ∈ G}^{δ} = { ∩_{i=1}^{∞} f^{−1}(H_{i}) : H_{i} ∈ G (i = 1, 2, …)} = {f^{−1}(∩_{i=1}^{∞}H_{i}) : H_{i} ∈ G (i = 1, 2, …)} = {f^{−1}(H) : H ∈ G ^{δ}}.

2. ⋃_{i∈I}{f^{−1}(H) : H ∈ G_{i}} = {A : ∃i ∈ I, amire A = f^{−1}(H) : H ∈ G_{i}} = {f^{−1}(H) : H ∈ ⋃_{i∈I}G_{i}}.

2.2.3. Lemma. Legyen f : Ω ↦.X és G ⊂ P (X) egy nemüres halmazrendszer. Ekkor σ{f^{−1}(H) : H ∈ G} = {f^{−1}(H) : H ∈ σ(G)}

Bizonyítás. Legyen H = {f^{−1}(H) : H ∈ G}. Az imént ismertetett transzfinit indukciós eljárás segítségével előállítjuk a H halmazrendszert tartalmazó legszűkebb σ-algebrát. Ehhez elkészítjük a P_{0}^{∗} = G, Q_{0}^{∗} = G^{C}, illetve ezzel párhuzamosan a P_{0} = H ={f^{−1}(H) : H ∈ P_{0}^{∗}} és Q_{0} = H^{C} = {(f^{−1}(H))^{C} : H ∈ G} = {f^{−1}(H^{C}) : H ∈ G} = {f^{−1}(H) : H ∈ G^{C}} = {f^{−1}(H) : H ∈ Q_{0}^{∗}} halmazokat, majd a fent látottakhoz hasonlóan a P_{α}^{∗}, Q_{α}^{∗} illetve a P_{α} és Q_{α} halmazokat. Legyen most α > 0 megszámlálható rendszám és tegyük fel, hogy minden β < α esetén teljesül a

P_{β}={f^{−1}(H) : H ∈ P_{β}^{∗}} és Q_{β} = {f^{−1}(H) : H ∈ Q_{β}^{∗}}

azonosság. Ekkor a 2.2.2. lemma alapján

P_{α} = ( ⋃_{β<α}Q_{β})^{δ} = ( ⋃_{β<α}{f^{−1}(H) : H ∈ Q_{β}^{∗}})^{δ} =

({f^{−1}(H) : H ∈ ⋃_{β<α}Q_{β}^{∗}})^{δ} = {f^{−1}(H) : H ∈ ( ⋃_{β<α}Q_{β}^{∗})^{δ}} = {f^{−1}(H) : H ∈ P_{α}^{∗}}

és teljesen hasonló módon Q_{α} = {f^{−1}(H) : H ∈ Q_{α}^{∗}} is igaz. Innen

σ{f^{−1}(H) : H ∈ G} = ⋃_{α<ω_{1}}P_{α} = ⋃_{α<ω_{1}} {f^{−1}(H) : H ∈ P_{α}^{∗}} =
{f^{−1}(H) : H ∈ ⋃_{α<ω_{1}}P_{α}^{∗}} = {f^{−1}(H) : H ∈ σ(G)}.

Ezzel bebizonyítottuk ezt a lemmát is.

Térjünk vissza most a 2.2.2 definícióra. Amennyiben a G halmazrendszer a B egy generátorrendszere, és rendre F^{ρ}, F^{B}, és F^{G} jelöli a σ{ρ}, σ{ρ(B) : B ∈ B}, illetve σ{ρ(B) : B ∈ G} σ-algebrákat, akkor az alábbi állítás fogalmazható meg.

2.2.1. Tétel. Tegyük fel, hogy a G generátorrendszer olyan, hogy S előállítható megszámlálható sok G-beli halmaz uniójaként. A fenti jelölések mellett:

F^{ρ} ≡ F^{B} ≡ F^{G}.

Bizonyítás. F^{G} = σ{ρ(B) : B ∈ G} = σ{π_{B} ∘ ρ : B ∈ G} = σ{ρ^{−1}(π_{B}^{−1}(C)) : B ∈ G, C ∈ B′} =

σ{ρ^{−1}(H) : H ∈ {π_{B}^{−1}(C) : B ∈ G, C ∈ B′}} = {ρ^{−1}(H) : H ∈ σ{π_{B}^{−1}(C) : B ∈ G, C ∈ B′}} a 2.2.3. lemma miatt. Megmutatjuk, hogy

{ρ^{−1}(H) : H ∈ σ{π_{B}^{−1}(C) : B ∈ G, C ∈ B′}} = {ρ^{−1}(H) : H ∈ σ{π_{B}^{−1}(C) : B ∈ σ(G), C ∈ B′}},

ehhez elég belátni, hogy σ{π_{B}^{−1}(C) : B ∈ G, C ∈ B′} = σ{π_{B}^{−1}(C) : B ∈ σ(G), C ∈ B′}, azaz, hogy σ{π_{B} : B ∈ G} = σ{π_{B} : B ∈ σ(G)}. Az világos, hogy

σ{π_{B} : B ∈ G} ⊂ σ{π_{B} : B ∈ σ(G)},  (2.5)

látni kell még a fordított irányú tartalmazást. Ehhez legyen az E =σ{π_{B} : B ∈ G} jelölés bevezetésével F^{∗} = {B ∈ σ(G) : π_{B} E-mérhető}. Ekkor nyilván G ⊂ F^{∗} ⊂ σ(G). Megmutatjuk, hogy F^{∗} σ-algebra, amiből F^{∗} = σ(G) adódik, azaz σ{π_{B} : B ∈ σ(G)} ⊂ E, mely (2.5)-tel együtt adja a kívánt állítást.

1. S előállítható megszámlálható sok G-beli halmaz uniójaként.

2. F^{∗} különbségzárt. Legyen A, B ∈ F^{∗}, A ⊂ B, ekkor π_{B\A}(ν) = ν(B\A) = (ν(B) − ν(A))χ(ν(A) < ∞) = (π_{B}(ν)−π_{A}(ν))χ(π_{A}(ν) < ∞) és mivel mérhető függvények különbsége illetve szorzata is mérhető, ezért B\A ∈ F^{∗}

3. F^{∗} zárt a megszámlálható diszjunkt unióra. Legyenek az A_{i} ∈ F^{∗}, (i = 1, 2, …) halmazok páronként diszjunktak, akkor π_{∪A_{i}}(ν) = ν(∪A_{i}) = ∑_{i} ν(A_{i}) = ∑_{i} π_{A_{i}}(ν) és megszámlálható sok mérhető függvény szintén mérhető, így ⋃_{i}A_{i} ∈ F^{∗}.

Tehát a 2.1.1 lemma miatt F^{∗} valóban σ-algebra. Folytatva a tétel bizonyítását, kimondhatjuk, hogy

σ{ρ^{−1}(H) : H ∈ {π_{B}^{−1}(C) : B ∈ G, C ∈ B′}} =
σ{ρ^{−1}(H) : H ∈ {π_{B}^{−1}(C) : B ∈ σ(G), C ∈ B′}} = F^{B}.

Hátra van még az F^{ρ} ≡ F^{B} azonosság igazolása. σ{ρ} = σ{ρ^{−1}(H) : H ∈ M_{R}} = σ{ρ^{−1}(H) : H ∈ {π_{B}^{−1}(C) : B ∈ σ(G), C ∈ B′}} definíció szerint, erről viszont az imént láttuk, hogy egyenlő F^{B}-vel. Ezzel a tételünket beláttuk.

Most térjünk rá arra az esetre, amikor (S, B) ≡ (ℝ^{+}, B_{0}˙), a G generátorrendszer pedig az alulról nyílt, felülről zárt intervallumok félgyűrűje. Értelmezzük a ρ_{t}, t ≥ 0 sztochasztikus folyamatot a következőképpen:

ρ_{t} = ρ((0, t]),

ha t > 0 és ρ_{0} legyen ρ_{0+}.

2.2.6. Definíció. Az imént bemutatott sztochasztikus folyamatot a C Cox folyamat középértékfolyamatának nevezzük.

Ekkor

σ(ρ) = σ{ρ(B) : B ∈ G} = σ{ρ((s, t]) : 0 < s ≤ t} =σ{ρ_{t} − ρ_{s} : 0 < s ≤ t} = : L.

ρ_{0} értelmezéséből adódóan könnyen látható, hogy L = σ{ρ_{t} − ρ_{s} : 0 ≤ s ≤ t}. Másrészt, ha L^{∗} = σ{ρ_{t} : 0 ≤ t}, akkor egyrészt L^{∗} ⊂ L, másrészt mivel ρ_{t} : 0 ≤ t L^{∗}-mérhető, azért ρ_{t} −ρ_{s} : 0 ≤ s ≤ t is, így L ⊂ L^{∗} is igaz. Tehát

σ(ρ) = σ{ρ_{t} : 0 ≤ t}.  (2.6)

Szintén természetes kérdésként vetődik fel, hogy mikor lesz egy Cox folyamat nem torlódó. Erre adunk most egy elégséges feltételt. A 2.1.4 állítás alapján elég lenne azt látni, hogy

ℙ(C(B) < ∞) = 1

minden B korlátos intervallumra. Fejezzük ki ezért ezt a valószínűséget (2.2) segítségével:

ℙ(C(B) < ∞) = ∑_{k=0}^{∞} ℙ(C(B) = k) = ∑_{k=0}^{∞} ∫_{M_{R}} P_{k}^{B}(ν)Q_{ρ}(dν) =
∫_{M_{R}} ∑_{k=0}^{∞} P_{k}^{B}(ν)Q_{ρ}(dν) = ∫_{M_{R}} χ(π_{B}(ν) < ∞)Q_{ρ}(dν) =
∫_{{π_{B}<∞}} 1Q_{ρ} = Q_{ρ}(π_{B} < ∞) = ℙ(ρ(B) < ∞).

Ha tehát ρ olyan, hogy 1 valószínűséggel véges minden korlátos intervallumra, akkor a C Cox folyamat nem torlódó. Ebben az esetben tehát az előző alfejezetben látottak miatt ekvivalens módon beszélhetünk magáról a folyamatról és annak számlálófolyamatáról. Ekkor (2.6)-ot figyelembe véve a C Cox folyamatot definiáló tulajdonságok a C_{t}, t ≥ 0 számlálófolyamatra az alábbi formában fogalmazhatóak meg.

1. C növekményei feltételesen függetlenek a σ{ρ_{t} : 0 ≤ t} σ-algebra mellett.

2. ℙ(C_{t} − C_{s} = k|ρ_{u}, 0 ≤ u) = \frac{(ρ_{t}−ρ_{s})^{k}}{k!} e^{−(ρ_{t}−ρ_{s})} minden 0 < s ≤ t esetén.

Ha pedig ρ m.b. monoton folytonos, akkor még

ℙ(C_{t} = k|ρ_{u}, 0 ≤ u) = \frac{ρ_{t}^{k}}{k!} e^{−ρ_{t}}

is igaz, természetesen az egyenlőséget majdnem biztosan értve.

Egy valós állapotterű Cox folyamat G középértékfüggvényére az intenzitásmértékre levezetett korábbi képlet alapján

G(t) = γ((0, t]) = ∫_{M_{R}} ν((0, t])Q_{ρ}(dν) = ∫_{M_{R}} G^{ν}(t)Q_{ρ}(dν)

adódik valahányszor 𝔼(ρ_{t}) < ∞, ahol G^{ν}(t) jelöli a ν intenzitásmértékű Poisson folyamat középértékfüggvényét, ρ_{t} pedig az imént bemutatott középértékfolyamatot.

2.2.7. Definíció. A C intenzitás-folyamata alatt az alábbi η_{t}, t > 0 sztochasztikus folyamatot értjük

η_{t} = \{ \begin{matrix}\frac{dρ}{dΛ}(t) & , ha ρ << Λ \\ 0 & egyébként \end{matrix} .

2.2.2. Megjegyzés. A szokásos η_{0} := η_{0+} értelmezéssel az intenzitásfolyamat is kiterjeszthető ℝ_{0}^{+}-ra.

Határozzuk most meg a C intenzitásfüggvényét. Tegyük fel, hogy a ρ véletlen mérték majdnem biztosan abszolút folytonos a Lebesgue-mértékre, azaz, hogy ℙ(ρ << Λ) = 1. Amennyiben a ν egészen σ-véges mérték abszolút folytonos a Λ Lebesgue mértékre, úgy jelöljük μ^{ν}-vel a \frac{dν}{dΛ} Radon-Nykodim deriváltat, ellenkező esetben μ^{ν} legyen azonosan 0. Ha a (ν, t) ↦ μ^{ν}(t) leképezés mérhető, akkor Fubini tételének értelmében fennáll a következő

γ(B) = ∫_{M_{R}} ν(B)Q_{ρ}(dν) = ∫_{{ν:ν<<Λ}} ν(B)Q_{ρ}(dν) =
∫_{{ν:ν<<Λ}} ∫_{B} μ^{ν}(t)Λ(dt)Q_{ρ}(dν) = ∫_{M_{R}} ∫_{B} μ^{ν}(t)Λ(dt)Q_{ρ}(dν) =

∫_{B} ∫_{M_{R}} μ^{ν}(t)Q_{ρ}(dν)Λ(dt),

vagyis a μ intenzitásfüggvény a

μ(t) = ∫_{M_{R}} μ^{ν}(t)Q_{ρ}(dν)

képlettel adható meg. Cox folyamatra példaként szolgál majd a második fejezet 6. alfejezetében bemutatásra kerülő Markov-modulált Poisson folyamat.

2.3. Néhány szó mátrix értékű függvények differenciálásáról

és a mátrix exponenciális függvényről

Ebben az alfejezetben néhány alapvető fogalmat és azonosságot ismertetünk mátrix értékű függvények differenciálásával kapcsolatban, majd ezután mátrixokra vonatkozó lineáris differenciálegyenleteket vizsgálunk, melyek megoldásának explicit alakban történő felírásához értelmezzük az e^{A} kifejezést, ahol a kitevőben egy A négyzetes mátrix áll. Ennek kapcsán vizsgáljuk majd az ún. mátrix exponenciális függvényt és ismertetünk egy numerikus eljárást is, amely segítségével hatékonyan számolható exp(A). Ezek alapvető elemeit képezik majd a későbbi tárgyalásunknak.

2.3.1. Definíció. Az A : ℝ ↦ ℝ^{m×n} mátrix értékű függvény deriváltja alatt az A′ : ℝ ↦ ℝ^{m×n} függvényt értjük, ahol (A′)_{ij} = A′_{ij}.

2.3.1. Lemma. Legyen adott az A : ℝ ↦ ℝ^{m×n} és B : ℝ ↦ ℝ^{n×p} mátrix értékű függvény, ekkor (AB)′ = A′B + AB′.

Bizonyítás.

[(AB)′]_{ij} = (AB)′_{ij} = (∑_{p=1}^{n}A_{ip}B_{pj})′ = ∑_{p=1}^{n} (A_{ip}B_{pj})′
= [(AB)′]_{ij} = (AB)′_{ij} = (∑_{p=1}^{n}A_{ip}B_{pj})′ = ∑_{p=1}^{n} (A_{ip}B_{pj})′
= ∑_{p=1}^{n}A′_{ip}B_{pj} + ∑_{p=1}^{n}A_{ip}B′_{pj} = ∑_{p=1}^{n} (A′)_{ip} B_{pj} + ∑_{p=1}^{n}A_{ip} (B′)_{pj}
= (A′B)_{ij} + (AB′)_{ij} = (A′B + AB′)_{ij}.

2.3.1. Következmény. Ha az A : ℝ ↦ ℝ^{m×n} és a C : ℝ ↦ ℝ^{p×q} mátrix értékű függvény konstans, akkor - tekintve, hogy A′ = 0 és C′ = 0 - a B : ℝ ↦ ℝ^{n×p} mátrix értékű függvény deriváltjára (AB)′ = AB′ és (BC)′ = B′C.

2.3.2. Lemma. Ha A : ℝ ↦ ℝ^{r×r} mátrix értékű függvény, akkor n ≥ 2 esetén

(A^{n})′ = ∑_{j=1}^{n−1} A^{j}A′A^{n−1−j}.

Bizonyítás. A 2.3.1 lemma alapján teljes indukcióval világos.

2.3.3. Lemma. Az A : ℝ ↦ ℝ^{r×r} mátrix értékű függvény inverzének deriváltjára

(A^{−1})′ = −A^{−1}A′A^{−1}.

Bizonyítás. A^{−1}(t)A(t) = I minden t ∈ ℝ esetén, ebből a 2.3.1 lemma alapján (A^{−1})′ A+A^{−1}A′ = 0, amiből átrendezés, majd A^{−1}-el való szorzás után adódik az állítás.

2.3.2. Definíció. Legyen A r × r-es valós mátrix, ekkor A exponenciálisa alatt az

exp(A) = ∑_{k=0}^{∞} \frac{A^{k}}{k!}

hatványsorral értelmezett r × r-es mátrixot értjük.

Az így értelmezett hatványsor mindig konvergens – a konvergenciát tetszőleges normára nézve definiálhatjuk, mivel véges dimenziós vektortérben minden norma ekvivalens – , fennáll ugyanis az alábbi lemma.

2.3.4. Lemma. Ha A_{0}, A_{1}, … m × n-es mátrixok, ∥·∥ tetszőleges mátrix norma és a ∑_{k=0}^{∞} ∥A_{k}∥ numerikus sor konvergens, akkor ∑_{k=0}^{∞} A_{k} is az.

Bizonyítás. A Cauchy kovergenciakritérium és a háromszög-egyenlőtlenség miatt

∥∑_{k=m+1}^{n}A_{k}∥ ≤ ∑_{k=m+1}^{n} ∥A_{k}∥ < ε,

ha m, n ≥ N(ε).

Innen már a ∑_{k=0}^{∞} \frac{A^{k}}{k!} hatványsor konvergenciája világos a ∑_{k=0}^{∞} \frac{∥A∥^{k}}{k!} sor konvergenciájából. Könnyen ellenőrizhetően igazak az alábbi lemmában felsorolt tulajdonságok.

2.3.5. Lemma. Az A ↦ exp(A) leképezés rendelkezik az alábbi tulajdonságokkal.

1. A csupa nulla mátrixra exp(0) = I.

2. exp(A) nemszinguláris, és inverze exp(−A).

3. exp(A) = [exp(\frac{A}{m})]^{m}.

4. exp(A + B) = exp(A) exp(B) akkor és csak akkor, ha A és B felcserélhetők.

2.3.3. Definíció. Legyen adott az A ∈ ℝ^{r×r} mátrix, ekkor az A által meghatározott mátrix exponenciális függvény alatt a t ↦ exp(At) összefüggéssel értelmezett mátrix értékű függvényt értjük.

2.3.6. Lemma. A mátrix exponenciális függvény rendelkezik az alábbi tulajdonságokkal.

1. exp[A(s + t)] = exp(As) exp(At).

2. t ↦ exp(At) inverze t ↦ exp(−At).

3. exp [(A + B)t] = exp(At) exp(Bt) akkor és csak akkor, ha A és B felcserélhetők.

4. \frac{∂}{∂t} exp(At) = A exp(At) = exp(At)A.

Bizonyítás. Az 1-3. állítások közvetlen következményei a (2.3.5) lemmának, míg a 4. állítás szintén egyszerűen adódik a exp(At) = ∑_{k=0}^{∞} \frac{(At)^{k}}{k!} hatványsor tagonkénti deriválásából.

A 2.3.6 lemma 4. tulajdonsága más szóval azt mondja, hogy ha A egy r × r-es mátrix és X(t) = exp(At), akkor az X mátrix értékű függvény kielégíti az X(0) = I kezdeti feltétel melletti

X′(t) = AX(t)

lineáris differenciál-egyenletet. A közönséges differenciál-egyenletek elméletéből jól ismert, hogy ennek az egyenletnek az említett kezdeti feltétel mellett ez az egyetlen megoldása létezik.

A későbbiek során szükségünk lesz a következő lemmára, ami egy bizonyos speciális alakú mátrixhoz tartozó exponenciális függvény kiszámításában nyújt segítséget.

2.3.7. Lemma. Legyenek adottak az r × r-es A és B mátrixok, majd definiáljuk a 2r × 2r-es C mátrixot a következő alakban:

C = ( \begin{matrix}A & B \\ 0 & A \end{matrix} ).

Ekkor

exp(Ct) = ( \begin{matrix}F (t) & G(t) \\ 0 & F (t) \end{matrix} ),  (2.7)

ahol F (t) = exp(At) és G(t) = ∫_{0}^{t} exp(A(t − s))B exp(As)ds.

Bizonyítás. A mátrix exponenciális függvény definíciójából világos, hogy exp(Ct) (2.7) alakú. Másrészt exp(Ct) az X′(t) = CX(t) differenciál-egyenlet egyértelmű megoldása az X(0) = I kezdeti feltétel mellett. Ezt felírva

( \begin{matrix}F′(t) & G′(t) \\ 0 & F′(t) \end{matrix} ) = ( \begin{matrix}A & B \\ 0 & A \end{matrix} ) ( \begin{matrix}F (t) & G(t) \\ 0 & F (t) \end{matrix} ) = ( \begin{matrix}AF (t) & AG(t) + BF (t) \\ 0 & AF (t) \end{matrix} ),

vagyis F (t) és G(t) az F (0) = I és G(0) = 0 kezdeti feltételek melletti

F′(t) = AF (t)
G′(t) = AG(t) + BF (t)

differenciál-egyenletrendszer egyértelmű megoldásaként adódik. Ebből F (t) = exp(At) nyilvánvaló, G(t) = ∫_{0}^{t} exp(A(t − s))B exp(As)ds esetén pedig

G′(t) = A exp(At) ∫_{0}^{t} exp(−As)B exp(As)ds + exp(At) exp(−At)B exp(At)
= AG(t) + BF (t).

ezzel bebizonyítottuk az állítást.

A következő segédtételünk ugyanezt a lemmát általánosítja, bár egy későbbi alkalmazásra való tekintettel kissé más megközelítésben mondjuk ki. Előtte vezessük be a következő jelöléseket. Legyen adott egy r × r-es A mátrix, és az I ={1, 2, …, r} indexhalmaz I illetve J részhalmazai. Ekkor A_{I,J}-vel fogjuk jelölni azt az |I| × |J|-es részmátrixot, amit az eredeti mátrixból kapunk a megadott indexhalmazok szerinti sorok és oszlopok kiválasztásával.

2.3.8. Lemma. Legyen A és B két r × r-es mátrix és tegyük fel, hogy az F_{n}(t) mátrix értékű függvényekből álló sorozat kielégíti a következő rekurziót:

F_{0}(0) = I és F_{n}(0) = 0, ha n ≥ 1,
F′_{0}(t) = F_{o}(t)A,
F′_{n}(t) = F_{n}(t)A + F_{n−1}(t)B, ha n ≥ 1.

Ekkor fenáll az alábbi két állítás.

1. F_{n}(t) explicit módon felírható az alábbi formában:

F_{0}(t) = exp(At)
F_{n}(t) = ∫_{0}^{t} F_{n−1}(s)B exp(A(t − s))ds, ha n ≥ 1.

2. Tekintsük az alábbi (n + 1)r × (n + 1)r-es mátrixot:

C = ( \begin{matrix}A & B & 0 & ⋯ & 0 \\ 0 & A & B & ⋯ & 0 \\ 0 & 0 & A & ... & 0 \\ ⋮ & ⋮ & ⋮ & ... & B \\ 0 & 0 & 0 & 0 & A \end{matrix} )

Legyen D(t) = exp(Ct) és 0 ≤ k ≤ n esetén definiáljuk az I_{k} = {kr + 1, kr + 2, …, kr + r} indexhalmazokat. Ekkor

F_{k}(t) = [D(t)]_{I_{0},I_{k}}.

Bizonyítás.

1. F_{0}(t) = exp(At) nyilvánvaló, ha pedig n ≥ 1, akkor tegyük fel, hogy F_{n−1}-et már meghatároztuk, ekkor az F′_{n}(t) = F_{n}(t)A + F_{n−1}(t)B egyenletnek egyértelműen létezik megoldása az F_{n}(0) = 0 kezdeti feltétel mellett, mégpedig F_{n}(t) = ∫_{0}^{t} F_{n−1}(s)B exp(A(t − s))ds, hiszen erre a függvényre

F′_{n}(t) = ∫_{0}^{t} F_{n−1}(s)B exp(−As)ds exp(At)A + F_{n−1}(t)B exp(−At) exp(At)
= F_{n}(t)A + F_{n−1}(t)B.

2. Legyen n rögzített. Ekkor az F_{0}(t), F_{1}(t), …, F_{n}(t) értékek egyértelműen meghatározhatók az első n + 1 rekurziós egyenletből. Vezessük be az

X(t) = ( \begin{matrix}F_{0}(t) & F_{1}(t) & F_{2}(t) & ⋯ & F_{n}(t) \\ 0 & F_{0}(t) & F_{1}(t) & ⋯ & F_{n−1}(t) \\ 0 & 0 & F_{0}(t) & ... & ⋮ \\ ⋮ & ⋮ & ⋮ & ... & F_{1}(t) \\ 0 & 0 & 0 & 0 & F_{0}(t) \end{matrix} )

új ismeretlent. Ekkor könnyen láthatóan a rekurziós egyenletek éppen azt mondják ki, hogy

X′(t) = X(t)C

és hogy X(0) = I. Ekkor tehát, X(t) = exp(Ct) és tekintve, hogy F_{k}(t) = [X(t)]_{I_{0},I_{k}} az állítás bebizonyíttatott.

2.3.1. Megjegyzés. A 2.3.7. lemmának más irányú átalalánosítása is létezik, lásd például [29] 1.Tételét.

A mátrix exponenciális kiszámítására rengeteg numerikus módszer létezik, a témakörről alapos áttekintést nyújt a [30] alatti dolgozat. Ezek közül mutatunk most itt be egyet, amit használni fogunk egy későbbi algoritmusunk során.

2.3.4. Definíció. Legyen A adott r × r-es négyzetes mátrix, ekkor exp(A) (p, q)-Padé approximációja alatt az

R_{pq}(A) = [D_{pq}(A)]^{−1}N_{pq}(A)

mennyiséget értjük, ahol

N_{pq}(A) = ∑_{j=0}^{p} \frac{(p + q − j)!p!}{(p + q)!j!(p − j)!}A^{j}

és

D_{pq}(A) = ∑_{j=0}^{p} \frac{(p + q − j)!q!}{(p + q)!j!(q − j)!}(−A)^{j}.

Van Loan [30]-ban megemlíti, hogy D_{pq}(A) regularitása biztosított, amennyiben p és q bizonyos értelemben elég nagyok vagy A sajátértékei negatívak. exp(A) Padé approximációval történő kiszámítása tehát lényegében egy megfelelő p és q érték meghatározásából, majd R_{pq}(A) kiszámításából áll. [30] részletesen kifejti, hogy milyen módon érdemes választani p-t és q-t a különböző pontossági szintek eléréséhez. A Matlab^{TM}-ban implementált "expm" eljárás a (6, 6) párhoz tartozó Padé approximációt használ. Először azonban átskálázza a mátrixot úgy, hogy ∞ normája \frac{1}{2}-nél kisebb legyen, majd a 2.3.5 lemma 3. tulajdonságát használva visszaállítja az eredeti skálázást.

2.3.1. Algoritmus (Padé-approximáció a mátrix exponenciális kiszámítására).

1. Bemenet: A négyzetes mátrix.

2. Válaszd meg az \frac{1}{2} ≤ f < 1 és e számokat úgy, hogy ∥A∥_{∞} = f2^{e} legyen. Ekkor s = 0 ∨ (e + 1) esetén \frac{∥A∥_{∞}}{2^{s}} < \frac{1}{2}.

3. E := R_{(6,6)}(A)

4. i = 1, 2, …, s esetén E := E^{2}.

5. Kimenet : E approximáció exp(A)-ra.

2.4. Véges állapotterű, folytonos idejű Markov folyamatok

2.4.1. Definíció. Azt mondjuk, hogy a ξ_{t} folyamat folytonos idejű és véges állapotterű Markov folyamat az F_{t}, t ≥ 0 filtrációra nézve, ha t ≥ 0 esetén ξ_{t} az E = {1, 2, …, r} értékek valamelyikét veszi fel és i ∈ E illetve t ≥ s esetén

ℙ(ξ_{t} = i|F_{s}) = ℙ(ξ_{t} = i|ξ_{s}).

A mi esetünkben az F_{t}, t ≥ 0 filtráció a ξ_{t} folyamat által generált filtráció lesz majd. Feltesszük, hogy az átmenet-valószínűségek stacionáriusak, azaz, hogy P_{ij}(t) = ℙ(ξ_{t+s} = j|ξ_{s} = i) = P_{ij}(t) = ℙ(ξ_{t} = j|ξ_{0} = i), ∀t, s ≥ 0-ra.

2.4.1. Lemma. A fenti feltételek mellett

1. P_{ij}(t) ≥ 0,

2. ∑_{j=1}^{r} P_{ij}(t) = 1, i, j ∈ E

3. P_{ij}(s + t) = ∑_{u=1}^{r} P_{iu}(s)P_{uj}(t), t, s ≥ 0-ra (Chapman-Kolmogorov összefüggés).

Bizonyítás. Az első két pont triviális, a harmadik a teljes valószínűség tételéből következik.

A továbbiakban feltesszük még, hogy P_{ij}(t) a nulla pontban jobbról differenciálható és P_{ij}(0+) = δ_{ij}. Legyen P = (P_{ij}). Ekkor a Chapman-Kolmogorov egyenletek eképpen írhatók fel mátrix alakban:

P (t + s) = P (t)P (s).

A deriválhatóság miatt létezik olyan Q mátrix, amire P (t) = exp(Qt).

2.4.2. Definíció. A Q mátrixot hívjuk a ξ_{t} Markov folyamat infinitezimális generátorának.

Mivel P (t) sorösszege 1, ezért Q sorösszege 0. Emiatt a modell r(r − 1) független paraméterrel rendelkezik. A Q mátrix elemei a főátlóbelieket kivéve mind nemnegatívak, tehát a főátlóban lévő elemek nempozitívak. Ezekre szokás bevezetni a q_{i} = −q_{ii} jelölést. Ezek segítségével a Markov folyamat fejlődése a következőképpen írható le. Ha a folyamat éppen az i-edik állapotban van, akkor ott töltött ideje q_{i} paraméterű exponenciális eloszlású, majd innen a \frac{q_{ij}}{q_{i}} valószínűségek szerint vált át valamelyik másik j ≠ i állapotba.

2.4.3. Definíció. A ξ_{t} Markov folyamat irreducibilis, ha tetszőleges állapotból el lehet jutni bármely másikba, azaz ha tetszőleges i, j ∈ E esetén létezik t ≥ 0, amire P_{ij}(t) > 0.

Definiáljuk egy adott i állapoton történő első áthaladás idejét a

T_{i} = inf{t > 0 : ξ_{t} = i, ξ_{t−} ≠ i}

összefügéssel. Legyen továbbá F_{ij} = ℙ(T_{j} < ∞|ξ_{0} = i).

2.4.4. Definíció. Azt mondjuk, hogy az i állapot pozitív ismétlődő állapot, ha F_{ii} < ∞ és 𝔼(T_{i}|ξ_{0} = i) < ∞.

2.4.5. Definíció. A ξ_{t} Markov folyamat ergodikus, ha minden állapot pozitív ismétlődő.

2.4.1. Tétel. Ha a ξ_{t}, t ≥ 0 Markov folyamat irreducibilis és legalább az egyik állapot pozitív ismétlődő, akkor minden állapot pozitív ismétlődő és F_{ij} = 1 minden i, j ∈ E-re.

2.4.6. Definíció. Az r dimenziós π = (π_{1}, … π_{r}) sorvektort valószínűszínűségi vektornak hívjuk, ha π_{i} ≥ 0 minen i-re és ∑_{i=1}^{r} π_{i} = 1.

2.4.7. Definíció. A p_{i}(t) = ℙ(ξ_{t} = i) függvényeket (i = 1, 2, …, r)állapot-valószínűség függvényeknek hívjuk.

Legyen p(t) = (p_{1}(t), …, p_{r}(t)). Ekkor nyilván p(t) minden t-re valószínűségi vektor és a teljes valószínűség tétele szerint

p_{i}(t) = ℙ(ξ_{t} = i) = ∑_{u=1}^{r}P_{ui}(t)p_{i} = p(0) [P (t)]^{(i)} ,

vagyis p(t) = p(0)P (t), amiből p′(t) = p(0)P′(t) = p(0)P (t)Q = p(t)Q, azaz p(t) időbeli fejlődése a

p′(t) = p(t)Q

differenciál-egyenlettel írható le.

2.4.8. Definíció. A p r dimenziós valószínűségi vektor a ξ_{t},t ≥ 0 Markov folyamat stacionárius eloszlása, ha ℙ(ξ_{0} = i) = p_{i} maga után vonja, hogy ℙ(ξ_{t} = i) = p_{i} minden t ≥ 0-re.

Ha létezik p stacionárius eloszlás, akkor arra nyilván p = pP (t) ∀t-re, vagy ezzel ekvivalensen

pQ = 0.  (2.8)

Megfordítva, ha a p valószínűségi vektor kielégíti a 2.8 egyenlőséget, akkor p a ξ_{t},t ≥ 0 Markov folyamat stacionárius eloszlása. Fennáll a következő tétel.

2.4.2. Tétel. Legyen a ξ_{t}, t ≥ 0 Markov folyamat irreducibilis és tegyük fel, hogy véges idő alatt 1 valószínűséggel csak véges sok állapotváltás történik. Ekkor ξ ergodikus akkor és csak akkor, ha létezik egy olyan p valószínűségi vektor, hogy p_{i} > 0 minden i ∈ E-re és p megoldása a

pQ = 0
p1=1

lineáris egyenletrendszernek. Ekkor p egyértelmű és

p_{j} = lim_{t↦∞} P_{ij}(t)

minden i, j ∈ E esetén.

A későbbi vizsgálódásaink során felmerülő véges állapotterű Markov folyamatokról mindig feltesszük majd, hogy irreducibilisek, stacionárius átmenetvalószínűségűek és azt, hogy egyértelműen létezik stacionárius eloszlás. Végezetül megemlítjük, hogy könnyen ellenőrízhetően a 2×2-es esetben a stacionárius eloszlás:

(p_{1}, p_{2}) = (\frac{q_{2}}{q_{1} + q_{2}} , \frac{q_{1}}{q_{1} + q_{2}} ) .

2.5. Kvázi-Markov folyamatok és Markov felújítási sorozatok

Szükségünk lesz a sztochasztikus folyamatok két újabb osztályának ismertetésére, melyekből később speciális esetként adódnak majd a számunkra érdekes Markov-modulált Poisson folyamatra vonatkozó következtetések. E két osztály pedig nevezetesen a kvázi-Markov folyamatok (SMP, az angol "semi-Markov process" elnevezésből) és a Markov felújítási sorozatok (MRS, az angol "Markov renewal sequence" elnevezésből) osztálya, habár tárgyalásunkból következni fog, hogy ez a két osztály tulajdonképpen egy és ugyanaz. Történeti kialakulását tekintve az SMP fogalma született meg előbb, Lévy[12], Smith[13], illetve Takács[27] egymástól függetlenül definiálták 1954ben.

Egy SMP egy olyan sztochasztikus folyamat, ami egy adott időpillanatban az előre rögzített, legfeljebb megszámlálható sok állapot valamelyikében van oly módon, hogy az egymás után meglátogatott állapotok Markov-láncot alkotnak, és egy adott állapotban a folyamat által eltöltött idő egy olyan valószínűségi változó, melynek eloszlása függhet az adott és a legközelebb meglátogatásra kerülő állapottól egyaránt. Ez szemléletesen egy olyan Markov láncként fogható fel, aminek az időparamétere véletlenszerűen transzformált. Mi a továbbiakban csak olyan SMP-kel foglalkozunk majd, melyek csak véges sok állapottal rendelkeznek.

Egy MRS pedig egy olyan valószínűségi változókból álló (X_{n}, Y_{n}) sorozat, ahol X_{n} egy SMP i. állapotváltásának idejét rögzíti, Y_{n} pedig az (i − 1). és az i. állapotváltás között eltelt időt. A szakirodalomban szokás még használni a Markov felújítási folyamat elnevezést, mi itt a sorozat szót használjuk, hangsúlyozva ezzel a folyamat diszkrét voltát. Az alfejezet további részeiben szabatosan is definiáljuk az előbbi szemléletes megfogalmazásokat, majd ismertetjük a terület néhány alapvető eredményét.

2.5.1. Definíció. Egy F : ℝ ↦ ℝ monoton növő, jobbról folytonos függvényt súlyfüggvénynek hívunk, ha F (∞) ≤ 1.

2.5.2. Definíció. Egy Q : ℝ ↦ ℝ^{r×r} mátrix értékű függvényt átmeneteloszlás-mátrixnak hívunk, ha a Q_{ij} elemek olyan súlyfüggvények, melyek kielégítik az alábbi feltételeket:

(i) Q_{ij}(t) = 0, ha t ≤ 0 és

(ii) Q_{i·}(t) = ∑_{j=1}^{r} Q_{ij}(∞) = 1, i = 1, …, r.

2.5.3. Definíció. Azt mondjuk, hogy a kétdimenziós (X_{n}, Y_{n}), n ≥ 0 sztochasztikus folyamat a π kiindulási valószínűségi vektor és a Q átmeneteloszlás-mátrix által meghatározott Markov felújítási sorozat, ha

1. Y_{0} = 0,

2. ℙ(X_{0} = i) = π_{i}, és

3. ℙ(X_{n} = i, Y_{n} ≤ y|X_{0}, X_{1}, Y_{1}, X_{2}, Y_{2}, …, X_{n−1}, Y_{n−1}) = ℙ(X_{n} = i, Y_{n} ≤ y|X_{n−1}) = Q_{X_{n−1},i}(y)

minden y ∈ ℝ és i = 1, …, r esetén.

2.5.1. Megjegyzés. Az így definiált MRS-re röviden a (π, Q) pár segítségével hivatkozunk majd. Itt és a továbbiakban valószínűségi változókra vonatkozó egyenlőségeken mindig majdnem biztosan fennálló egyenlőséget értünk.

Erre a modellre úgy tekinthetünk, mint ami egy olyan rendszer időbeli fejlődését írja le, ami minden t időpillanatban az E = {1, …, r} állapotok valamelyikében van és két állapot között véletlenszerűen (Y_{n} szerint) időzik. Ez a véletlenszerűen eltöltött idő pedig az aktuális és az éppen látogatásra kerülő állapottól függ. Hogy jobban megértsük a modellt, a következő példákban rámutatunk, hogy miként illeszthető be ebbe a tárgyalásba néhány igen jól ismert sztochasztikus folyamat.

2.5.1. Példa. Ha r = 1, akkor (X, Y) egy felújítási folyamatot határoz meg, ha az n. felújítás ideje ∑_{i=0}^{n} Y_{i}.

2.5.2. Példa. Ha Y_{n} = 1 minden n ≥ 1-re, akkor (X, Y) egy diszkrét idejű Markov láncot határoz meg.

2.5.3. Példa. Ha Y_{n} exponenciális eloszlású olyan paraméterrel, ami X_{n−1}-től függ, akkor egy folytonos idejű Markov láncot kapunk.

Vezessük be a T_{n} = ∑_{i=0}^{n} Y_{i} jelölést, ekkor T_{n} az n-edik állapotváltás idejeként interpretálható.

2.5.2. Megjegyzés. Megjegyezzük, hogy itt az "állapotváltás" szót használtuk, ami azonban nem feltétlenül jelent a szó szoros értelmében vett állapotváltást, tekintve, hogy egyelőre nem zártuk ki a Q_{ii} > 0 lehetőséget.

2.5.1. Lemma. (X_{n}, T_{n}), n ≥ 0 kielégíti az alábbi egyenlőségeket:

ℙ(X_{n} = i, T_{n} ≤ t|X_{0}, X_{1}, T_{1}, X_{2}, T_{2}, …, X_{n−1}, T_{n−1}) = Q_{X_{n−1},i}(t − T_{n−1})  (2.9)

és

ℙ(X_{n} = i|X_{0}, X_{1}, T_{1}, X_{2}, T_{2}, …, X_{n−1}, T_{n−1}) = Q_{X_{n−1},i}(∞).  (2.10)

Bizonyítás. (2.9) könnyen látható a

ℙ(X_{n} = i, Y_{n} ≤ t − T_{n−1}|X_{0}, X_{1}, T_{1}, X_{2}, T_{2}, …, X_{n−1}, T_{n−1})

átrendezett alakból és abból a tényből, hogy a σ-algebra, amelyre a feltételt vesszük megegyezik a 2.5.3 definícióban szereplő σ-algebrával. (2.10) pedig a Lebesgue monoton konvergencia tétel egyszerű következménye.

Vezessünk még be néhány jelölést. Legyen Q̃ = Q(∞). A 2.5.2 definíció alapján világos, hogy Q̃ sztochasztikus mátrix. Jelölje U_{c}(x) a c pontra koncentrált valószínűségi változó eloszlásfüggvényét. Definiáljuk az F mátrixértékű függvényt az alábbi képlettel:

F_{ij} := \{ \begin{matrix}\frac{Q_{ij}}{Q̃_{ij}} & , ha & Q̃_{ij} > 0 \\ U_{1} & , ha & Q̃_{ij} = 0 \end{matrix} .

Legyen továbbá F = (F_{ij}). Ezek segítségével könnyen beláthatók az alábbi állítások.

2.5.2. Lemma. A fenti jelölések mellett fennállnak az alábbi összefüggések

1. ℙ(Y_{n} ≤ y|X_{0}, X_{1}, …, X_{n−1}) = Q_{X_{n−1},·}(y)

2. ℙ(Y_{n} ≤ y|X_{0}, X_{1}, …, X_{n}) = F_{X_{n−1},X_{n}}(y)

3. ℙ(X_{n} = i|X_{0}, X_{1}, …, X_{n−1}) = Q̃_{X_{n−1},i} 4. ℙ(Y_{n_{1}} ≤ y_{1}, Y_{n_{2}} ≤ y_{2}, …, Y_{n_{k}} ≤ y_{k}|X_{n}, n ≥ 0)

= ℙ(Y_{n_{1}} ≤ y_{1}, Y_{n_{2}} ≤ y_{2}, …, Y_{n_{k}} ≤ y_{k}|X_{n_{1}}, X_{n_{k}}, …, X_{n_{k}}) = ∏_{i=1}^{k}F_{X_{n_{i}−1},X_{n_{i}}} (y_{i}).

Vegyük észre, hogy a 3. állítás éppen azt jelenti, hogy X_{0}, X_{1}, … Markov lánc, melynek átmenetvalószínűség mátrixa Q̃.

Értelmezzük az N_{t}, t ≥ 0 egész értékű sztochasztikus folyamatot az N_{t} = sup{n ≥ 0 : T_{n} ≤ t} képlettel és. Ezenkívül legyen i = 1, …, r esetén N_{t}^{i} = #{k : X_{k} = i, k = 1, … N_{t}}, vagyis N_{t}^{i} azt adja meg, hogy a (0, t] intervallumban N_{t} hányszor került az i állapotba. Világos, hogy minden t-re N_{t} = ∑_{i=1}^{r} N_{t}^{i}.

2.5.4. Definíció. Az N_{t}, t ≥ 0 sztochasztikus folyamatot az MRP számlálófolyamatának, míg az N_{t}^{i}, t ≥ 0 sztochasztikus folyamatokat az MRP-hez tartozó i. típusú felújítási folyamatnak nevezzük.

Könnyen ellenőrizhető, hogy az elnevezés jogos, vagyis N_{t}^{i} valóban felújítási folyamat. Az iménti definíciók figyelembe vételével bevezethetünk egy másik sztochasztikus folyamatot, amely minden időpillanatban az MRP aktuális állapotának az értékét veszi fel.

2.5.5. Definíció. Azt mondjuk, hogy a Z_{t}, t ≥ 0 folyamat a π k.v.v. és a Q átmeneteloszlás-mátrix által meghatározott kvázi-Markov folyamat, ha Z_{t} = X_{N_{t}}.

A 2.2 abrán látható egy kvázi-Markov folyamat lehetséges realizációja.

12345 0=T_{0} T_{1}  T_{2}  T_{3} T_{4} T_{5}  T_{6}

2.2. ábra. Egy SMP lehetséges realizációja

Bizonyítás nélkül közöljük az alábbi tételeket, a részletek megtalálhatók R. Pyke [18] illetve [19] cikkében.

2.5.1. Tétel. ℙ(N_{t} < ∞, ∀t ≥ 0) = 1.

Legyen most az MRS véges állapotterű, ekkor kimondható az alábbi tétel.

2.5.2. Tétel. Definiáljuk minden t ≥ 0 és n ≥ 0 esetén a

P_{ij}(n, t) = ℙ(N_{t} = n, Z_{t} = j|Z_{0} = i)

mennyiségeket. Ekkor

P_{ij}(0, t) = δ_{ij}U_{0}(t)(1 − Q_{i·}(t))

és n = 1, 2, … esetén

P_{ij}(n, t) = ∑_{k=1}^{r} P_{kj}(n − 1, t) ∗ Q_{ik}(t).

2.5.3. Megjegyzés. Itt ∗ konvolúciót jelöl, azaz

(f ∗ g) (t) = \{ \begin{matrix}∫_{0−}^{t} f(t − y)dg(y) & , ha & t ≥ 0 \\ 0 & , ha & t < 0 \end{matrix} .

Ezután a konvolúció művelete mátrixértékű függvényekre is definiálható:

K ∗ L = ∑_{k=1}^{r} K_{ik} ∗ L_{kj},

mely formálisan éppen ugyanúgy van definiálva, mint a hagyományos mátrixszorzás azzal a különbséggel, hogy itt a szorzás operátor konvolúcióra van cserélve. Legyen továbbá K^{(∗0)} = I és n = 1, 2, … esetén definiáljuk a K^{(∗n)} mátrixértékű függvényt a K mátrixértékű függvény önmagával vett n-szeres konvolúciójaként, azaz legyen K^{(∗n)} = K^{(∗(n−1))} ∗ K.

A 2.5.2 tétel így fogalmazható meg tömören mátrix jelölésekkel:

P (0, t) = I − H(t)  (2.11)

és

P (n, ·) = ∑_{n=0}^{∞} Q^{(∗n)} ∗ (I − H),  (2.12)

ahol I_{ij} = δ_{ij}U_{0}(t) és H_{ij}(t) = δ_{ij}Q_{i·}(t).

2.6. A Markov-modulált Poisson folyamat, mint speciális Cox

folyamat

A Markov-modulált Poisson folyamat, vagy röviden MMPP(az angol Markov-modulated Poisson Process elnevezés alapján) egy olyan általános Poisson folyamatként képzelhető el, amelyben a μ(t) intenzitásfüggvény nem determinisztikus, hanem maga is egy sztochasztikus folyamat, mégpedig egy véges állapotterű folytonos idejű Markov-folyamat által meghatározott az alábbiak szerint: μ(t) = λ_{i}, ha a Markov-folyamat éppen az i-edik állapotban van.

Hogy precízzé tegyük definíciónkat, legyen ξ_{t} egy véges állapotterű folytonos idejű Markov folyamat, melynek állapottere az egyszerűség kedvéért legyen E = {1, 2, … r}. Tekintsük továbbá azt a λ : E ↦ ℝ^{+} fügvényt, amelyre λ(i) = λ_{i} minden i ∈ E esetén az előre rögzített λ_{1}, λ_{2}, …, λ_{r} ∈ ℝ^{+} számokkal. Az MMPP-t egy olyan Cox folyamatként fogjuk definiálni, melynek η intenzitásfolyamatára η_{t} = λ(ξ_{t}), t > 0. Legyen tehát G = {(s, t] : 0 < s ≤ t}, majd definiáljuk a ρ véletlen mértéket a következőképpen: és ρ((s, t]) := ∫_{s+}^{t} η_{t}dt. Ezzel egy egészen σ-véges mértéket értelmeztünk, mely egyértelműen terjed ki a generált σ-algebrára, vagyis B_{0}-ra. Ha most a kiterjesztett mértéket továbbra is ρ jelöli, úgy ezzel megadtunk egy B_{0}-n értelmezett egészen σ-véges mértéket. Ezzel már definiálható a Markov-modulált Poisson folyamat fogalma.

2.6.1. Definíció. Legyen a ρ véletlen mérték az előző bekezdésben leírtak szerint definiált, ekkor az általa meghatározott Cox folyamatot Markov-modulált Poisson folyamatnak (MMPP) hívjuk.

Az előző alfejezetekre visszatekintve az adódik, hogy egy MMPP paraméterezhető a mögöttes Markov-folyamat Q infinitezimális generátormátrixával, illetve a λ = [λ_{1}, λ_{2}, …, λ_{r}]^{T} vektorral. Ez összesen r(r − 1) + r = r^{2} paramétert jelent. Később ezen paraméterek becslésére mutatunk majd be egy eljárást. Szokás még bevezetni a Λ = diag(λ) jelölést, mellyel a paramétertér az alábbi formában írható le:

Θ = {(Q, Λ) : Q irreducibilis és Λ = diag(λ), ahol λ_{i} > 0∀i}.

Igazából a mögöttes Markov folyamat kezdeti eloszlása is paraméter, ezt azonban nem fogjuk becsülni, hanem Q függvényeként választjuk majd meg. A Q és Λ paramétereket illetően fenáll az alábbi lemma.

2.6.1. Lemma. Amennyiben Q irreducibilis, úgy a Λ = 0 triviális esetet leszámítva Q−Λ nemszinguláris.

Ahhoz, hogy egy MMPP-t teljesen leírjunk, szükséges, hogy megadjuk a mögöttes Markov folyamat t = 0 pillanatbeli állapotát. Jelöljük ezt π-vel. Az alábbiakban ismertetjük a (Q, Λ, π) paraméterekkel rendelkező MMPP fejlődését, ahol tehát π és Q rendre a mögöttes Markov folyamat kezdeti eloszlásának vektora illetve infinitezimális generátora. Ehhez előbb belátjuk az alábbi lemmát.

2.6.2. Lemma. Legyenek ξ_{1} és ξ_{2} független exponenciális eloszlású valószínűségi változók λ_{1} és λ_{2} paraméterekkel. Legyen továbbá

U = min(ξ_{1}, ξ_{2})

és

N = \{ \begin{matrix}1 & , ha & ξ_{1} < ξ_{2} \\ 2 & , ha & ξ_{1} ≥ ξ_{2} \end{matrix} .

Ekkor fenáll a következő három állítás.

2.3. ábra. Egy MMPP lehetséges realizációja

1. ℙ(N = 1) = \frac{λ_{1}}{λ_{1}+λ_{2}} és ℙ(N = 2) = \frac{λ_{2}}{λ_{1}+λ_{2}} .

2. U exponenciális eloszlású λ_{1} + λ_{2} paraméterrel.

3. N és U függetlenek. Bizonyítás.

1. Mivel ξ_{1} és ξ_{2} függetlenek, azért ξ_{1} − ξ_{2} eloszlásfüggvénye felírható ξ_{1} és −ξ_{2} eloszlásfüggvé-

nyeinek konvolúciójaként:

ℙ(ξ_{1} < ξ_{2}) = F_{ξ_{1}−ξ_{2}}(0) = F_{ξ_{1}} ∗ F_{−ξ_{2}}(0) = ∫_{−∞}^{∞} F_{ξ_{1}}(−y)dF_{−ξ_{2}}(y)

= ∫_{−∞}^{0} F_{ξ_{1}}(−y)dF_{−ξ_{2}}(y) = ∫_{−∞}^{0}(1 − e^{−λ_{1}y})d(e^{λ_{2}y})

= ∫_{−∞}^{0}(1 − e^{−λ_{1}y})λ_{2}e^{λ_{2}y}dy = 1 − λ_{2} \frac{1}{λ_{1} + λ_{2}} = \frac{λ_{1}}{λ_{1} + λ_{2}} .

2. Ismét ξ_{1} és ξ_{2} függetlenségéből adódóan ℙ(U ≥ x) = ℙ(ξ_{1} ≥ x, ξ_{2} ≥ x) = ℙ(ξ_{1} ≥ x)ℙ(ξ_{2} ≥ x) = e^{−λ_{1}x}e^{−λ_{2}x} = e^{−(λ_{1}+λ_{2})x}, ha x ≥ 0.

3. Az f_{Y} (y) = f_{h(Y)}(h(y))|J (h(y))| összefüggés alkalmazásával — ahol |J (h(·))| jelöli a h függvény Jacobi-determinánsát — és ξ_{1} illetve ξ_{2} függetlenségének kihasználásával felírható, hogy

f_{(ξ_{1}−ξ_{2},ξ_{2})}(s, t) = f_{(ξ_{1},ξ_{2})}(s + t, t) = f_{ξ_{1}}(s + t)f_{ξ_{2}}(t).

Emiatt

ℙ(N = 2, U ≥ x) = ℙ(ξ_{1} − ξ_{2} ≥ 0, ξ_{2} ≥ x) = ∫_{x}^{∞}∫_{0}^{∞}f_{(ξ_{1}−ξ_{2},ξ_{2})}(s, t)dsdt

= ∫_{x}^{∞}∫_{0}^{∞}f_{ξ_{1}}(s + t)f_{ξ_{2}}(t)dsdt = ∫_{x}^{∞}f_{ξ_{2}}(t) ∫_{0}^{∞}f_{ξ_{1}}(s + t)dsdt

= ∫_{x}^{∞}f_{ξ_{2}}(t) ∫_{t}^{∞}f_{ξ_{1}}(u)dudt = ∫_{x}^{∞}f_{ξ_{2}}(t)e^{−λ_{1}t}dt = λ_{2} ∫_{x}^{∞}e^{−(λ_{1}+λ_{2})t}dt

= \frac{λ_{2}}{λ_{1} + λ_{2}} e^{−(λ_{1}+λ_{2})x} = ℙ(N = 2)ℙ(U ≥ x),

ami bizonyítja a függetlenséget.

Ezek után egyaránt nevezzük eseménynek azt, ha ξ_{t} állapotot vált t-ben vagy ha t érkezési ideje az MMPP-nek. Ekkor az MMPP fejlődése az alábbiak szerint jellemezhető. Ha s-ben bekövetkezik egy esemény, akkor feltéve, hogy ξ_{s+} = i, a következő eseményig eltelt idő exponenciális eloszlású q_{i} + λ_{i} paraméterrel és az eltelt időtől függetlenül az esemény \frac{q_{ij}}{q_{j}+λ_{i}} valószínűséggel állapotváltás vagy \frac{λ_{i}}{q_{j}+λ_{i}} valószínűséggel ugrás. Ez azt jelenti, hogy (N_{t}, ξ_{t}) egy megszámlálható állapotterű Markov folyamat, ezért a

P_{ij}(n, t) = ℙ(N_{t} = n, ξ_{t} = j|ξ_{0} = i)  (2.13)

átmenet-valószínűségekre n ≥ 1 esetén felírhatók az alábbi előre haladó Chapman-Kolmogorov egyenletek:

P′_{ij}(n, t) = _{∑l≠j}P_{il}(n, t)q_{lj} + P_{ij}(n, t)(q_{j} − λ_{j}) + P_{ij}(n − 1, t)λ_{j},

ami mátrix jelöléssel így írható tömören:

P′(n, t) = P (n, t)(Q − Λ) + P (n − 1, t)Λ,

ahol nyilván P (0, 0) = I és P (n, 0) = I, ha n ≥ 1. Hasonlóan n = 0-ra

P′(n, t) = P (n, t)(Q − Λ)

adódik. Ekkor a P^{∗}(z, t) = ∑_{n=0}^{∞} P (n, t)z^{n} generátorfüggvény kielégíti a

\frac{∂}{∂t}P^{∗}(z, t) = P^{∗}(z, t)(Q − Λ) + zP^{∗}(z, t)Λ

egyenlőséget, továbbá P^{∗}(z, 0) = I, amiből már következik, hogy P^{∗}(z, t) = exp [(Q − Λ + zΛ)t]. Ennek segítségével meghatározzuk N_{t} várható értékét minden rögzített t-re. Ha az M(t) mátrix az

M(t) = [\frac{∂}{∂z} P^{∗}(z, t)]_{z=1}

képlettel definiált, akkor világos, hogy

M_{ij}(t) = ∑_{n=0}^{∞} nP_{ij}(n, t),

ami azt jelenti, hogy

E(N_{t}) = ∑_{n=0}^{∞} nP (N_{t} = n) = ∑_{n=0}^{∞} ∑_{j=1}^{r} nP (N_{t} = n, ξ_{t} = j)
= ∑_{n=0}^{∞} ∑_{j=1}^{r} ∑_{i=1}^{r} nP_{ij}(n, t)ℙ(ξ_{0} = i)
= ∑_{j=1}^{r} ∑_{i=1}^{r} M_{ij}(t)π_{i} = πM(t)1.

A mátrix differenciálás szabályait figyelembe véve

M(t) = ∑_{n=1}^{∞} \frac{t^{n}}{n!} ∑_{p=0}^{n−1} Q^{p}ΛQ^{n−1−p}.

Kis számolással az is belátható, hogy ekkor

M(t)1 = 1πλt+( exp(Qt) − I)(Q + 1π)^{−1}λ,

amiből azt kapjuk, hogy

E(N_{t}) = π1πλt+π( exp(Qt) − I)(Q + 1π)^{−1}λ.

Ha most π speciálisan a ξ_{t} Markov folyamat stacionárius eloszlása, akkor π exp(Qt) = π és így

E(N_{t}) = π1πλt = πλt.

2.6.1. Megjegyzés. A P (n, t)-re levezetett összefüggések kiszámolásához felhasználhattuk volna a 2.11 és 2.12 alatti összefüggéseket is.

Tegyük fel most, hogy T rögzített és hogy a T időpontig megfigyelt legutolsó állapotváltás a s < T időpontban történik, akkor az exponenciális eloszlás örökifjú tulajdonsága miatt a T utáni legelső esemény q_{i} + λ_{i} paraméterű exponenciális eloszlás szerint következik be, feltéve, hogy ξ_{t} T -ben az i. állapotban volt. Ez azt jelenti, hogy az N_{t}^{∗} = N_{T −t} − N_{T} összefüggéssel egy MMPP-t értelmezünk, mely esetén ξ_{t}^{∗} = ξ_{T +t} a mögöttes Markov folyamat a π_{i}^{∗} = ℙ(ξ_{0}^{∗} = i) = ℙ(ξ_{T} = i) kezdeti valószínűségekkel. Tehát míg N_{t}^{∗} csupán a π^{∗} kezdeti eloszlásvektorban különbözik N_{t}-től. Ha π a Q stacionárius eloszlása, akkor még π^{∗} = π is igaz. Ez alapján az 2.13 alatti mennyiségek segítségével meghatározhatjuk a ℙ(N_{T +t} − N_{T} = n), n ≥ 0 valószínűségeket.

ℙ(N_{T +t} − N_{T} = n) = ℙ(N_{t}^{∗} = n) = ∑_{j=1}^{r} ℙ(N_{t}^{∗} = n, ξ_{t}^{∗} = j)
= ∑_{i=1}^{r} ∑_{j=1}^{r} ℙ(N_{t}^{∗} = n, ξ_{t}^{∗} = j|ξ_{0}^{∗} = i)ℙ(ξ_{0}^{∗} = i)
= ∑_{i=1}^{r} ∑_{j=1}^{r} P_{ij}(n, t)ℙ(ξ_{T} = i) = π^{∗}P (n, t)1.

A Markov-folyamatokról szóló fejezetben láttuk, hogy π^{∗} = π exp(Qt), tehát

ℙ(N_{T +t} − N_{T} = n) = π exp(Qt)P (n, t)1,

ha pedig π a mögöttes Markov folyamat stacionárius eloszlása, akkor

ℙ(N_{T +t} − N_{T} = n) = πP (n, t)1.

A következőkben bemutatjuk, hogy miként gondulhatunk az MMPP-re úgy, mint egy MRS-re. Ennek érdekében legyen X_{n}, n ≥ 0 a mögöttes Markov folyamat állapota az n. megfigyeléskor (X_{0} legyen a t = 0 időpontbeli állapot). Legyen továbbá Y_{0} = 0 és n ≥ 1 esetén jelölje Y_{n} az n. és az (n − 1). érkezés között eltelt időt. Ekkor (X_{n}, Y_{n}), n ≥ 0 sorozat, ennek megokolása megtalálható például [6]-ben. Fenáll továbbá az alábbi tétel, melynek bizonyítása megtalálható [17]-ben ezért annak közlését most mellőzzük.

2.6.1. Tétel. Az imént definiált (X_{n}, Y_{n}), n ≥ 0 Markov felújítási sorozat átmeneteloszlás-mátrixa:

Q(x) = ∫_{0}^{x} exp[(Q − Λ)u]duΛ
= {I − exp[(Q − Λ)x]} (Λ − Q)^{−1}Λ
= {I − exp[(Q − Λ)x]} Q(∞).

2.6.2. Megjegyzés. Ha az időhorizont kezdete nem felel meg egy érkezési időpontnak, akkor a ℙ(ξ_{1} = j, Y_{1} ≤ y|ξ_{0} = i) valószínűségeket külön definiálni kell.

A Q̃ = Q(∞) = (Λ−Q)^{−1}Λ mátrix sztochasztikus, ahogy azt láttuk az 5. alfejezetben. Ráadásul Q̃ azon Markov lánc átmenet-valószínűség mátrixa, ami rögzíti a folyamat állapotát az i. érkezéskor (beágyazott Markov lánc). Ha pedig minden λ_{i} pozitív - ahogy azt mi végig fel is tesszük az általunk vizsgált MMPP-k esetén - akkor Q̃ egyben irreducibilis is. Ha π jelöli a Q mátrixhoz tartozó stacionárius eloszlást, akkor könnyű látni, hogy a

π̃ = \frac{πΛ}{πλ}

valószínűségi vektor kielégíti a π̃ = π̃Q̃ lineáris egyenletrendszert.

Térjünk rá a mögöttes Markov folyamat kiindulási vektorának megválasztásának kérdésére. Az MMPP intervallum-stacionárius verzióját kapjuk, ha a kiindulási valószínűségi vektort π̃-nak választjuk. Az ún. környezet-stacionárius verzióhoz jutunk, ha a kiindulási valószínűségi vektort a Q mátrix π stacionárius vektorának választjuk. Ez a Markov felújítási sorozatként történő leírásban azt jelenti, hogy a ℙ(X_{0} = i) valószínűségeket π_{i}-re állítjuk és aztán a ℙ(ξ_{1} = j, Y_{1} ≤ y|ξ_{0} = i) valószínűségeket úgy definiáljuk, hogy a belőlük alkotott mátrix értékű függvény az y helyen éppen Q(y)-nal legyen egyenlő. Ekkor az időhorizont kezdete ugyan nem egy érkezési időpont, a mögöttes ξ Markov folyamat azonban stacionárius lesz. Mi ez utóbbi verziót fogjuk majd a későbbiekben használni.

Végezetül kimondunk egy tételt, mely a későbbiekben fontos szerepet játszik majd. A bizonyítását illetően lásd [31] vagy [17].

2.6.2. Tétel. Tegyük fel, hogy adott egy (Y_{1}, …, Y_{n}) minta egy ϑ = (Q, Λ) paraméterű MMPP-ből származó köztes időkből, ekkor (Y_{1}, …, Y_{n}) együttes sűrűségfüggvénye

f_{(Y_{1},...,Y_{n})}^{ϑ}(y_{1}, …, y_{n}) = π^{ϑ} (∏_{k=1}^{n} f^{ϑ}(y_{k})) 1,

ahol

f^{ϑ}(t) = exp [(Q − Λ)t] Λ.  (2.14)

Az MMPP paramétereinek becslése

az EM algoritmus segítségével

A fejezet első alfejezete egy rövid összefoglaló az EM algoritmusról, mint paraméterbecslési eljárásról. A 2. alfejezet a Rydén [23] által prezentált EM algoritmust hivatott bemutatni, míg a 3. alfejezet az algoritmus implementálásával foglalkozik, figyelembe véve a W.J.J. Roberts et. al. által közölt új eredményket. Az alfejezet saját eredménye a kiinduló becslés meghatározására szolgáló k-means klaszterezést használó eljárás.

3.1. Az EM algoritmus

Az úgynevezett EM (expectation maximization) algoritmus egy széles körben használható iteratív paraméterbecslési eljárás olyan feladatok megoldására, ahol a minta egy része nem megfigyelhető. Ez a szituáció lép fel például akkor, ha (X, Y) egy véletlen vektor az f_{ϑ} együttes sűrűségfüggvénnyel és a ϑ ∈ θ paramétert szeretnénk becsülni, viszont csak az Y valószínűségi változót tudjuk megfigyelni. Az alábbiakban megmutatjuk, hogy hogyan kezelhetők az ilyen esetek az EM algoritmus segítségével.

3.1.1. Definíció. Legyen (Ω, A, ℙ_{ϑ}), (ϑ ∈ θ) statisztikai mező, (X, B) mérhető tér és X : Ω ↦ X valószínűségi változó. Tegyük fel, hogy a statisztikai mező dominált, azaz létezik Λ|B mérték, hogy ℙ_{ϑ} ∘ X^{−1} << Λ minden ϑ ∈ θ-ra. Legyen x ∈ X egy megfigyelés az X valószínűségi változóból. Ekkor az x megfigyelés loglikelihood-függvénye alatt a ϑ ↦ log f_{ϑ}(x) függvényt értjük, ahol f_{ϑ} = \frac{dΛ}{d(ℙ_{ϑ}∘X^{−1})}.

3.1.2. Definíció. Legyen x ∈ X egy megfigyelés az X valószínűségi változóból. Az előző definíció jelölései mellett azt mondjuk, hogy ϑ̂ a ϑ ∈ θ paraméter maximum-likelihood (ML) becslése, ha maximalizálja a ϑ ↦ log f_{ϑ}(x) loglikelihood-függvényt.

Ha most (X, Y) egy véletlen vektor és ϑ ∈ θ a becsülendő paraméter, de csak Y -t tudjuk megfigyelni, akkor nem tudjuk kiszámolni az ML becslést, mivel a nem megfigyelhető adat nem áll rendelkezésünkre. Ehelyett keresünk egy olyan függvényt, ami bizonyos értelemben jól közelíti log f_{ϑ}(X, Y)-t, de csak Y -tól függ. Ehhez tekintsük az alábbi lemmát.

3.1.1. Lemma. Legyen (Ω, A, ℙ) egy valószínűségi mező és tekintsük az L^{2}(Ω)-beli valós értékű X illetve Y valószínűségi változókat. Definiáljuk az F függvényhalmazt az alábbi módon.

F = {f : ℝ ↦ ℝ mérhető függvény, 𝔼f(Y)^{2} < ∞}

Ekkor

min_{f∈F} ∥X − f(Y)∥_{2} = ∥X − 𝔼(X|Y)∥_{2} .

- ahol ∥·∥_{2} L^{2}-normát jelöl, - és a minimum a f(y) = 𝔼(X|Y = y) függvény esetén vétetik fel.

Ha tehát abból indulunk ki, hogy már van egy θ_{0} feltételezésem arról, hogy mi lehet a valódi paraméter, akkor a log f_{ϑ}(X, Y) függvény legjobb közelítése a fenti lemma szerinti értelemben l(ϑ_{0}, ϑ, y) = 𝔼_{ϑ_{0}}( log f_{ϑ}(X, Y)|Y = y), ami már csak a y megfigyeléstól függ, így ennek már meg tudjuk határozni a ϑ-ban vett maximumát. Ha ez a maximum a ϑ_{1}-ben vétetik fel, akkor azt reméljük, ϑ_{1} már közelebb áll az ML becsléshez, mint ϑ_{0} és ϑ_{1}-re megismételhetjük az iménti eljárást. Tehát az algoritmus az alábbi módon implementálható, az egyszerűség kedvéért megállási kritérium nélkül.

3.1.1. Algoritmus (EM algoritmus). 1. bemenő adatok:ϑ_{0} kiindulási becslés, y megfigyelés a Y valószínűségi változóból.

2. n. iteráció:

a) Határozd meg az l(ϑ_{n−1}, ϑ, y) = 𝔼_{ϑ_{n−1}}( log f_{ϑ}(X, Y)|Y = y) függvényt! (ún. E-lépés)

b) ϑ_{n} := arg max_{ϑ} l(ϑ_{n−1}, ϑ, y). (ún. M-lépés).

3. kimenő adat: {ϑ_{n}} végtelen sorozat.

A gyakorlatban általában az E és M lépések közül az egyik nehezen kivitelezhető, míg a másik könnyen, hogy mikor melyik, az az adott probléma sajátosságaitól függ. Megmutatható, hogy l(ϑ_{0}, y) ≤ l(ϑ_{1}, y) ≤ …, azaz, hogy minden egyes iterációs lépés növeli az Y -ból vett y minta loglikelihood-függvényét. Az is igazolható továbbá, hogy bizonyos meghatározott regularitási feltételek mellett (lásd Wu [32]) ϑ_{n} a ϑ ↦ l(ϑ, y)-ból számított ML-becsléshez konvergál.

3.1.1. Megjegyzés. Bár jelen dolgozatban nem kerül tárgyalásra, csupán nagy vonalakban megemlítjük, hogy létezik az EM algoritmusnak egy általánosabb változata az ún. GEM(generalized EM). A GEM esetében nem követeljük meg az EM iterációban leírt M lépés teljesülését, csupán az a törekvésünk, hogy

l(ϑ_{n}, y) ≤ l(ϑ_{n+1}, y)

teljesüljön. Világos, hogy minden EM algoritmus egyben GEM is.

3.2. Egy EM algoritmus Markov-modulált Poisson folyama-

tokra

A következő fejezetben ismertetünk egy EM algoritmust, melyet először T. Rydén [23] mutatott be 1996-os cikkében. Tételezzük fel, hogy van egy n elemű megfigyelésünk az y_{1}, y_{2}, …, y_{n} köztes időkből. Ez nyilván nem teljes adat, hiszen a mögöttes Markov folyamatról nem kapunk kielégítő információt. Ebben a helyzetben tehát természetesen adódó lehetőségként merül fel az EM algoritmus alkalmazása, a kérdés csupán az, hogy mit válasszunk nem megfigyelhető adatnak. Itt elsősorban két lehetőség kínálja magát:

1. Választható hiányzó adatnak az {X_{k}} beágyazott Markov-lánc, vagy

2. a mögöttes {X_{t}} Markov-folyamat egész trajektóriája.

Rydén korábbi munkájában vizsgálta az első esetet és megmutatta, hogy bár ekkor at E-lépés egyszerű, az M-lépés ezzel ellentétben nem explicit és csak numerikus módszerek alkalmazása segítségével kivitelezhető. Majd 1996-os cikkében vizsgálta a második esetet, ahol éppen az E-lépés a bonyolultabb, az M-lépés viszont - ahogy azt nemsokára látni fogjuk - explicit módon kifejezhető. Jelen fejezetben ezt az algoritmust fogjuk ismertetni. Ekkor tehát a megfigyelésünk egy Y = (Y_{1}, Y_{2}, …, Y_{n}) véletlen elem az ((ℝ^{+})^{n} , B) mérhető térből (ahol B a Borel-halmazok σalgebrája). A hiányzó (vagy nem megfigyelhető) adat pedig az X Markov folyamat, ami szintén egy véletlen elem egy alkalmasan választott mérhető térben. Először vezessünk be néhány statisztikát és jelölést az Y megfigyelés függvényében.

1. t_{0} := 0 és k = 1, 2, …, n esetén t_{k} := ∑_{i=1}^{k} y_{i}, azaz t_{k} nem más, mint az i-edik érkezés időpontja, a t_{0} = 0 konvencióval. Legyen továbbá T = t_{n}.

2. Tegyük fel, hogy a mögöttes Markov folyamat ugrásai a 0 < u_{1} < u_{2} < … < u_{m} < T időpontokban történnek, majd itt is állapodjunk meg az u_{0} = 0 és u_{m+1} = T jelölésekben.

3. 1 ≤ k ≤ m+1-re I_{k} := [u_{k−1}, u_{k}), azaz I_{k} a mögöttes Markov folyamat (k −1)-edik és k-adik ugrása közti időintervallumot jelöli.

4. Legyen továbbá s_{k} = X_{u} (ahol u ∈ I_{k} tetszőleges), és z_{k} = #{k : t_{k} ∈ I_{k}, (k ≠ 0)}, azaz s_{k} a Markov-folyamat állapota a (k − 1)-edik és k-adik ugrása között, míg z_{k} az érkezések száma a Markov-folyamat (k − 1)-edik és k-adik ugrása között.

5. A t ≤ T ideig bekövetkezett érkezések számát pedig jelöljük N_{t}-vel:

N_{t} := #{k : 0 < k ≤ n, t_{k} ≤ t}.

3.2.1. Megjegyzés. Világos, hogy (y_{1}, …, y_{n})-ből megkonstruálható N_{t} és ugyanez fordítva is igaz.

Ezek után - ahogy az [23]-ban bemutatásra kerül - az (Y, X) minta szerinti teljes L_{C} likelihoodfüggvény a következőképpen írható fel:

L_{C}(Q, Λ) = π_{s_{1}} ∏_{k=1}^{m} (q_{s_{k}}e^{−q_{s_{k}} Δu_{k}} \frac{q_{s_{k},s_{k+1}}}{q_{s_{k}}} ) e^{−q_{s_{m+1}} Δu_{m+1}}

× ∏_{k=1}^{m+1} (\frac{(λ_{s_{k}}Δu_{k})^{z_{k}}}{z_{k}!} e^{−λ_{s_{k}} Δu_{k}}) \frac{z_{k}!}{(Δu_{k})^{z_{k}}}

= π_{s_{1}} ∏_{k=1}^{m} (e^{−q_{s_{k}} Δu_{k}}q_{s_{k},s_{k+1}}) e^{−q_{s_{m+1}} Δu_{m+1}} ∏_{k=1}^{m+1} (λ_{s_{k}}^{z_{k}}e^{−λ_{s_{k}} Δu_{k}}) .

Ebből a loglikelihood-függvény:

l_{C}(Q, Λ) = log π_{s_{1}} − ∑_{k=1}^{m}(q_{s_{k}}Δu_{k} + log q_{s_{k},s_{k+1}}) − q_{s_{m+1}}Δu_{m+1}  (3.1)

+ ∑_{k=1}^{m+1}(z_{k} log λ_{s_{k}} − λ_{s_{k}}Δu_{k}).

Világos, hogy log π_{s_{1}} = ∑_{i=1}^{r} χ(X_{0} = i) log π_{i}. Határozzuk meg most q_{i} együtthatóját! Könnyen láthatóan ez éppen ∑_{{k:1≤k≤m+1,s_{k}=i}} Δu_{k}. Erre bevezethetjük a D_{i} jelölést. D_{i} tehát nem más, mint az az idő, amit az {X_{t}} Markov folyamat az i állapotban tölt a [0, T] intervallum alatt. Észrevetjük, hogy D_{i} felírható az alábbi, kissé kompaktabb formában is:

D_{i} = ∫_{0}^{T} χ(X_{t} = i).

Nyilván λ_{i} együtthatója is ugyanez lesz. A log q_{ij}, (i ≠ j) tag együtthatója pedig így fejezhető ki: #{k : 1 ≤ k ≤ m, s_{k} = i, s_{k+1} = j}. Ez tehát nem más, mint az {X_{t}} Markov folyamat (i, j)- állapotváltásainak a száma a [0, T] intervallum alatt. Erre bevezetjük az m_{ij} jelölést. Másképpen kifejezve:

m_{ij} = #{t : 0 < t ≤ T, X_{t−} = i, X_{t} = j}.

Hátravan még log λ_{i} együtthatójának a meghatározása. Ha ezt n_{i} jelöli, akkor

n_{i} = ∑_{{k:1≤k≤m+1,s_{k}=i}} z_{k} = ∫_{0}^{T} χ(X_{t_{k}} = i).

Tehát n_{i} azon események száma, melyek az alatt az idő alatt következtek be, míg {X_{t}} az i állapotban volt. Ezekkel tehát (3.1) az alábbi alakban írható fel:

l_{C}(Q, Λ) = ∑_{i=1}^{r} χ(X_{0} = i) log π_{i} − ∑_{i=1}^{r} D_{i}q_{i}

+ ∑_{i=1}^{r} ∑_{j=1, j≠i}^{r} m_{ij} log q_{ij} + ∑_{i=1}^{r}(n_{i} log λ_{i} − λ_{i}q_{i}).

Tegyük fel most, hogy ϑ_{0} rögzített és legyen ℙ_{0} = ℙ_{ϑ_{0}}. Ekkor

𝔼_{0}(l_{C}(Q, Λ)|Y = y) = ∑_{i=1}^{r} χ(X_{0} = i) log π̂_{i} − ∑_{i=1}^{r} D̂_{i}q_{i}  (3.2)

+ ∑_{i=1}^{r} ∑_{j=1, j≠i}^{r} m̂_{ij} log q_{ij} + ∑_{i=1}^{r}(n̂_{i} log λ_{i} − λ_{i}q_{i}),

ahol

1. π̂_{i} = ℙ_{0}(X_{0} = i|Y = y) , 2. m̂_{ij} = 𝔼_{0}(m_{ij}|Y = y),

3. D̂_{i} = 𝔼_{0}(D_{i}|Y = y) = ∫_{0}^{T} ℙ_{0}(X_{t} = i|Y = y) és végül

4. n̂_{i} = 𝔼_{0}(n_{i}|Y = y) = ∑_{k=1}^{n} ℙ_{0}(X_{t_{k}} = i|Y = y).

A (3.2) összegben lévő első tagot elhagyjuk, mert a minta méretének növekedésével ennek mérete elhanyagolható a többihez képest. Ezek után egyszerűen ellenőrizhető, hogy a q_{i} = ∑_{j≠i} q_{ij} feltételek mellett a fennmaradó kifejezés egyértelmű maximuma a

q̂_{ij} = \frac{m̂_{ij}}{D̂_{i}} , i ≠ j és a λ̂_{i} = \frac{n̂_{i}}{D̂_{i}} , i = 1, …, r

helyeken vétetik fel. Ez lesz az algoritmus M-lépése. Az E-lépés meghatározása pedig m̂_{ij}, D̂_{i} és n̂_{i} feltételes várható értékek meghatározásából áll, ami a fentinél jóval bonyolultabb feladat. Vezessük be a következő jelöléseket. Legyen

α_{t}(i, ϑ) = π^{ϑ} (∏_{k=1}^{N_{t}} f^{ϑ}(y_{k})) \overline{F}_{ϑ}(t − N_{t})1_{i}  (3.3)

és

β_{t}(i, ϑ) = 1_{j}^{⊤}f^{ϑ}(t_{N_{t}+1} − t) (∏_{k=N_{t}+2}^{n} f^{ϑ}(y_{k})) 1,  (3.4)

ahol f^{ϑ} a (2.14) alatt definiált és

\overline{F}_{ϑ}(t) = exp [(Q − Λ)t] .  (3.5)

Itt azzal a konvencióval élünk, hogy ∏_{p=k}^{n}(·) = I, amennyiben k > n. Ezzel egyidejűleg legyen

L_{ϑ}^{∗}(k) = \{ \begin{matrix}π^{ϑ} ∏_{p=1}^{k}f^{ϑ}(y_{p}) & , k = 1, 2, …, n \\ π^{ϑ} & , k = 0 \end{matrix}

és

R_{ϑ}^{∗}(k) = ∏_{p=k}^{n}f^{ϑ}(y_{p})1, k = 1, 2, …

Könnyen láthatóan L′_{ϑ^{0}}(k) kielégíti az

L_{ϑ}^{∗}(0) = π^{ϑ}, L_{ϑ}^{∗}(k) = L_{ϑ}^{∗}(k − 1)f^{ϑ}(y_{k}), k = 1, 2, …, n,  (3.6)

R_{ϑ}^{∗}(k) pedig az

R_{ϑ}^{∗}(n + 1) = 1, R_{ϑ}^{∗}(k) = f^{ϑ}(y_{k})R_{ϑ}^{∗}(k + 1), k = n, …, 1  (3.7)

rekurziót. Ryden [23]-ben megmutatta, hogy ekkor

m̂_{ij} = \frac{q_{ij}^{0}}{L(ϑ^{0}, y)} ∫_{0}^{T} α_{t}(i, ϑ^{0})β_{t}(j, ϑ^{0})dt,  (3.8)

továbbá, hogy

D̂_{i} = \frac{1}{L(ϑ^{0}, y)} ∫_{0}^{T} α_{t}(i, ϑ^{0})β_{t}(i, ϑ^{0})dt  (3.9)

és, hogy

n̂_{i} = \frac{π^{0}}{L(ϑ^{0}, y)}L_{ϑ^{0}}^{∗}(k)R_{ϑ^{0}}^{∗}(k + 1),  (3.10)

ahol L(ϑ, y) az y mintához tartozó likelihood-függvény értéke a ϑ helyen, ami a 2.6.2 tétel alapján π^{ϑ} (∏_{k=1}^{n} f^{ϑ}(y_{k})) 1. Szokás a (3.3), illetve (3.4) alatti valószínűséget előre tekintő, illetve hátráló sűrűségnek nevezni. Hasonlóan a (3.6) és (3.7) alatti rekurziókat előre tekintő, illetve hátráló rekurzióknak. A további átalakítások során a képletek egyszerűsítése érdekében elhagyjuk a ϑ-tól való függést, hiszen a továbbiakban úgyis minden kiértékelés ϑ^{0} mellett folyik. Világos, hogy

α_{t}(i) = L^{∗}(N_{t})\overline{F} (t − t_{N_{t}})1_{i} és β_{t}(i) = 1_{i}^{⊤}f(t_{N_{t}+1} − t)R′(N_{t} + 2),

tehát m̂_{ij}, D̂_{i} és n̂_{i} kifejezései az alábbiak szerint alakíthatók tovább, az L = L(ϑ^{0}, y) egyszerűsítő jelölés bevezetése mellett.

\frac{L}{q_{ij}} m̂_{ij} = ∫_{0}^{T} L^{∗}(N_{t})\overline{F} (t − t_{N_{t}})1_{i}1_{j}^{⊤}f(t_{N_{t}+1} − t)R^{∗}(N_{t} + 2)dt

= ∑_{k=1}^{n} ∫_{t_{k−1}}^{t_{k}−} L^{∗}(N_{t})\overline{F} (t − t_{N_{t}})1_{i}1_{j}^{⊤}f(t_{N_{t}+1} − t)R^{∗}(N_{t} + 2)dt,

ami pedig - tekintve, hogy a [t_{k−1}, t_{k}) intervallumon N_{t} = k − 1 -, így írható tovább:

= ∑_{k=1}^{n} ∫_{t_{k−1}}^{t_{k}−} L^{∗}(k − 1)\overline{F} (t − t_{k−1})1_{i}1_{j}^{⊤}f(t_{k} − t)R^{∗}(k + 1)dt

= ∑_{k=1}^{n} L^{∗}(k − 1) ∫_{t_{k−1}}^{t_{k}−} \overline{F} (t − t_{k−1})1_{i}1_{j}^{⊤}f(t_{k} − t)dtR^{∗}(k + 1)

Ugyanezt elvégezve D̂_{i}-re és n̂_{i}-ra, azt kapjuk, hogy

LD̂_{i} = ∑_{k=1}^{n} ∫_{t_{k−1}}^{t_{k}−} L^{∗}(k − 1)\overline{F} (t − t_{k−1})1_{i}1_{i}^{⊤}f(t_{k} − t)R^{∗}(k + 1)dt,

illetve

Ln̂_{i} = ∑_{k=1}^{n} L^{∗}(k)1_{i}1_{i}^{⊤}R^{∗}(k + 1) = ∑_{k=1}^{n} (L^{∗})^{i}(k)R_{i}^{∗}(k + 1).

Ezzel megadtuk az algoritmus M, illetve E lépéseit, az implementációval kapcsolatos további gyakorlati problémákat a következő alfejezetben tárgyaljuk.

3.3. Az algoritmus implementálása

Az előző fejezetben képletszerűen felírtuk az algoritmus E, illetve M lépését, most rátérünk az algoritmus implementálásának tárgyalására. Ehhez kicsit át kell alakítanunk a képleteinket, ugyanis a jelenlegivel a következő probléma merül fel. Megmutatható, hogy az előretekintő illetve hátráló sűrűségek exponenciális sebességgel tartanak 0-hoz, amennyiben n ↦ ∞. Emiatt az előretekintő illetve hátráló rekurziók nem feltétlenül adnak pontos eredményt, ha a valós kalkulációk során olyan számítógépet használunk, amelynek a számítási korlátait meghaladják az esetlegesen felmerülő kiugróan kicsi vagy nagy értékek. W.J.J. Roberts, Y. Ephraim és E. Dieguez 2006-os eredményükben - lásd [21] - bemutattak egy skálázási eljárást, melynek alkalmazása megóv minket a túlcsordulástól. Először ezt ismertetjük. Ennek érdekében tetszőleges y vektor esetén k = 1, …, dim y -ra y^{k} jelölje az y első k koordinátájából álló vektort, míg y^{0} legyen konstans 0. Mint ismeretes, két valószínűségi változó - monjuk X és Y - együttes és az egyiknek a másikra vonatkozó feltételes sűrűségfüggvénye között - amennyiben azok léteznek - az alábbi kapcsolat áll fenn:

f_{(X,Y)}(x, y) = f_{X|Y} (x|y)f_{Y} (y).

Ezt n-szer alkalmazva - kissé pongyola jelöléssel, elhagyva a ϑ-tól való függést - fennáll az L(y) = ∏_{l=1}^{n}f_{Y^{l}|Y^{l−1}}(y_{l}|y^{l−1}) azonosság, amit a c_{l} = f_{Y^{l}|Y^{l−1}}(y_{l}|y^{l−1}) jelöléssel együtt egyszerüen L = ∏_{l=1}^{n}c_{l} alakba írhatunk. Ezzel az m̂_{ij}-re vonatkozó korábbi alfejezetben szereplő előállítás az alábbiak szerint alakul át.

m̂_{ij} = ∑_{k=1}^{n} \frac{q_{ij}^{0}}{c_{k}} (π∏_{l=1}^{k−1} \frac{f(y_{l})}{c_{l}}) ∫_{t_{k−1}}^{t_{k}−} \overline{F} (t − t_{k−1})1_{i}1_{j}^{⊤}f(t_{k} − t)dt (∏_{l=k+1}^{n} \frac{f(y_{l})}{c_{l}}1).

Itt újabb jelöléseket vezetünk be a két zárójeles kifejezésre, legyen

L_{ϑ}(k) = π^{ϑ}∏_{l=1}^{k−1} \frac{f_{ϑ}(y_{l})}{c_{l}(ϑ)} és R_{ϑ}(k) = ∏_{l=k}^{n} \frac{f_{ϑ}(y_{l})}{c_{l}(ϑ)}1.  (3.11)

Ezenkívűl állapodjunk még meg abban, hogy L_{ϑ}(0) = π^{ϑ}, illetve R_{ϑ}(n + 1) = 1. A továbbiakban itt sem írjuk majd ki mindig a ϑ-tól való függést, hogy a kiértékelés milyen ϑ mellett történik, az a kontextusból úgyis kiderül. Mivel egyrészt f_{Y^{k}}(y^{k}) = ∏_{l=1}^{k}c_{l}, másrészt f_{Y^{k}}(y^{k}) = π∏_{l=1}^{k}f(y_{l})1, azért

c_{k} = π∏_{l=1}^{k−1} \frac{f(y_{l})}{c_{l}}f(y_{k})1 =L(k − 1)f(y_{k})1.  (3.12)

Ezzel (3.11) és (3.12) alapján a következő rekurziók írhatók fel:

L(k) = \frac{L(k − 1)f(y_{k})}{L(k − 1)f(y_{k})1}, k = 1, …, n

és

R(k) = \frac{f(y_{k})R(k + 1)}{L(k − 1)f(y_{k})1}, k = n, …, 1.

Megmutatható (lásd [21]), hogy az ily módon átskálázott rekurziók esetén egyrészt

[L(k)]_{j} = ℙ_{0}(X_{t_{k}} = j|Y^{k} = y^{k}),

tehát [L(k)]_{j} legfeljebb 1, másrészt

[R(k)]_{j} < \frac{1}{[L(k − 1)]_{j}} ,

azaz [R(k)]_{j} felülről korlátos és [R(k)]_{j}[L(k − 1)]_{j} < 1.

Ezek után tehát felírható az

m̂_{ij} = ∑_{k=1}^{n} \frac{q_{ij}^{0}}{c_{k}} ∫_{t_{k−1}}^{t_{k}−} L(k − 1)\overline{F} (t − t_{k−1})1_{i}1_{j}^{⊤}f(t_{k} − t)R(k + 1)dt

= ∑_{k=1}^{n} \frac{q_{ij}^{0}}{c_{k}} ∫_{t_{k−1}}^{t_{k}−} [1_{i}^{⊤}\overline{F} (t − t_{k−1})^{⊤}L(k − 1)^{⊤}] [R(k + 1)^{⊤}f(t_{k} − t)^{⊤}1_{j}] dt

= ∑_{k=1}^{n} \frac{q_{ij}^{0}}{c_{k}} ∫_{t_{k−1}}^{t_{k}−} [\overline{F} (t − t_{k−1})^{⊤}L(k − 1)^{⊤}R(k + 1)^{⊤}f(t_{k} − t)^{⊤}]_{ij} dt

= ∑_{k=1}^{n} \frac{q_{ij}^{0}}{c_{k}} ∫_{t_{k−1}}^{t_{k}−} \{[f(t_{k} − t)R(k + 1)L(k − 1)\overline{F} (t − t_{k−1})]^{⊤}\}_{ij} dt  (3.13)

összefüggés, valamint amennyiben ⊙ jelöli az elemenkénti mátrix szorzást és n̂ jelöli azt az r × 1-es vektort, aminek i. eleme n̂_{i}, akkor

n̂ = ∑_{k=1}^{n} L(k)^{⊤} ⊙ R(k + 1),

továbbá világos, hogy

D̂_{i} = \frac{m̂_{ii}}{q_{ii}^{0}} .

Most már csak az maradt hátra, hogy megmutassuk, miként lehet hatékonyan kiértékelni a (3.13) előállításban szereplő integrálokat. Ennek érdekében jelölje m̂ azt az r × r-es mátrixot, amelyben az i. sor j. eleme m̂_{ij}. Ekkor, könnyen láthatóan

m̂ = Q^{0} ⊙ ∑_{k=1}^{n} \frac{1}{c_{k}} ∫_{t_{k−1}}^{t_{k}−} [f(t_{k} − t)R(k + 1)L(k − 1)\overline{F} (t − t_{k−1})]^{⊤} dt.

Ebből m̂ transzponáltjára a következő adódik:

m^{⊤} = (Q^{0})^{⊤} ⊙ ∑_{k=1}^{n} \frac{1}{c_{k}} ∫_{t_{k−1}}^{t_{k}−} f(t_{k} − t)R(k + 1)L(k − 1)\overline{F} (t − t_{k−1})dt,

ami az I_{k} = ∫_{t_{k−1}}^{t_{k}−} f(t_{k} −t)R(k+1)L(k−1)\overline{F} (t−t_{k−1})dt jelölés bevezetésévél tovább egyszerűsödik:

m̂ = Q^{0} ⊙ ∑_{k=1}^{n} \frac{I_{k}^{⊤}}{c_{k}}  (3.14)

Itt figyelembe véve a (2.14) és (3.5) előállításokat és alkalmazva az s = t−t_{k−1} integrál-transzformációt, azt kapjuk, hogy

I_{k} = ∫_{0}^{y_{k}} exp((Q^{0} − Λ^{0})(y_{k} − y))Λ^{0}R(k + 1)L(k − 1) exp((Q^{0} − Λ^{0})y)dy.

Ennek kiszámításához alkalmazzuk a 2.3.7 lemmát az A = Q^{0} − Λ^{0} és B = Λ^{0}R(k + 1)L(k − 1) szereposztással, legyen ennek érdekében

S_{k} = ( \begin{matrix}Q^{0} − Λ^{0} & Λ^{0}R(k + 1)L(k − 1) \\ 0 & Q^{0} − Λ^{0} \end{matrix} )

Ekkor tehát a 2.3.7 lemmából az következik, hogy I_{k} nem más, mint az exp(S_{k}y_{k}) mátrix jobb felső r × r-es blokkja. Így tehát I_{k} kiszámítását visszavezettük egy magasabb dimenziós mátrix exponenciálisának meghatározására. Ezt pedig az 2.3 fejezetben ismertetett Padé-féle diagonális approximáció módszerének segítségével tesszük meg, megjegyezve, hogy a mátrix exponenciális kiszámítására alkalmas "expm" Matlab^{TM} függvény éppen ezt az algoritmust használja.

Most már minden adva van ahhoz, hogy megfogalmazhassuk az algoritmusunkat.

3.3.1. Algoritmus (Rydén EM algoritmusa). 1. Bemenő adat: y = (y_{1}, y_{2}, …, y_{n}) megfigyelés egymás utáni köztes időkből.

2. Kezdeti ϑ^{0} = (Q^{0}, Λ^{0}) becslés meghatározása a minta alapján.

3. Tegyük fel, hogy ϑ^{l−1}-et már meghatároztuk, ekkor az l. iterációs lépés:

(a) π = π(Q^{l−1}, Λ^{l−1}) kezdeti eloszlás meghatározása.

(b) L(0) := π és k = 1, …, n-ig L(k) := \frac{L(k−1)f(y_{k})}{L(k−1)f(y_{k})1}. (Itt, és a további iterációs lépésekben is f-et ϑ^{l−1} mellett kell kiértékelni.)

(c) R(n + 1) := 1 (ahol 1 jelöli az r × 1-es csupa 1 vektort) és k = 1, …, n-ig R(k) := \frac{f(y_{k})R(k+1)}{L(k−1)f(y_{k})1}.

(d) k = 1, …, n-ig c_{k} := L(k − 1)f(y_{k})1.

(e) k = 1, …, n-ig

S_{k} := ( \begin{matrix}Q^{l−1} − Λ^{l−1} & Λ^{l−1}R(k + 1)L(k − 1) \\ 0 & Q^{l−1} − Λ^{l−1} \end{matrix} )
R_{k} := exp(S_{k}y_{k}) (Padé mátrix approximáció) és
(I_{k})_{ij} := (R_{k})_{i,r+j} .

(f)

m̂ := Q^{l−1} ⊙ ∑_{k=1}^{n} \frac{I_{k}^{⊤}}{c_{k}}

(g)

n̂ := ∑_{k=1}^{n} L(k)^{⊤} ⊙ R(k + 1)

(h) Ekkor

l(ϑ^{(l−1)}) = l(Q^{(l−1)}, Λ^{(l−1)}) = ∑_{k=1}^{n} log c_{k}

és az új ϑ^{l} becslés pedig:

q_{ij}^{l} = q_{ii}^{l−1} \frac{m̂_{ij}}{m̂_{ii}} , i, j = 1, …, r, i ≠ j és λ_{i}^{l} = q_{ii}^{l−1} \frac{n̂_{i}}{m̂_{ii}} , i = 1, …, r.

4. Kimenő adat: {ϑ^{l}} becsléssorozat.

Nyitva maradt még az algoritmus 2. pontjának, azaz a kiindulási becslés kezdeti értékének megválasztásának a kérdése. Erre Rydén [23]-ban javasol egy eljárást, aminek segítségével kellően közel lehet kerülni a kezdeti paraméterekhez. Ennek azonban az a hátránya, hogy több fázisból áll és használja például a rejtett Markov modell paraméreinek a becslésére vonatkozó EM algoritmust, ami önmagában egy, az imént bemutatott algoritmushoz hasonló kaliberű procedúra. Bár a rejtett Markov modell paramétereire vonatkozó becslés kiszámítására létezik egy eljárás a Matlab^{TM} "Statistics Toolbox" csomagjában, most azonban mégsem ezt az utat választjuk, hanem bemutatunk egy egyszerűbb intuitív eljárást, mely jóval pontatlanabb, ám egyszerűségéből adódóan rendkívül gyorsan számolható és legalább a paraméterek egymáshoz való nagyságrendi viszonya eltalálható vele.

Egy MMPP esetén a λ_{i}, (i = 1, …, r) intenzitáshányadok véletszerűen váltakoznak. Ha a mögöttes Markov folyamat éppen az i. állapotban van, akkor az ez alatt bekövetkező érkezések \frac{1}{λ_{i}}

köztes idő sorszáma 1  2  3  4  5  6
köztes idő  0.0042 0.0118 0.0106 0.0380 0.0254 0.0200
valódi állapot  1  1  1  2  2  1
becsült állapot  1  1  1  2  2  1

3.1. táblázat. A kiinduló becslés pontossága

köztes idő sorszáma 7  8  9  10  11  12
köztes idő  0.0050 0.0064 0.0035 0.0284 0.0105 0.0413
valódi állapot  1  1  1  1  1  1
becsült állapot  1  1  1  2  1  2

3.2. táblázat. A kiinduló becslés pontossága (folytatás)

várható értkű és \frac{1}{λ_{i}} szórású exponenciális eloszlású valószínűségi változók, ezért úgy képzeljük, hogy az ez alatt az idő alatt bekövetkező érkezések közti időkből vett megfigyeléseink \frac{1}{λ_{i}} környéken kell, hogy szóródjanak. Ez alapján megpróbáljuk r csoportra osztani az Y megfigyelésünket, mégpedig úgy, hogy az egyes csoprtokon belüli elemek az euklideszi távolságra nézve közel legyenek egymáshoz. Ezt pedig az euklideszi távolságra vonatkozó k-means klaszterező eljárással - egy magyar nyelvű bevezető található például az [1] alatti jegyzetben - tesszük meg, megjegyezve, hogy ez az eljárás már szintén előre programozottan rendelkezésünkre áll a Matlab^{TM}-ban. Ha elkészültek a klasztereink, akkor azt feltételezzük, hogy az egy klaszteren belüli megfigyelsések f.a.e. exponenciális valószínűségi változók. Ekkor a paraméter maximum-likelihood becslése a mintaátlag reciproka, amit könnyen megkaphatunk, hiszen az egy klaszteren belüli elemek mintaátlaga éppen a klaszterközéppont. Összefoglalva, legyen adott az Y = (Y_{1}, Y_{2}, …, Y_{n}) köztes időkből vett minta, és tegyük fel, hogy a klaszterezés az {1, 2, …, n} halmaz alábbi partícióját eredményezi:

{1, 2, …, n} = ⋃_{i=1}^{r}K_{i}.

Ekkor a λ_{i}-k becslésére a

λ̂_{i} = \frac{n}{_{∑j∈K_{i}}Y_{j}}, i = 1, 2, …, r

összefüggés adódik. Ennek az eljárásnak a kritikus pontja az, hogy a köztes idők közül valójában nem sorolható mindegyik köztes idő tisztán valamelyik állapothoz, ugyanis van olyan, ami egy adott állapotban indul, de már egy másik állapotban ér véget. Ettől a nehézségtől most az egyszerűség kedvéért eltekintünk. Illusztráció gyanánt álljon itt a következő példa.

3.3.1. Példa. 9000 elemű mintát szimuláltunk a q_{12} = 1 és q_{21} = 10, illetve λ_{1} = 100 és λ_{2} = 10 paraméterekkel rendelkező két állapotú MMPP-ből. A minta első néhány értékét a (3.1) és (3.2) táblázatokban foglaltuk össze. A "valódi állapot" mező mutatja a mögöttes Markov-folyamat állapotát a megfelelő érkezés bekövetkeztekor, a "becsült állapot" mező pedig a klaszterezés során kapott javaslatot.

Látjuk tehát, hogy a klaszterező eljárás több ízben is hibázik a köztes idők klasszifikálása során. A klaszterközéppontokra c_{1} = 0.0075 illetve c_{2} = 0.041 adódik, melyből a λ̂_{1} = \frac{1}{c_{1}} = 133.6267 illetve a λ̂_{2} = \frac{1}{c_{2}} = 24.3681 kiinduló becsléseket kapjuk.

A Q generátormátrix becslésénél visszautalunk a korábbi fejezetekben ismertetett összefüggésekre, nevezetesen arra, hogy amennyiben a mögöttes Markov folyamat éppen az i. állapotban van, úgy az ott töltött ideje q_{i} paraméterű exponenciális eloszlású valószínűségi változó, majd ezek után a \frac{q_{ij}}{q_{j}} valószínűségek szerint ugrik át egy másik j ≠ i állapotba. Ebből kiindulva a klaszterezett minta alapján megpróbáljuk reprodukálni az egyes állapotokban eltöltött időket.

3.3.2. Példa. Például a (3.1) és (3.2) táblázatokban összefoglalt minta esetén a klaszterezés szerint a Markov folyamat kezdetben az 1. állapotban volt, majd a 4. érkezési idő már a 2. állapot szerinti intenzitással történt. Ez azt jelenti, hogy az állapotváltás valahol T_{3} és T_{4} között volt. Mivel a 2. állapot szerinti érkezésre az állapotváltástól számítva várhatóan \frac{1}{λ_{2}} időt kell várni, ezért a Markov folyamat 1. állapotban eltöltött idejét - úgy értve, hogy az első állapotváltásig - a ξ_{11} = T_{4} − \frac{1}{λ̂_{2}} kifejezéssel becsüljük. Ezután megkeressük a következő csupa 1-esből álló szakaszt. Látható, hogy ez a 6. megfigyeléstől kezdődik és a 9. megfigyelésnél ér véget. Tehát itt azt feltételezzük, hogy az állapotváltás T_{9} és T_{10} között volt. Így az ebben az állapotban töltött időt most a ξ_{12} = T_{10}− \frac{1}{λ̂_{2}} −T_{5}− (y_{6} − \frac{1}{λ̂_{1}} ) kifejezéssel becsülhetjük. Ezzel egy ξ_{11}, ξ_{12}, … mintát kapunk q_{1} paraméterű exponenciális eloszlásból, majd ezután ismét ezek mintaátlagának reciprokával becsüljük a q_{1} paramétert.

Mindezt precízebben megfogalmazva, jelölje K_{i} ∈ {1, 2, …, r} az Y_{i} megfigyelés klaszterének sorszámát. Legyen H_{ij} a j.-edik maximális méretű csupa i-ből állő blokk elemeinek az indexe Kban (i = 1, …, r és H_{ij} = ∅, ha már nincs több ilyen blokk). Legyen n_{i} = #{j : H_{ij} ≠ ∅}. Ezzel tehát ismét partícionáltuk az {1, 2, …, n} halmazt:

{1, 2, …, n} = ⋃_{i=1}^{r}⋃_{j=1}^{∞}H_{ij}.

H_{ij} legkisebb eleme legyen m_{ij}, legnagyobb eleme pedig M_{ij}. Ekkor legyen

ξ̂_{ij} = (T_{M_{ij}+1} − T_{m_{ij}−1} − \frac{1}{λ̂_{K_{M_{ij} +1}}} ) − (Y_{m_{ij}} − \frac{1}{λ̂_{K_{M_{ij}}}} )

= T_{M_{ij}+1} − T_{m_{ij}} + \frac{1}{λ̂_{K_{M_{ij}}}} − \frac{1}{λ̂_{K_{M_{ij} +1}}} , j = 1, …, n_{i}

Ebből már megkapható q_{i} becslése:

q̂_{i} = \frac{n_{i}}{∑{n_{i} \atop j=1}ξ̂_{ij}} .

Ezzel megbecsültük a Q mátrix főátlóját, hátravan még a többi elem becslése. Ehhez rögzített i ∈ {1, …, r} esetén minden j ≠ i-hez készítsük el az A_{ij} = #{0 ≤ k ≤ n − 1 : K_{k} = i, K_{k+1} = j}. Itt tehát A_{ij} adja meg a klaszterezés alapján várt (i, j) állapotváltások számát. Legyen

A_{i·} = ∑_{j=1,j≠i}^{r} A_{ij},

vagyis A_{i·} azt adja meg, hogy összesen hány állapotváltás volt az i állapotból. Ekkor a korábban említettek miatt j ≠ i esetén q_{ij} természetes becsléseként

q̂_{ij} = q̂_{i} \frac{A_{ij}}{A_{i·}}

adódik. i = 1, …, r esetén legyen A_{ii} = −A_{i·} és jelöljük B-vel azt az r × r-es mátrixot, ahol az i-edik sor j-edik eleme \frac{A_{ij}}{A_{i}}. Ekkor tehát a Q mátrix becslése

Q̂ = CB,

ahol C =diag(q), vagyis a q elemeiből képzett diagonális mátrix.

3.3.3. Példa. A 3.3.1 példa esetén ekkor q̂_{12} = 5.2394 illetve q̂_{21} = 14.2649 adódik kiinduló becslésként.

Összefoglalásképpen, bár a módszer nem szolgáltat túl nagy pontosságot és érzékeny arra, hogy például a λ_{i} paraméterek mennyire jól szeparálhatók egymástól - tehát kellő pontosságú eredmény inkább akkor várható, ha előre tudjuk, hogy a λ_{i} intenzitáshányadok között nincsenek nagyon közeliek - , arra azonban mindenképpen alkalmas, hogy gyorsan jussunk intuitív kezdeti becslésekhez. A következőkben ismertetjük e becslés kiszámításának konkrét implementációját.

3.3.2. Algoritmus (Kiinduló becslés keresése). 1. Bemenő adatok: y n dimenziós vektor, r (az illesztendő modell rendje)

2. Klaszterezés k-means eljárás segítségével, kimenetek:

(a) k n dimenziós vektor, k_{i} jelenti y_{i} klaszterének sorszámát (1-től r-ig terjedhet)

(b) c = (c_{1}, …, c_{r}) vektor, amely tartalmazza a klaszterközpontokat

3. 1. kimenet: a λ_{i} paraméterek becslése: λ̂_{i} := \frac{1}{c_{i}} , i = 1, …, r.

4. Levágjuk az utolsó blokkot a mintából (az ugyanis egy csonka megfigyelést adna), n := dim y

5. inicializálás:

(a) S legyen r × 1-es vektor S(i) := \{ \begin{matrix}0 & , ha & i ≠ k_{1} \\ y_{1} & , ha & i = k_{1} \end{matrix}

(b) N legyen r × 1-es vektor, csupa nulla elemekkel.

6. Első iteráció: i = 2, …, n-ig:

(a) Ha k_{i} = k_{i−1}, akkor S(k_{i−1}) := S(k_{i−1}) + y_{i}, különben

(b) S(k_{i−1}) := S(k_{i−1}) + y_{i} − c_{k_{i}}; S(k_{i}) := S(k_{i}) + c_{k_{i}}; N(k_{i−1}) := N(k_{i−1}) + 1.

7. Utólagos módosítás: S(k_{n}) := S(k_{n}) − c_{k_{n}}.

8. A q_{i} paraméterek becslése q̂_{i} := \frac{S(i)}{N(i)}, i = 1, …, r.

9. Az A mátrix inicializálása r × r-es csupa nulla mátrixként, majd A feltöltése, 2. iteráció: i = 1, …, n − 1-ig:

(a) Ha k_{i} ≠ k_{i+1}, akkor A_{k_{i},k_{i+1}} := A_{k_{i},k_{i+1}} + 1.

10. i = 1, …, r-ig A_{i·} := ∑_{j=1,j≠i}^{r} A_{ij} és A_{ii} := −A_{i·}.

11. B :=diag(\frac{1}{A_{1}} , …, \frac{1}{A_{r}})A és C :=diag(q_{1}, …, q_{r}).

12. 2. kimenet: Q becslése: Q̂ = CB

Ezzel tehát befejeztük az implementáció kérdésének tárgyalását, mielőtt azonban rátérnénk az EM algoritmus gyakorlati alkalmazására, megemlítjük még annak egy fontos tulajdonságát. Nevezetesen azt, hogy megőrzi a nullákat, vagyis, ha a ϑ^{0} kezdeti becslés esetén q_{ij}^{0} = 0 valamilyen i-re és j-re, akkor a következő iteráció során kiszámolt becslés is rendelkezni fog ezzel a tulajdonsággal, azaz q_{ij}^{1} = 0 szintén fenáll. Ez könnyen látható abból a tényből, hogy q_{ij}^{1} = \frac{m̂_{ij}}{D̂_{i}} és a (3.14) előállítás alapján m̂_{ij} = q_{ij}^{0} (∑_{k=1}^{n} \frac{I_{k}^{⊤}}{c_{k}} )_{ij}. Ez a tulajdonság igen hasznos lehet bizonyos speciális struktúrájú MMPP-k becslésekor, nevezetesen amikor tudjuk, hogy a Q mátrix bizonyos elemei nullával egyenlők.

Markov-modulált Poisson folyamat a

gyakorlatban

Jelen fejezet az előző fejezetben bemutatott algoritmus alkalmazása egy gyakorlati példára. Az 1. alfejezetben ismertetjük, hogy milyen eljárást használtunk MMPP-ből vett minta szimulálására. A 2. alfejezet fontos önálló eredménye a 4.2.1 lemma, - ami a 2.3.8 lemmán alapul -, enélkül ugyanis igencsak nehézkes lenne az alfejezet során bemutatásra kerülő ℙ(ζ_{k+1} < K) farokvalószínűségek numerikus számolása. A 3. alfejezben közöljük vizsgálódásaink numerikus eredményeit.

4.1. Minta generálása adott paraméterű MMPP-ből

Ebben az alfejezetben bemutatjuk, hogy miként tudunk mintát generálni egy megadott paraméterű Markov modulált Poisson folyamatból vett köztes időkből. E célból legyen tehát a ϑ = (Q, Λ) paraméter rögzített. Kihasználva a mögöttes Markov folyamat fejlődését leíró tulajdonságokat és azt, hogy egy λ intenzitású Poisson folyamat köztes idői független, λ paraméterű exponenciális eloszlásúak, az alábbi módon járhatunk el. Véletlenszerűen kiválasztjuk a Markov folyamat kezdeti állapotát a π stacionárius eloszlás szerint. Ekkor a Markov folyamat a megadott állapottól függő - nevezetesen q_{i} paraméterű, amennyiben i ez az állapot - exponenciális eloszlás szerinti időt tölt itt el, majd szintén a megadott állapottól függően átvált egy másik állapotba. Rögzítsük ezt az időt és jelöljük τ_{1}-gyel. Ez alatt az időintervallum alatt a folyamatunk szintén az adott állapottól függő - nevezetesen λ_{i} intenzitáshányadú, amennyiben i ez az állapot - homogén Poisson folyamatként viselkedik. Ez azt jelenti, hogy a megfigyelések az adott paraméter szerinti, f.a.a. exponenciális eloszlású valószínűségi változók szerint érkeznek. Generáljunk tehát ezekből egy mintát egészen addig, amíg a mintaelemek összege még kisebb, mint τ_{1}. Ezután határozzuk meg a \frac{q_{ij}}{q_{i}} valószínűségek szerint a Markov folyamat új állapotát, majd ismét szimuláljuk az itt eltöltött τ_{2} időt, mely alatt ismét konstans a Poisson folyamat intenzitáshányada. Ezt folytassuk egészen addig, amíg elő nem állt a kellő elemszámú minta. Ugyanez a következő módon fogalmazható meg egzakt módon.

4.1.1. Algoritmus (MMPP szimuláció). 1. Bemenet: ϑ = (Q, Λ) paraméter, n mintaelemszám.

2. π stacionárius eloszlás meghatározása a πQ = 0, π1 = 1 egyenletrendszer megoldásával.

3. Válasszunk egy s_{0} számot véletlenszerűen az E = {1, 2, …, r} halmazból π szerint. (i valószínűsége π_{i}).

4. A Markov folyamat aktuális állapotát az s változóban tároljuk, s inicializálása:

s := s_{0}.

.

5. Inicializálás: y := 0, t :=vélexp(q_{s}) és a :=vélexp(λ_{s}).

6. Ciklus: amíg dim y < n + 1 :

(a) Ciklus: amíg (a < t) és (dim y < n + 1)

i. y := (y, a − ∑_{i=1}^{dim y} y_{i}). ii. a := a+vélexp(λ_{s}).

(b) a := t.

(c) e ∈ E\s véletlen állapot választása a \frac{q_{sj}}{q_{s}} valószínűségek szerint

(d) s := e.

(e) t := t+vélexp(q_{s}) és a := a+vélexp(λ_{s}).

7. Ciklus: i = 1, …, n-re y_{i} := y_{i+1}

8. kimenet: n hosszú y sorvektor MMPP köztes időkből.

Itt vélexp(λ) alatt egy λ paraméterű exponencális eloszlás szerint generált véletlen számot értünk. Ha van egy olyan generátorunk, ami egyenletes eloszlású U véletlen számot állít elő a (0, 1) intervallumból, akkor tetszőleges abszolút folytonos F eloszlásból vett véletlen minta generálható vele, ugyanis F^{−1}(U) eloszlása éppen F lesz. Mi a Matlab beépített "randexp" függvényét használtuk exponenciális eloszlású minta generálására.

4.2. Egy biztosításmatematikai feladat

A Markov-modulált Poisson folyamatnak rengeteg gyakorlati alkalmazása van, használata igen közkedvelt az olyan sorbanállási illetve kiszolgálási modellek esetén, ahol az egyes kiszolgálási idők intenzitása az időben véletlenszerűen váltakozik. Jelen dolgozat a továbbiakban a következő problémát vizsgálja. Tételezzük fel, hogy adott egy biztosító társaság, akinek a biztosítottjaival meghatározott káresemények történhetnek. Most nem koncentrálunk a kárnagyságra, csupán a károk száma érdekel bennünket. Ha egy adott pillanattól kezdve (ezt választjuk origónak) elkezdjük megfigyelni a károk bekövetkezésének idejét, akkor ezek meghatároznak egy konfigurációt az ℝ^{+} halmazon. Ha n(B) jelöli a konfigurációból az ℝ^{+} egy B Borel-halmazába eső pontok számát, akkor ezáltal egy N pontfolyamat n realizációját kapjuk. A gyakorlatban persze a B Borel-halmaz nyilván többnyire egy intervallum, azaz a biztosítót az érdekli, hogy egy meghatározott összefüggő időperiódusban hány kár következik be. Persze ennek nem kell feltétlenül így lennie, hiszen érdekelheti például a biztosítót a csak tavasszal bekövetkező károk statisztikája, ami máris egy nem összefüggő halmazt eredményez, bár ekkor is vissza tudjuk vezetni a problémát diszjunkt intervallumok uniójára.

Az alapszituáció a következő: a társaság a kárfolyamatot MMPP-vel próbálja modellezni. Képzeljük el például, hogy a biztosító portfóliója gépjármű felelősség-biztosításokból áll, ekkor gépkocsibalesetek bekövetkezési idejét kísérjük figyelemmel. Köztudott, hogy esős időszak esetén több káresemény történik, azon kívűl kevesebb. Az, hogy mikor van esős időszak, nem determinisztikus, feltehető például, hogy az esős illetve napos időszakok egy kétállapotú Markov folyamat szerint váltakoznak. Világos, hogy ez esetben természetes feltevésként adódik az MMPP, mint kárfolyamat. Nyilván ennél több állapot is elképzelhető, például az előbbi esetnél maradva, több időjárástól függő kategóriát is bevezethetnénk.

Egyéb modellre is gondolhatunk, például egy egészségbiztosításokból álló portfólió esetén lehet a nagyobb intenzitással járó véletlentől függő időszak egy járvány által okozott. A továbbiakban nem kívánunk különösebb hangsúlyt fektetni a konkrét interpretációra, csupán magára a modellre fókuszálunk, azaz a kárfolyamatról csak mint matematikai objektumról beszélünk és ezt modellezzük egy MMPP-vel. A feladatunk a következő lesz. Feltesszük, hogy a biztosító több egymást követő egyenlő időintervallumon keresztül (pl. nap, hónap, év) rögzíti a kárstatisztikát.

A cél előrejelzést adni az éppen soron következő időintervallumban bekövetkező károk ζ számára, illetve az előre rögzített α megbízhatósági szintekhez tartozó K(ζ) = inf{K ∈ ℕ_{0} : ℙ(ζ < K) ≥ 1 − α} kvantilis értékekre. Előbbi egyfajta legjobb tudásunk szerint elvárható becslést (biztosítási szaknyelven Best Estimate) szolgáltat, míg az utóbbi érték arról informál, hogy mi az a szám, aminél nagyobb kárszám már egyhez közeli valószínűséggel nem következik be, tehát ez tulajdonképpen egy bizonyos értelemben a legnagyobb lehetséges veszteséget jelenti.

Tegyük fel, hogy az említett egyenlő nagyságú időperiódusok hossza h, az egyszerűség kedvéért hívjuk ezeket napoknak. Tegyük fel, hogy k darab napból álló megfigyeléssel rendelkezünk. Az l. nap alatt egészen pontosan az ((l − 1) h, lh] intervallumot értjük, erre a továbbiakban, mint I_{l} is hivatkozunk majd. Jelölje az l. nap alatt bekövetkező károk számát ζ_{l}, ekkor ζ_{l} = N(((l − 1) h, lh]). Mi az 𝔼(ζ_{k+1}) illetve az inf{K ∈ ℕ_{0} : ℙ(ζ_{k+1} < K) ≥ 1 − α} értékeket akarjuk meghatározni az MMPP paramétereire kapott becslés alapján. Itt két különböző esetet fogunk vizsgálni.

Az egyik esetben feltesszük, hogy a biztosítónak nemcsak a napi kárszámok, hanem a pontos káridőpontok is rendelkezésére állnak az egész (0, kh] intervallumon. Ekkor az előző fejezetben ismertetett EM algoritmus segítségével megbecsüljük az MMPP paramétereit, majd az így kapott paraméterekkel kiszámoljuk a minket érdeklő előrejelzéseket.

4.1. ábra. Napi kárszámok egy MMPP egyik realizációja esetén

A másik esetben azt tesszük fel, hogy csak a napi (ζ_{1}, ζ_{2}, …, ζ_{k}) kárszámokból álló megfigyelés áll a rendelkezésünkre. Ekkor itt is az lenne a cél, hogy az ezen mintából számított egzakt MLbecslést felírva megbecsüljük az MMPP paramétereit, majd kiszámoljuk a megfelelő előrejelzéseket. A likelihood-függvény explicit alakban való megadásának bonyolultsága miatt azonban ehelyett egyszerűsített közelítő modelleket alkalmazunk, név szerint az alábbiakat tekintjük majd:

1. A napi kárszámok f.a.e. Poisson valószínűségi változók.

2. A napi kárszámok f.a.e. keverék Poisson valószínűségi változók.

A további célunk az lesz, hogy numerikus példákat megvizsgálva megállapítsuk, hogy egyrészt mennyire tér el a valóságtól a teljes trajektória ismeretében számított előrejelzés, másrészt, hogy mennyire tér el a valóságtól illetve a teljes trajektória ismeretében számított becsléstől az egyszerűsített és csupán napi kárszámokat felhasználó előrejelzés.

A továbbiakban feltesszük, hogy a kárfolyamat által meghatározott MMPP valós paraméterei Q és Λ, és, hogy a π kezdeti eloszlás a Q stacionárius eloszlása. Az előző fejezetben látottak alapján

𝔼(ζ_{k+1}) = 𝔼(N_{kh+h} − N_{kh}) = 𝔼(N_{kh+h}) − 𝔼(N_{kh})
= πλ(kh + h) − πλkh = πλh,

így az előrejelzésünk a kárszámra πλh, ha ismerjük a valós paramétereket.

Ha ζ tetszőleges diszkrét valószínűségi változó, akkor

ℙ(ζ < K) = ∑_{n=0}^{K−1} ℙ(ζ = n),

ezért az 1 − α megbízhatósági szinthez tartozó kvantilis kiszámítása az alábbi módon algoritmizálható.

4.2.1. Algoritmus. 1. Bemenet: ζ v.v., 1 − α.

2. Inicializálás: K := 1.

3. Ciklus: Amíg ℙ(ζ < K) < 1 − α

(a) K := K + 1

4. Kimenet: K kvantilis.

Látjuk, hogy az algoritmus egyetlen kérdéses pontja, hogy milyen könnyen tudjuk számolni a ℙ(ζ < K) mennyiségeket. Ha ezt ζ_{k+1}-re akarjuk alkalmazni, akkor tehát a ℙ(ζ_{k+1} < K) mennyiségek számítására kell hatékony módszert adnunk.

ℙ(ζ_{k+1} < K) = ∑_{n=0}^{K−1} ℙ(ζ_{k+1} = n) = ∑_{n=0}^{K−1} ℙ(N_{kh+h} − N_{kh} = n)
= ∑_{n=0}^{K−1} πP (n, h)1 = π (∑_{n=0}^{K−1} P (n, h)) 1,

tehát a célunk elérése érdekében az S(K, h) = ∑_{n=0}^{K−1} P (n, h) kiszámítása szükségeltetik. Erre a következő állítás szolgáltat megoldást.

4.2.1. Lemma. Tekintsük az alábbi nr × nr-es mátrixot:

C = ( \begin{matrix}Q − Λ & Λ & 0 & ⋯ & 0 \\ 0 & Q − Λ & Λ & ⋯ & 0 \\ 0 & 0 & Q − Λ & ... & 0 \\ ⋮ & ⋮ & ⋮ & ... & Λ \\ 0 & 0 & 0 & Q − Λ & Q − Λ \end{matrix} )

Legyen D(h) = exp(Ch), ekkor S(K, h) nem más, mint a D(h) mátrixból a legelső r sorának megtartásával és a többi sorának elhagyásával képzett E(h) mátrix, így S(K, h)1 éppen ezen mátrix sorösszegeiből képzett oszlopvektor.

Bizonyítás. A második fejezet MMPP-kről szóló alfejezete alapján tudjuk, hogy P (n, h), n = 0, 1, …, K − 1 kielégíti az alábbi differenciál-egyenletrendszert:

P (0, 0) = I és P (n, 0) = 0, ha n ≥ 1,
P′(0, t) = P (0, t)(Q − Λ),
P′(n, t) = P (n, t)(Q − Λ) + P (n − 1, t)Λ, ha n ≥ 1.

Emiatt a 2.3.8 lemma szerint P (n, h) = [D(h)]_{I_{0},I_{n}}, amiből már egyszerűen következik az állítás.

4.2.1. Megjegyzés. A mátrix exponenciális kiszámítása itt is a már korábban megismert Padé approximáció segítségével történik.

Ezzel tehát megmutattuk, hogy hogyan számolhatjuk ki a kérdéses előrejelzéseket abban az esetben, ha az összes érkezési időpont ismert. Most meghatározzuk ugyanezeket az egyszerűsített modellek esetén. Először tekintsük azt az esetet, amikor azt feltételezzük, hogy a ζ_{1}, ζ_{2}, …, ζ_{k} kárszámok f.a.e. exponenciálisak. Jelölje a paramétert μ. Az világos, hogy ekkor 𝔼(ζ_{k+1}) = μ. A kvantilis meghatározása a 4.2.1 algoritmus alapján semmilyen nehézséget nem okoz, tekintve, hogy a ℙ(ζ_{k+1} = n) = \frac{μ^{n}}{n!} e^{−μ} valószínűségek explicit alakban adottak. Már csak a paraméter becslése maradt hátra, ez azonban szintén könnyen elvégezhető, hiszen ekkor a likelihood-függvény

L(μ) = ∏_{i=1}^{k} \frac{μ^{ζ_{i}}}{ζ_{i}!} e^{−μ},

vagyis a loglikelihood-függvény

l(μ) = log μ ∑_{i=1}^{k} ζ_{i} − kμ − c = k(\overline{ζ} log μ − μ) − c,

ahol \overline{ζ} = \frac{1}{k} ∑_{i=1}^{k} ζ_{i} és c alkalmas ζ_{1}, ζ_{2}, …, ζ_{k}-tól függő konstans. Ennek a kifejezésnek pedig μ = \overline{ζ} esetén van maximuma, tehát a paraméter ML becslése nem más, mint μ̂ = \overline{ζ}.

Valamivel nehezebb dolgunk van, ha a ζ_{1}, ζ_{2}, …, ζ_{k} kárszámokat független, azonos eloszlásúaknak feltételezzük ugyan, de nem tisztán Poisson, hanem kevert Poisson eloszlást tételezünk fel. Itt kevert Poisson eloszlás alatt a következőt értjük.

4.2.1. Definíció. Legyen τ : Ω ↦ {μ_{1}, μ_{2}, …, μ_{r}} valószínűségi változó és vezessük be a p_{i} = ℙ(τ = μ_{i}) jelöléseket. Azt mondjuk, hogy a ζ valószínűségi változó a τ keverő valószínűségi változó (vagy a p keverő eloszlás) szerinti kevert Poisson eloszlás, ha

ℙ(ζ = i|τ) = \frac{τ^{i}}{i!} e^{−τ}.

Durván fogalmazva ez azt jelenti, hogy a ζ valószínűségi változó p_{i} valószínűséggel μ_{i} paraméterű Poisson eloszlású.

4.2.2. Megjegyzés. Természetesen nem szükséges, hogy τ csak véges sok értéket vegyen fel, lehet akár folytonos eloszlás is, mi azonban most csak ezt az esetet fogjuk használni.

Ebben az esetben tehát a paraméterek a μ = (μ_{1}, μ_{2}, …, μ_{r}) és a p = (p_{1}, p_{2}, …, p_{r}) vektorok. Először ezek segítségével megadjuk a ℙ(ζ_{k+1} = i) illetve az 𝔼(ζ_{k+1}) értékeket. A teljes valószínűség tételének értelmében

ℙ(ζ_{k+1} = i) = ∑_{j=1}^{r} ℙ(ζ = i|τ = μ_{j})ℙ(τ = μ_{j}) = ∑_{j=1}^{r} \frac{μ_{j}^{i}}{i!} e^{−μ_{j}}p_{j},

a várható érték pedig a teljes várható érték tételének figyelembe vételével

𝔼(ζ_{k+1}) = ∑_{j=1}^{r} 𝔼(ζ_{k+1}|τ = μ_{j})ℙ(τ = μ_{j}) = ∑_{j=1}^{r} p_{j}μ_{j}.

Tehát a kvantilisre és a kárszámra vonatkozó előrejelzéseinket most is könnyen el tudjuk végezni, ha a μ és p paraméterek becslése már ismert. A likelihood-függvény most

L(μ, p) = ∏_{i=1}^{k} (∑_{j=1}^{r} \frac{μ_{j}^{ζ_{i}}}{ζ_{i}!} e^{−μ_{j}}p_{j}),

vagyis a loglikelihood

l(μ, p) = ∑_{i=1}^{k} log (∑_{j=1}^{r} \frac{μ_{j}^{ζ_{i}}}{ζ_{i}!} e^{−μ_{j}}p_{j})

alakot ölt. Ez esetben a likelihood egyenletek

0 = \frac{∂}{∂μ_{s}} l(μ, p) = ∑_{i=1}^{k} (∑_{j=1}^{r} \frac{μ_{j}^{ζ_{i}}}{ζ_{i}!} e^{−μ_{j}}p_{j})^{−1} \frac{μ_{s}^{ζ_{i}}}{ζ_{i}!} e^{−μ_{s}}p_{s} (\frac{ζ_{i}}{μ_{s}} − 1) ,

0 = \frac{∂}{∂p_{s}} l(μ, p) = ∑_{i=1}^{k} (∑_{j=1}^{r} \frac{μ_{j}^{ζ_{i}}}{ζ_{i}!} e^{−μ_{j}}p_{j})^{−1} \frac{μ_{s}^{ζ_{i}}}{ζ_{i}!} e^{−μ_{s}},
0 = ∑_{i=1}^{r} p_{i} − 1, s = 1, …, r

alakúak, amivel nyilván ekvivalens a

0 = \frac{∂}{∂μ_{s}} l(μ, p) = ∑_{i=1}^{k} (∑_{j=1}^{r} \frac{μ_{j}^{ζ_{i}}}{ζ_{i}!} e^{−μ_{j}}p_{j})^{−1} \frac{μ_{s}^{ζ_{i}−1}}{(ζ_{i} − 1)!},

0 = \frac{∂}{∂p_{s}} l(μ, p) = ∑_{i=1}^{k} (∑_{j=1}^{r} \frac{μ_{j}^{ζ_{i}}}{ζ_{i}!} e^{−μ_{j}}p_{j})^{−1} \frac{μ_{s}^{ζ_{i}}}{ζ_{i}!} ,
0 = ∑_{i=1}^{r} p_{i} − 1, s = 1, …, r

rendszer. Ebben az esetben az ML becslés egzakt felírása tehát nem lehetséges, a megoldás olyan numerikus módszerek alkalmazása árán érhető csak el, amely túlmutat a jelen dolgozat keretein. Ehelyett momentum módszert alkalmazunk a paraméterek becslésére, vagyis felírunk annyi összefüggést a momentumok és a paraméterek között, amelyekből már egyértelműen kifejezhetőek a paraméterek, majd az így adódó kifejezésekben a momentumokat a tapasztalati megfelelőjükkel helyettesítjük. Mivel egy η paraméterű ζ Poisson-eloszlás generátorfüggvénye

G(z) = e^{η(z−1)},

azért az eloszlás f_{k} = 𝔼(ζ(ζ − 1) · . . . · (ζ − k + 1)) — az f_{0} = 1 megállapodással — faktoriális momentumai könnyen meghatározhatók deriválás segítségével

f_{k} = \frac{∂}{∂z^{k}} [G(z)]_{z=1} = η^{k}.

Kevés utánagondolással látható, hogy a mi esetünkben a momentum módszerrel felírt összefüggések ekvivalensek a faktoriális momentumokra felírt egyenletekkel. Az i-edik faktoriális momentum:

f_{i} = 𝔼(ζ_{k+1}(ζ_{k+1} − 1) · . . . · (ζ_{k+1} − i + 1))
= ∑_{j=1}^{r} 𝔼(ζ_{k+1}(ζ_{k+1} − 1) · . . . · (ζ_{k+1} − i + 1)|τ = μ_{j})ℙ(τ = μ_{j})
= ∑_{j=1}^{r} p_{j}μ_{j}^{i},

tehát a

∑_{j=1}^{r} p_{j}μ_{j}^{i} = f_{i}, i = 0, … 2r − 1

egyenletrendszert kell megoldanunk. Legyen ϑ = (μ_{1}, …, μ_{r}, p_{1}, …, p_{r}) és tegyük fel, hogy a megoldás ϑ = F (f_{0}, …, f_{2r−1}), ekkor a ϑ becslése

ϑ̂ = F (f̂_{0}, …, f̂_{2r−1}),

ahol f̂_{0}, …, f̂_{2r−1} az f_{0}, …, f_{2r−1} faktoriális momentumok ζ_{1}, … ζ_{k} függvényében számított tapasztalati becslése. Az alábbiakban megmutatjuk, hogy miként számolható a ϑ̂ becslés abban az esetben, amikor r = 2. Ekkor a következő egyenletrendszert kell megoldanunk:

p_{1} + p_{2} = 1
p_{1}μ_{1} + p_{2}μ_{2} = f̂_{1}
p_{1}μ_{1}^{2} + p_{2}μ_{2}^{2} = f̂_{2}
p_{1}μ_{1}^{3} + p_{2}μ_{2}^{3} = f̂_{3},

ahol

f̂_{1} = \frac{1}{k} ∑_{i=1}^{k} ζ_{i} = \overline{ζ},
f̂_{2} = \frac{1}{k} ∑_{i=1}^{k} ζ_{i}(ζ_{i} − 1) = \frac{1}{k} ∑_{i=1}^{k} ζ_{i}^{2} − \frac{1}{k} ∑_{i=1}^{k} ζ_{i},

és

f̂_{3} = \frac{1}{k} ∑_{i=1}^{k} ζ_{i}(ζ_{i} − 1)(ζ_{i} − 2) = \frac{1}{k} ∑_{i=1}^{k} ζ_{i}^{3} − \frac{3}{k} ∑_{i=1}^{k} ζ_{i}^{2} − \frac{1}{k} ∑_{i=1}^{k} ζ_{i}.

Az egyenletrendszer első két egyenletéből kapjuk, hogy

p_{1} = \frac{f̂_{1} − μ_{2}}{μ_{1} − μ_{2}} és p_{2} = \frac{μ_{1} − f̂_{1}}{μ_{1} − μ_{2}} .  (4.1)

Írjunk, mindenütt 1 − p_{1}-et p_{2} helyére, ekkor a rendszer így módosul:

p_{1}(μ_{1} − μ_{2}) + μ_{2} = f̂_{1}  (4.2)
p_{1}(μ_{1}^{2} − μ_{2}^{2}) + μ_{2}^{2} = f̂_{2}
p_{1}(μ_{1}^{3} − μ_{2}^{3}) + μ_{2}^{3} = f̂_{3},

A második egyenletbe p_{1} kifejezését behelyettesítve (f̂_{1} − μ_{2})(μ_{1} + μ_{2}) + μ_{2}^{2} = f̂_{2}, azaz

(f̂_{1} − μ_{2})μ_{1} = f̂_{2} − f̂_{1}μ_{2}  (4.3)

adódik. A harmadik egyenlet bal oldala pedig

p_{1}(μ_{1}^{3} − μ_{2}^{3}) + μ_{2}^{3} = p_{1}(μ_{1} − μ_{2})(μ_{1}^{2} + μ_{1}μ_{2} + μ_{2}^{2}) + μ_{2}^{3}
= (f̂_{1} − μ_{2})(μ_{1}^{2} + μ_{1}μ_{2} + μ_{2}^{2}) + μ_{2}^{3}
= (f̂_{1} − μ_{2})μ_{1}^{2} + (f̂_{2} − f̂_{1}μ_{2})μ_{2} + f̂_{1}μ_{2}^{2}
= (f̂_{1} − μ_{2})μ_{1}^{2} + f̂_{2}μ_{2},

vagyis

(f̂_{1} − μ_{2})μ_{1}^{2} + f̂_{2}μ_{2} = f̂_{3}.

Szorozva az egyenlet mindkét oldalát (f̂_{1} − μ_{2})-vel, azt kapjuk, hogy

(f̂_{2} − f̂_{1}μ_{2})^{2} + (f̂_{1} − μ_{2})f̂_{2}μ_{2} = f̂_{3}(f̂_{1} − μ_{2}),

majd ezt rendezve végül az

(f̂_{1}^{2} − f̂_{2})μ_{2}^{2} + (f̂_{3} − f̂_{1}f̂_{2})μ_{2} + f̂_{2}^{2} − f̂_{1}f̂_{3} = 0  (4.4)

másodfokú egyenlethez jutunk. Ha μ_{2} kielégíti ezt az egyenletet, akkor μ_{1} is, hiszen az (4.3) összefüggés miatt

μ_{1} + μ_{2} = \frac{f̂_{2} − f̂_{1}μ_{2}}{(f̂_{1} − μ_{2})} + μ_{2} = \frac{f̂_{2} − f̂_{1}μ_{2} + (f̂_{1} − μ_{2})μ_{2}}{(f̂_{1} − μ_{2})} = \frac{f̂_{2} − μ_{2}^{2}}{(f̂_{1} − μ_{2})} = −\frac{f̂_{3} − f̂_{1}f̂_{2}}{(f̂_{1}^{2} − f̂_{2})} ,

és

μ_{1}μ_{2} = \frac{(f̂_{2} − f̂_{1}μ_{2})μ_{2}}{(f̂_{1} − μ_{2})} = \frac{f̂_{2}^{2} − f̂_{1}f̂_{3}}{(f̂_{1}^{2} − f̂_{2})} ,

ahol az utolsó egyenlőség mindkét átalakítás-sorozat esetén csupán (4.4) átrendezése. Ez azt jelenti, hogy μ_{1} és μ_{2} kielégítik az

(f̂_{1}^{2} − f̂_{2})x^{2} + (f̂_{3} − f̂_{1}f̂_{2})x + f̂_{2}^{2} − f̂_{1}f̂_{3} = 0  (4.5)

másodfokú egyenlethez tartozó Viète-formulákat, vagyis a másik gyök csakis μ_{1} lehet. Összefoglalva tehát, μ_{1} és μ_{2} becslései a (4.5) másodfokú egyenlet két gyökeként adódnak, ami után már p_{1} és p_{2} a (4.1) kifejezések alapján egyértelműen meghatározott.

4.3. Numerikus eredmények

Numerikus példákon teszteltük az EM algoritmus, illetve a korábbi alfejezetben leírt becslési modellek alkalmazásának hatékonyságát. A most következő szakaszban ezen vizsgálódásaink eredményéről számolunk be.

Ryden a már korábban említett 1996-os cikkében két különböző 2 rendű MMPP esetén közölt numerikus eredményeket az általa prezentált algoritmus segítségével. Ezek rendre az általa 2A, illetve 2B-vel jelölt esetek. Az 4.1 táblázat mutatja az ezeknek megfelelő paramétereket. Ezek jól láthatóan úgy lettek megválasztva, hogy az egyik esetben λ_{1} és λ_{2} jól elkülöníthető legyen egymástól, míg a másikban kevésbé.

\begin{matrix}Eset & q_{12} & q_{21} & λ_{1} & λ_{2} \\ 2A & 10 & 1 & 100 & 10 \\ 2B & 5 & 2 & 100 & 50 \end{matrix}

4.1. táblázat. A Rydén által vizsgált 2A és 2B eset

Mi újraimplementáltuk az algoritmust a [21]-ben ajánlott javításoknak - és a 3.3 alfejezetben leírtaknak - megfelelően, tehát alkalmunk nyílik összevetni az eredményeket a Rydén által kapottakkal. Az összehasonlítás persze nem teljesen korrekt, mivel más trajektórákból készítettük a becslést és más kiindulási becslést használtunk. Mi különösen az előző alfejezetben felvetett problémára fektettünk hangsúlyt, ezért a mintagenerálásnál a következő módon jártunk el. A korábban tárgyalt azonos hosszúsági időperiódusok (amiket az egyszerűség kedvéért napoknak kereszteltünk el) hosszát egységnyinek választottuk, majd feltettük, hogy k = 365. Persze k egész más érték is lehetne - ügyelve arra, hogy a köztes időkből vett minta kellően nagy maradjon -, mi most önkényesen választottuk ezt az értéket, azt az elképzelt esetet szemléltetve, amikor a biztosítónak egy évnyi megfigyelése van napi kárszámokból (vagy a jobbik esetben a kárfolyamat teljes trajektóriájából). Ez az eljárás a 2A esetben 6828 dimenziós, míg a 2B esetben egy jóval nagyobb 23483 dimenziós köztes időkből vett mintát eredményezett. Rydén a becsléseit 5000 elemű mintából számította. A 4.2 táblázat az első néhány nap szimulált kárszámait mutatja.

Nap 1 2 3 4 5 6 7 8 9 10
2A 18 14 16 9 10 9 41 17 24 30
2B 68 72 56 71 77 72 65 56 47 74

4.2. táblázat. Szimulált kárszámok a 2A és 2B esetben

Bár mindkét esetben rosszabb kezdeti értékről indítottuk az EM algoritmust, mint Rydén, az általunk kapott becslés nagyjából ugyanolyan jól közelítette a valós paramétereket, sőt, némelyik paraméter esetén jobb is volt annál. A kapott becsléseket a 4.3 táblázatban foglaltuk össze.

A továbbiakban e két példán keresztül bemutatjuk, hogy milyen értékeket kaptunk az előző alfejezetben ismertetett előrejelzésekre a szintén ott tárgyalt teljes trajektóriából számoló (tehát a Rydén-féle EM algoritmust használó), illetve az egyszerűsített, f.a.e. Poisson illetve kevert Poisson kárszámot feltételező modellek esetén. Először térjünk rá az 𝔼(ζ_{k+1}) előrejelzésre (emlékezve, hogy esetünkben k = 365). Az ide vonatkozó eredmények az 4.4 táblázatban láthatók. A "Valós paraméter" fejlécű oszlopban az igazi paraméterek szerint számított érték van, míg az "MMPP", "Tiszta Poisson", illetve "Kevert Poisson" fejlécű oszlopokban rendre a már korábban említett módon becsült paraméterekből számított értékek állnak. Itt persze a gyakorlati életben nyilván csak a

Rydén becslése  Új becslés
Eset q_{12} q_{21} λ_{1} λ_{2} q_{12} q_{21} λ_{1} λ_{2}
2A 9.61 0.921 96.0 9.63 10.733 1.095 99.278 10.102
2B 5.94 2.26 105 50.8 4.930 2.148 99.457 49.041

4.3. táblázat. Rydén becslése és a mostani becslés

valamilyen módon egészekre kerekített értékek az érdekesek — tekintve, hogy a kárszám csak egész értékű lehet —, most azonban kiírtuk a tizedesjegyeket is, hogy jobban összehasonlíthatóak legyenek a különböző modellek.

Előrejelzés Valós paraméter MMPP Tiszta Poisson Kevert Poisson
2A  18.182  18.360 18.348  18.348
2B  64.286  64.341 64.337  64.337

4.4. táblázat. Az 𝔼(ζ_{k+1}) előrejelzések a három különböző modell esetén

Azt látjuk tehát, hogy a 3 különböző modell esetén a becsült előrejelzés csak minimálisan tér el a valós paraméterek ismeretében számítottól, a tisztán Poisson és kevert Poisson modellek esetén pedig teljesen megegyezik összhangban az előző alfejezet végén tett megjegyzésünkkel. Ebből az a sejtés adódik, hogy amennyiben az 𝔼(ζ_{k+1}) előrejelzésre vagyunk kíváncsiak, úgy elegendő csupán a napi kárszámok ismerete és az egyszerűsített Poisson modell használata.

Vizsgáljuk meg, hogy milyen becslések adódtak a kvantilisekre az egyes modellekből. Kétféle megbízhatósági szint - 0.95 és 0.99 - esetére határoztuk meg a becsléseket. Az eredmények az 4.5 és 4.6 táblákban láthatók.

Szint=0.95 Valós érték MMPP Tiszta Poisson Kevert Poisson
2A  42  42  27  43
2B  90  91  79  90

4.5. táblázat. A 0.95-os megbízhatósági szinthez tartozó kvantilisbecslések

Látható, hogy a teljes trajektória ismeretében számított becslés nagyon jól közelít. A tisztán Poisson modell mindkét megbízhatósági szint esetén teljesen rossz eredményt ad, ráadásul számunkra kedvezőtlen irányban, ugyanis alábecsüli a kvantilis értéket. Ezzel szemben a kevert Poisson modell 0.95-os megbízhatósági szinten meglepően jó eredményt ad. A magasabb megbízhatósági szinten azonban már ez a modell is jóval nagyobbat téved. Ez alapján az a sejtés fogalmazható meg, hogy a kvantilis becslése esetén az egyszerűsített Poisson modell alkalmazása teljesen elvethető, ugyanis az közel sem szolgáltatja az általunk elvárt eredményt. A kevert Poisson modell azonban alkalmas tűnik arra, hogy egyszerű eljárással viszonylag pontos értéket kapjunk vele.

Az imént említett 2 fő példán kívűl még további 52 kísérletet végeztünk. A valós paramétereket véletlenszerűen választottuk, ügyelve arra, hogy az így keletkező minta „élethű” legyen, azaz a napi kárszámok alakulása reális értéket tükrözzön. Az egyes futási eredmények számszerű értékeit a függelékben található táblázatokban foglaltuk össze. Az alábbiakban ismertetjük a tesztelés során összegyűlt tapasztalatainkat.

A Ryden-féle becslés pontosságáról azt tapasztaltuk, hogy az euklideszi norma szerint vett átlagos eltérés a valós illetve a becsült paraméterek között 3.04, ami 3.26%-a a valós paraméterek átlagos euklideszi normájának. A futások részeletes jegyzőkönyve megtalálható az A.1, A.2, A.3

Szint=0.99 Valós érték MMPP Tiszta Poisson Kevert Poisson
2A  58  56  30  49
2B  102  102  85  98

4.6. táblázat. A 0.99-os megbízhatósági szinthez tartozó kvantilisbecslések

és A.4 táblázatokban. Az 𝔼(ζ_{k+1}) értékekre kapott becslések esetén a közelítő Poisson modellből származó becslés egyetlen eset kivételével mindig jobb eredményt adott, mint a jóval kifinomultabb MMPP modellből a teljes trajektória ismeretében számított becslés. Ez utóbbi esetén az eltérések 0.02% és 7.07% százalék között mozogtak, az átlagos eltérés pedig 1,84% volt. Ugyanezen százalékok a Poisson modell esetén 0.000023%, 0.2% és 0.05%. A 7.07% százalékos eltérést produkáló MMPP futás esetén a Poisson futás eredményének eltérése csupán 0.02% volt, vagyis még mindig legfeljebb akkora, mint az MMPP legjobb esetében. A részletes futási eredményeket az A.5 és az A.6 táblázat tartalmazza. Összességében elmondható tehát, hogy a Poisson közelítés igen jól teljesít és a megvizsgált esetekben a teljes trajektóriából származó plusz információ sem tudott javítani az így kapott becslésen.

A 95%-os megbízhatósági szinthez tartozó kvantlisek vizsgálatakor az MMPP szerinti becslések legnagyobb tévedése 5 volt, de a legtöbb esetben pontosan eltalálták a valós értéket. A hibák a 4.7 táblázatban összefoglalt eloszlás szerint alakultak.

Eltérés  0  1  2  3  4  5
A futások számának %-ában 46.3% 31.5% 11.1% 7.4% 0% 3.7%

4.7. táblázat. Az MMPP becslésből számított kvantilisek eltéréseinek a gyakorisága

Elmondható tehát, hogy az esetek 96%-ában az így becsült kvantilisek legfeljebb 3-mal tértek el a tényleges kvantilistől. A kevert Poisson modellből számított kvantilisbecslések csupán 3 futás esetén volt az eltérés nagyobb, mint 5, ezek értéke azonban kiugróan magas volt: két ízben 7, illetve egyszer 14. Ebben az esetben szintén teljesült az, hogy 0 eltérésnek a legnagyobb a gyakorisága, összehasonlításképpen közöljük az 4.7 táblázatnak megfelelő összegzést a kevert Poisson esetben.

\begin{matrix}Eltérés & 0 & 1 & 2 & 3 & 4 & 5 & 5< \\ % & 38.9% & 25.9% & 16.7% & 5.6% & 5.6% & 1.9% & 5.9% \end{matrix}

4.8. táblázat. A kevert Poisson becslésből számított kvantilisek eltéréseinek a gyakorisága

Tehát a kevert Poisson modell pontossága elmarad a szofisztikáltabb MMPP modelltől, bár az esetek 87%-ában nem volt nagyobb a hiba, mint 3, ami azért még mindig kellően magas arány. Persze hátrányként jegyezhető meg, hogy azokban az esetekben, amikor ez a hiba meghaladta az 5-öt, a tévedés mértéke nagyon nagy volt. Ezt valamelyest enyhíti az a tény, hogy a kiugróan nagy hibák pozitív irányúak voltak, azaz nem becsülték alul a lehetséges legnagyobb kárt.

Ezekkel ellentétben - ahogy az a korábbi 2 fő példa esetén már körvonalazódott - a tisztán Poisson közelítés a megvizsgált esetek közül egyszer sem találta el pontosan a valós kvantilis értéket, és a hibák szórása is igen nagy volt, a legkisebb érték 2, a legnagyobb 39 volt. Ezen túlmenően az eltérés mindig negatív irányú volt, azaz alábecsültük a valós kvantilis értékeket, aminek a hatása végzetes lehet a biztosítótársaság számára, ha eszerint képzi a minimális szavatoló tőke szükségletét. Összehasonlításképpen megemlítjük, hogy itt csak az esetek 24%-ában volt a hiba mértéke legfeljebb 3, ami lényegesen elmarad az előbbi két modell eredményétől. Itt az átlagos eltérés 10.79 volt. A futások részletes jegyzőkönyvét az A.7 és az A.8 táblázat tartalmazza.

Végezetül közöljük ugyanezeket az adatokat a 99%-os megbízhatósági szint esetére is. Az MMPP modell esetén az eltérési gyakoriságok ebben az esetben a 4.9 táblázatban összefoglaltak szerint alakultak.

Eltérés 0  1  2  3  4  5  6  7
%  37% 37% 14.8% 3.7% 3.7% 1.9% 0% 1.9%

4.9. táblázat. MMPP eltérések gyakorisága 0.99-es mb. szinten

Ez már nem olyan jó eredmény, mint az alacsonyabb megbízhatósági szint esetén, de azért az esetek túlnyomó részében itt is csak 0 vagy 1 hiba volt. Legfeljebb 3 hiba volt az esetek 93%-ában. Sokkal szembetűnőbb a romlás a kevert Poisson modell eredményében, itt ugyanis meglehetősen nagy tévedések is előfordultak, a legrosszabb esetben 27-tel kisebb volt a kevert Poisson becslésből számított kvantilis, mint a valós érték. Itt csak egyetlen egy olyan eset volt, amikor a becsült kvantilis pontos volt és csak valamivel több, mint az esetek felében (52%) volt a tévedés mértéke legfeljebb 3. Szemléltetésképpen közöljük a gyakoriságok táblázatát.

\begin{matrix}Eltérés & 0 & 1 & 2 & 3 & 4 & 5 & 5< \\ % & 1.9% & 24.1% & 4.8% & 9.3% & 7.4% & 11.1% & 31.5% \end{matrix}

4.10. táblázat. Kevert Poisson gyakoriságok 0.99-es mb. szinten

A tisztán Poisson modell ismét nem került közel a valós értékekhez, az átlagos eltérés a valós és becsült értékek között 18.07 volt, a legkisebb eltérés 4, a legnagyobb 65. Ennél a legnagyobb hibát eredményező futásnál az MMPP becslés hibája 0 volt, míg a kevert Poisson-é a meglehetősen nagy 27 érték. Továbbra is fenáll, hogy a tisztán Poisson modell végig alábecsülte az értékeket. A további részletek megtalálhatók az A.9 és az A.10 táblázatokban.

Összességében elmondható, hogy a tapasztalatok alapján az 𝔼(ζ_{k+1}) becslései közül a legjobb eredményt az egyszerűsített Poisson modell alkalmazása adja, vagyis a napi kárszámok mintaátlaga. Ezzel szemben a Poisson modell alkalmazása téves eredményeket szolgáltat mindkét vizsgált megbízhatósági szint melletti kvantilisbecslésekre. A kevert Poisson becslés az esetek túlnyomó részében megfelelő hibahatáron belül mozgott a 95%-os megbízhatósági szinten, az egyenletesen jó közelítésben azonban elvétve akadt egy-két kiugróan rossz érték, ami miatt ezen modell alkalmazása magában rejt némi bizonytalansági tényezőt. Ráadásul a megbízhatósági szint növelésével jelentősen romlottak a hibastatisztikák. Tehát a kísérletek alapján azt sejtjük, hogy a kevert Poisson közelítés alkalmazásának van létjogosultsága az alacsonyabb megbízhatósági szint esetén, 0.99%-os megbízhatósági szinten azonban már csak korlátozottan javasolt a használata. Úgy tűnik, hogy a kvantilisbecslések már jóval érzékenyebbek arra, hogy az egész trajektóriából számítjuk a becslést, vagy csupán napi kárszámokra szorítkozva, tehát amennyiben rendelkezésünkre áll az egész trajektóriából - vagyis a köztes időkből - vett minta, úgy akkor kapjuk a legmegfelelőbb eredményt, ha az - egyszerűsített modellekhez képest jóval bonyolultabb - EM algoritmus által megbecsült paraméterekből számítjuk becsléseinket.

Záró megjegyzések

5.1. Összegzés

A dolgozat során ismertettük a Markov-modulált Poisson folyamat fő tulajdonságait. Áttekintettük a paramétereinek becslésére szolgáló Rydén által bevezetett EM algoritmust, majd módosítottuk és implementáltuk azt a Matlab programozási környezetben a W.J.J. Roberts et. al. által eszközölt új skálázási technikának és a mátrix exponenciális függvényekhez kapcsolódó integrálok kiértékelésében elért javításoknak megfelelően. Így már egy egyszerűen programozható, hatékony algoritmushoz jutottunk. Kiegészítettük az algoritmust egy k-means klaszterezésen alapuló, kiindulási becslés előállítására alkalmas eljárással. Ezután numerikus példákat vizsgáltunk egy biztosításmatematikai feladaton keresztül, melynek során C. Van Loan egyik lemmájának analógiájára konstruáltunk egy algoritmust az MMPP modell esetén egy rögzített hosszúságú intervallumba eső megfigyelések számának az eloszlása által meghatározott, különböző megbízhatósági szintekhez tartozó kvantilis értékek numerikus számítására.

5.2. További nyitott kérdések

Bár a numerikus példák alapján a Rydén-féle EM algoritmus efficiens eszköznek bizonyult az MMPP-k becslésére, továbbra is nyitva maradt kérdés az algoritmus konvergencia-tulajdonságainak egzakt vizsgálata, mellyel se [23], se [21] nem foglalkozik. Annak ellenére, hogy az általunk implementált EM algoritmus működik tetszőleges rendű MMPP esetén, mi mégis csak 2-állapotú eseteket vizsgáltunk. További vizsgálódásokra volna szükség a magasabb dimenziós esetben is. A numerikus eseteken túlmenően szükség volna a 4. fejezet 2. alfejezetében leírt kvantilisbecslések korrekt matematikai érzékenységvizsgálatára, meg kellene vizsgálni, hogy mekkorát változik a kvantilis, a paraméterek kis változtatásának hatására. Amikor a napi kárszámokból álló (ζ_{1}, …, ζ_{k}) elemű mintából becsültünk, akkor közelítő módszereket alkalmaztunk. Ehelyett korrektebb lenne a (ζ_{1}, …, ζ_{k}) együttes sűrűségfüggvényének felírása, s abból a loglikelihood-becslés meghatározása. Ez azonban a következő problémát veti fel. Az együttes sűrűségfüggvény:

f_{(ζ_{1},...,ζ_{k})}(z_{1}, …, z_{k}) = ℙ(ζ_{1} = z_{1}, …, ζ_{k} = z_{k})
= 𝔼(ℙ(ζ_{1} = z_{1}, …, ζ_{k} = z_{k}|ρ_{u}, u ≥ 0)),

ahol ρ_{t} az MMPP intenzitásfolyamata. Egységnyi hosszúságú napokkal számolva ekkor ζ_{i} = N((i − 1, i]), ahol N jelöli a szóban forgó MMPP-t. Az előző átalakításban kihasználva a növekmények feltételes függetlenségét, azt kapjuk, hogy

f_{(ζ_{1},...,ζ_{k})}(z_{1}, …, z_{k}) = 𝔼 (∏_{j=1}^{k} ℙ(N((j − 1, j]) = z_{j}|ρ_{u}, u ≥ 0))

= 𝔼 (∏_{j=1}^{k} \frac{(ρ_{j} − ρ_{j−1})^{z_{j}}}{z_{j}!} e^{−(ρ_{j}−ρ_{j−1})})

= 𝔼 (∏_{j=1}^{k} \frac{(ρ_{j} − ρ_{j−1})^{z_{j}}}{z_{j}!} e^{−(ρ_{j}−ρ_{j−1})})

= 𝔼 (e^{−ρ_{k}} ∏_{j=1}^{k} \frac{(ρ_{j} − ρ_{j−1})^{z_{j}}}{z_{j}!} ) ,

ahol, mint tudjuk ρ_{j} − ρ_{j−1} = ∫_{j−1+}^{j} λ(ξ_{s})ds, ahol ξ_{s} a mögöttes Markov folyamat. Ahhoz, hogy meghatározzuk a fenti várható értéket, a (ρ_{1}, ρ_{2}, …, ρ_{k}) változók együttes eloszlását kellene felírni, ami nem tűnik túl könnyű feladatnak.

5.3. Köszönetnyilvánítás

Szeretném hálámat kifejezni témavezetőmnek, Arató Miklósnak, aki számomra e kutatási területet kijelölte és útmutatásával, illetve hasznos tanácsaival segített munkám sikeres elvégzésében. További köszönet illeti Fodor Erzsébetet, aki időt és fáradtságot nem kímélve gondosan átnézte a dolgozatomat - nagy szerepet vállalva ezzel a felmerülő nyelvtani hibák kiküszöbölésében -, valamint felbecsülhetetlen értékű számítógép-használati lehetőség biztosításával is hozzájárult a bemutatott futási eredmények elkészültéhez.

A. Függelék

A futási eredmények összefoglaló táblázatai

A.1. A valós paraméterek és az EM algoritmus által kapott

becslés

Valós paraméterek  Kiinduló becslés
Futás Mintaelemszám q_{12} q_{21} λ_{1} λ_{1} q_{12} q_{21} λ_{1} λ_{1}
2A  6828 10 1 100 10 5.641 3.407 38.287 4.880
2B  23483 5 2 100 50 17.699 17.088 107.055 21.921
1  8100 2.2 9.7 26.3 3.8 4.422 5.408 34.548 6.666
2  4766 10 1 50 10 4.158 3.335 23.160 4.294
3  5953 5 1 50 10 4.360 3.405 29.812 4.702
4  5488 5 1 40 10 3.813 3.365 26.159 4.504
5  5524 3 1 30 10 4.067 3.594 26.246 4.793
7  5934 2 1 30 10 3.140 3.370 26.691 4.511
8  6537 50 5 100 10 6.179 4.179 34.821 5.494
9  6371 55 5 100 10 6.290 4.207 34.036 5.627
10  6343 60 5 100 10 6.021 4.154 33.853 5.449
11  5946 65 5 100 10 5.085 4.020 29.886 5.191
12  5232 70 4 90 10 4.760 3.828 25.965 4.879
13  5476 60 4 90 10 5.839 3.868 29.340 5.092
14  5020 65 4 80 10 5.227 3.717 25.586 4.835
15  7096 60 5 80 15 6.724 5.285 34.796 6.927
16  7940 10 1 100 15 6.789 4.769 41.527 6.478
17  10431 10 1 100 20 7.888 6.571 50.927 8.860
18  7478 10 1 125 10 5.494 3.253 43.500 4.771
19  8053 10 1 150 10 7.281 3.456 50.903 5.163
20  8479 10 1 175 10 5.529 3.233 50.802 4.721
21  10232 10 1 200 10 6.314 3.254 64.513 4.895

A.1. táblázat. A becsülendő paraméterek és a kiinduló becslések

Valós paraméterek  Kiinduló becslés
Futás Mintaelemszám q_{12} q_{21} λ_{1} λ_{1} q_{12} q_{21} λ_{1} λ_{1}
22  9318 10 2 100 10 5.769 3.754 52.542 5.428
23  5499 15 1 100 10 5.461 3.451 30.048 4.795
24  5072 20 1 100 10 4.285 3.144 26.379 4.263
25  9452 5 1 100 10 4.445 3.332 49.220 4.797
26  10439 5 1 100 15 6.503 4.909 53.498 7.051
27  7688 5 1 75 10 4.502 3.354 40.639 4.678
28  10495 5 1 125 10 5.620 3.330 61.084 5.091
29  6039 10 1 90 10 4.965 3.342 32.660 4.597
30  6071 10 1 80 10 5.366 3.424 33.209 4.809
31  5635 10 1 70 10 4.698 3.403 29.230 4.612
32  5390 10 1 60 10 5.465 3.540 29.061 4.907
33  5103 30 2 150 10 4.935 3.302 27.013 4.494
34  5927 20 1 150 10 5.800 3.207 35.229 4.677
35  9080 25 2 145 15 7.689 5.241 48.364 7.129
36  9643 15 1 130 20 8.311 6.451 47.957 8.527
37  5853 15 3 70 5 4.638 2.689 35.391 3.629
38  9252 10 3 80 10 5.163 3.949 48.224 5.419
39  10609 15 3 100 15 7.974 5.659 56.416 7.710
40  10819 10 3 80 15 6.287 5.279 53.672 7.267
41  12366 10 5 70 15 6.865 6.454 58.884 8.677
42  13107 10 7 60 20 8.105 8.088 59.801 10.791
43  10754 20 6 60 20 9.341 7.561 52.550 9.916
44  9763 25 6 70 15 7.795 6.306 48.484 8.194
45  8971 30 6 70 15 7.831 6.125 44.815 7.994
46  8512 30 5 70 15 7.140 5.702 42.502 7.504
47  6125 35 5 65 10 5.763 4.241 31.444 5.450
48  7478 30 4 65 15 7.466 5.464 37.856 7.234
49  10550 20 4 75 20 8.152 7.046 50.771 9.254
50  10225 20 5 60 20 7.940 7.241 48.296 9.340
51  9541 20 3 65 20 7.991 6.963 45.252 8.880
52  10804 20 2 70 25 9.080 8.278 49.781 10.558

A.2. táblázat. A becsülendő paraméterek és a kiinduló becslések (folytatás)

Valós paraméterek  EM becslés
Futás q_{12} q_{21} λ_{1} λ_{1} q_{12} q_{21} λ_{1} λ_{1}
2A 10 1 100 10 10.733 1.095 99.278 10.102
2B 5 2 100 50 4.930 2.148 99.457 49.041
1 2.2 9.7 26.3 3.8 2.053 8.348 26.310 5.461
2 10 1 50 10 10.570 1.041 45.980 9.826
3 5 1 50 10 5.644 1.074 49.898 9.944
4 5 1 40 10 5.739 1.061 41.569 10.121
5 3 1 30 10 2.874 1.053 29.123 10.022
6 1 1 20 10 1.335 1.045 20.449 10.513
7 2 1 30 10 1.966 0.997 29.579 9.502
8 50 5 100 10 55.794 5.315 102.427 9.862
9 55 5 100 10 62.465 6.101 95.444 9.844
10 60 5 100 10 62.678 5.652 103.602 9.608
11 65 5 100 10 74.598 5.004 106.530 10.246
12 70 4 90 10 70.758 5.352 78.119 9.512
13 60 4 90 10 60.706 4.184 91.077 9.761
14 65 4 80 10 63.104 3.092 87.467 10.145
15 60 5 80 15 58.291 5.076 71.170 14.951
16 10 1 100 15 12.554 1.083 103.653 14.704
17 10 1 100 20 9.074 1.012 99.493 20.694
18 10 1 125 10 9.725 0.987 124.496 9.947
19 10 1 150 10 9.946 0.958 147.670 10.009
20 10 1 175 10 11.092 0.962 177.055 9.888
21 10 1 200 10 10.218 1.102 197.928 9.711
22 10 2 100 10 9.646 1.943 102.143 10.077
23 15 1 100 10 17.092 0.971 101.196 10.176
24 20 1 100 10 20.551 0.996 100.011 9.750
25 5 1 100 10 4.970 1.053 100.234 10.202
26 5 1 100 15 5.142 0.992 98.231 15.178
27 5 1 75 10 5.107 1.076 74.611 9.788
28 5 1 125 10 4.945 0.957 126.361 9.930
29 10 1 90 10 10.085 0.927 88.062 9.978
30 10 1 80 10 10.406 1.042 82.846 10.006
31 10 1 70 10 11.414 1.028 72.281 10.341
32 10 1 60 10 11.212 1.041 63.726 10.231
33 30 2 150 10 19.278 0.532 151.836 10.177
34 20 1 150 10 20.643 1.058 144.530 9.650
35 25 2 145 15 25.910 2.091 146.129 15.091
36 15 1 130 20 17.087 0.992 129.806 20.418
37 15 3 70 5 15.710 3.063 71.107 5.304
38 10 3 80 10 10.183 2.884 80.160 9.829
39 15 3 100 15 15.094 3.122 97.627 14.890
40 10 3 80 15 9.896 2.999 80.088 14.363

A.3. táblázat. A becsülendő paraméterek és az EM algoritmus által kapott becslés

Valós paraméterek  EM becslés
Futás q_{12} q_{21} λ_{1} λ_{1} q_{12} q_{21} λ_{1} λ_{1}
41 10 5 70 15 10.470 5.269 70.948 15.232
42 10 7 60 20 9.749 7.175 59.275 18.720
43 20 6 60 20 22.396 8.116 59.556 18.557
44 25 6 70 15 24.131 6.585 69.926 14.973
45 30 6 70 15 30.715 5.565 74.287 15.577
46 30 5 70 15 29.833 5.355 70.803 14.802
47 35 5 65 10 31.168 4.005 64.055 10.718
48 30 4 65 15 31.987 3.900 66.252 14.914
49 20 4 75 20 20.601 4.224 74.448 19.563
50 20 5 60 20 18.984 5.438 57.875 19.463
51 20 3 65 20 19.159 2.857 64.783 20.379
52 20 2 70 25 17.048 2.370 62.820 24.989

A.4. táblázat. A becsülendő paraméterek és az EM algoritmus által kapott becslés (folytatás)

A.2. A várható kárszám becslése

Becsült értékek

Futás Valós érték MMPP modell Poisson modell Kevert Poisson modell
2A  18.182  18.360  18.348  18.348
2B  64.286  64.341  64.337  64.337
1  22.140  22.195  22.192  22.192
2  13.636  13.066  13.058  13.058
3  16.667  16.330  16.310  16.310
4  15.000  15.030  15.036  15.036
5  15.000  15.144  15.134  15.134
6  15.000  14.874  14.863  14.863
7  16.667  16.259  16.258  16.258
8  18.182  17.913  17.910  17.910
9  17.500  17.461  17.455  17.455
10  16.923  17.383  17.378  17.378
11  16.429  16.298  16.290  16.290
12  14.324  14.337  14.334  14.334
13  15.000  15.003  15.003  15.003
14  14.058  13.757  13.753  13.753
15  20.000  19.454  19.441  19.441
16  22.727  21.770  21.753  21.753
17  27.273  28.599  28.578  28.578
18  20.455  20.498  20.488  20.488
19  22.727  22.099  22.063  22.063

A.5. táblázat. Az 𝔼(ζ_{k+1}) előrejelzések valós és becsült értékei

Becsült értékek

Futás Valós érték MMPP modell Poisson modell Kevert Poisson modell
20  25.000  23.234  23.230  23.230
21  27.273  28.041  28.033  28.033
22  25.000  25.512  25.529  25.529
23  15.625  15.069  15.066  15.066
24  14.286  13.924  13.896  13.896
25  25.000  25.939  25.893  25.893
26  29.167  28.607  28.600  28.600
27  20.833  21.071  21.063  21.063
28  29.167  28.806  28.753  28.753
29  17.273  16.551  16.545  16.545
30  16.364  16.638  16.633  16.633
31  15.455  15.457  15.438  15.438
32  14.545  14.775  14.767  14.767
33  13.415  13.984  13.981  13.981
34  16.667  16.225  16.238  16.238
35  24.630  24.879  24.877  24.877
36  26.875  26.419  26.419  26.419
37  15.833  16.040  16.036  16.036
38  26.154  25.350  25.348  25.348
39  29.167  29.070  29.066  29.066
40  30.000  29.648  29.641  29.641
41  33.333  33.884  33.879  33.879
42  36.471  35.913  35.910  35.910
43  29.231  29.463  29.463  29.463
44  25.645  26.754  26.748  26.748
45  24.167  24.583  24.578  24.578
46  22.857  23.325  23.321  23.321
47  16.875  16.791  16.781  16.781
48  20.882  20.494  20.488  20.488
49  29.167  28.902  28.904  28.904
50  28.000  28.016  28.014  28.014
51  25.870  26.142  26.140  26.140
52  29.091  29.605  29.600  29.600

A.6. táblázat. Az 𝔼(ζ_{k+1}) előrejelzések valós és becsült értékei (folytatás)

A.3. Kvantilisbecslések

Megbízhatósági szint: 95%  Becsült értékek
Futás  Valós érték MMPP modell Poisson modell Kevert Poisson modell
2A  42  42  27  43
2B  90  91  79  90
1  33  33  31  32
2  26  24  20  25
3  35  34  24  34
4  30  30  23  30
5  29  29  23  29
6  25  25  22  25
7  32  32  24  32
8  31  30  26  31
9  30  29  26  30
10  29  29  25  29
11  28  27  24  27
12  24  24  22  24
13  25  25  23  25
14  28  23  21  24
15  30  29  28  30
16  46  43  31  50
17  49  52  39  51
18  50  51  29  54
19  59  57  31  59
20  67  62  32  72
21  75  76  38  89
22  53  55  35  53
23  34  31  23  32
24  29  28  21  30
25  64  65  36  71
26  66  64  39  66
27  50  50  30  49
28  78  78  39  77
29  39  37  24  40
30  36  36  25  39
31  32  32  23  32
32  29  29  22  31
33  30  31  21  37
34  38  36  24  40
35  46  46  34  49
36  49  46  36  49
37  33  33  24  33

A.7. táblázat. A 95%-os megbízhatósági szinthez tartozó kvantilisbecslések

Megbízhatósági szint: 95%  Becsült értékek
Futás  Valós érték MMPP modell Poisson modell Kevert Poisson modell
38  50  49  35  48
39  52  51  39  51
40  53  53  40  50
41  54  54  45  52
42  53  52  47  51
43  43  43  40  44
44  40  41  37  40
45  37  38  34  38
46  36  36  33  38
47  28  28  25  26
48  33  32  29  33
49  44  44  39  43
50  41  41  38  41
51  39  39  36  40
52  42  42  40  42

A.8. táblázat. A 95%-os megbízhatósági szinthez tartozó kvantilisbecslések (folytatás)

Megbízhatósági szint: 99%  Becsült értékek
Futás  Valós érték MMPP modell Poisson modell Kevert Poisson modell
2A  58  56  30  49
2B  102  102  85  98
1  37  37  35  36
2  34  31  23  29
3  46  44  27  40
4  38  38  26  35
5  35  35  26  33
6  29  29  26  29
7  38  37  27  37
8  37  37  29  35
9  36  35  29  35
10  34  35  29  33
11  33  32  27  31
12  29  28  25  27
13  30  31  26  29
14  28  28  24  29
15  35  34  31  31
16  61  56  34  58
17  63  67  43  58
18  70  71  33  62
19  82  80  35  66
20  95  88  36  84

A.9. táblázat. A 99%-os megbízhatósági szinthez tartozó kvantilisbecslések

Megbízhatósági szint: 99%  Becsült értékek
Futás  Valós érték MMPP modell Poisson modell Kevert Poisson modell
21  107  107  42  80
22  68  70  39  59
23  45  42  26  38
24  38  38  24  35
25  85  86  39  80
26  86  84  43  73
27  65  65  33  50
28  105  106  43  85
29  53  51  28  47
30  48  49  28  45
31  43  42  26  38
32  38  38  25  36
33  43  45  24  45
34  52  50  27  47
35  58  58  38  56
36  63  59  40  56
37  43  43  27  37
38  62  61  39  54
39  64  63  43  57
40  64  64  44  56
41  63  63  49  58
42  60  59  52  57
43  49  49  44  50
44  47  48  40  45
45  44  44  38  43
46  42  43  36  43
47  34  34  28  30
48  38  37  33  39
49  52  51  43  48
50  47  47  42  46
51  45  46  40  46
52  48  48  44  47

A.10. táblázat. A 99%-os megbízhatósági szinthez tartozó kvantilisbecslések (folytatás)

Irodalomjegyzék

[1] F. Bodon. Adatbányászati algoritmusok. http://www.cs.bme.hu/ bodon/magyar/adatbanyaszat/ tanulmany/adatbanyaszat.pdf, 2006.

[2] E. Çinlar. Markov renewal theory. Adv. Appl. Prob., 1:123–187, 1969.

[3] D.R. Cox. Some statistical methods connected with series of events (with discussion). J. R. Statist. Soc., B 17(6):129., 1955.

[4] L. Czách. Mértékelmélet. kézirat, 2001.

[5] D. Daley, D.J. és Vere-Jones. An introduction to the theory of point processes. Springer Series in Statistics. New York: Springer.

[6] Wolfgang Fischer and Kathleen S. Meier-Hellstern. The markov-modulated poisson process (mmpp) cookbook. Perform. Eval., 18(2):149–171, 1993.

[7] Rice-J.A. Fredkin, D.R. Maximum likelihood estimation and identification directly from single-channel recordings. Proc. R. Soc. Lond., 249:125–132, 1992.

[8] L.A. Freed, D. S.-Shepp. A poisson process whose rate is a hidden markov process. Adv. Appl. Prob., 14:21–36, 1982.

[9] M. Laczkovich. Valós függvénytan. ELTE egyetemi jegyzet, Budapest, 1995.

[10] M.L. Leroux, B.G. Puterman. Maximum-penalized-likelihood estimation for independent and markovdependent mixture models. Biometrics, 48:545–558, 1992.

[11] Rabiner L.R. és Sondhi M.M. Levinson, S.E. An introduction to the application of the theory of probabilistic functions of a markov process in automatic speech recognition. Bell Syst. Tech. J., 62:1035–1074, 1983.

[12] P. Lévy. Systèmes semi-markovians à au plus une infinité dénombramble d’états possibles. Proc. Int. Congr. Math., 2:294., 1954.

[13] P. Lévy. Regenerative stochastic processes. Proc. Roy. Soc. (London), 232:6–31, 1955.

[14] Gy. Mihaletzky. Kockázati folyamatok. ELTE Eötvös kiadó, 2001.

[15] U. Naor, P. és Yechiali. Queueing problmes with heterogenous arrivals and service. Oper. Res., 19.

[16] M.F. Neuts. A queue subject to extroneous phase changes. Adv. Appl. Probab., 3.

[17] M.F. Neuts. Structured Stochastic Matrices of M/G/1-type and their Applications. Marcel Dekker, New York, NY, 1989.

[18] R. Pyke. Markov renewal processes: definitions and preliminary properties. Annals of Mathematical Statistics, 32:1231–1242, 1961.

[19] R. Pyke. Markov renewal processes with finitely many states. Annals of Mathematical Statistics, 32:1242–1259, 1961.

[20] R.-D. Reiss. A course on point processes. Springer-Verlag, 1993.

[21] Ephraim Y. és Dieguez E. Roberts, W.J.J. On rydén’s em algorithm for estimating mmpps. IEEE SIGNAL PROCESSING LETTERS, 13(6):373–376, 2006.

[22] T. Rydén. Parameter estimation for markov-modulated poisson processes. Stochastic Models, 10(4), 1994.

[23] T. Rydén. An em algorithm for estimation in markov-modulated poisson processes. Comput. Stat. Data Anal., 21(4):431–447, 1996.

[24] G.A.F. Seber. Multivariate Observations. Wiley, New York, 1984.

[25] H. Spath. Cluster Dissection and Analysis: Theory, FORTRAN Programs, Examples. Halsted Press, New York, 1985.

[26] Matlab szoftver. http://www.mathworks.com.

[27] L. Takács. Some investigations concerning recurrent stochastic processes of a certain type. Magyar Tud. Akad. Mat. Kutató Int. Közl., 3:115–128, 1954.

[28] S. Karlin H. M. Taylor. Sztochasztikus folyamatok. Gondolat, Budapest, 1985.

[29] C. Van Loan. Computing integrals involving the matrix exponential. IEEE transactions on automatic control, 23(3):395–404, 1978.

[30] C. Van Loan, C.-Moller. Nineteen dubious ways to compute the exponential of a matrix, twenty-five years later. SIAM Review, 45:3–49, 2003.

[31] D. Vere-Jones. The e-m algorithm for continuous-time point process examples. kézirat, 2005.

[32] C.F.J. Wu. On the convergence properties of the em algorithm. Ann. Statist., 11:95–103, 1992.