Slovník | Vyhledávání | Mapa webu
 
Analýza a modelování dynamických biologických datVybrané kapitoly z matematického modelování Zobecnění Poissonova procesu Simulace nehomogenního Poissonova procesu

Logo Matematická biologie

Simulace nehomogenního Poissonova procesu

V této části si ukážeme praktický způsob, jak generovat trajektorii Poissonova procesu , s danou intenzitou měnící se v čase. Protože se intenzita v čase mění, nelze přímo použít způsob, kdy se vygenerují délky intervalů z exponenciálního rozdělení a jejich kumulativní součty pak určují okamžiky událostí v Poissonově procesu. Jednou z možností je využití definičních vlastností (i*) a (ii*). Na časové ose postupujeme po velmi malých krocích délky . V každém kroku, v čase , je pravděpodobnost příchodu události rovna a považujeme ji za konstantní během krátké doby . Vygenerujeme realizaci z alternativního rozdělení s pravděpodobností úspěchu ,
a v případě nastání úspěchu (tj. realizace je rovna 1) zapíšeme nový okamžik příchodu události. Tento postup je ale zdlouhavý a nepřesný, z důvodu nutnosti velmi jemného dělení časové osy.

V praxi se používá jiný, přesnější a svým přístupem zajímavý přístup k simulaci trajektorie nehomogenního Poissonova procesu. Pro proměnnou intenzitu určíme maximální hodnotu funkce na zkoumaném intervalu a označíme ji ; např. pro periodickou intenzitu (4) je . Následně, způsobem popsaným v části o homogenním Poissonově procesu, vygenerujeme trajektorii Poissonova procesu pro konstantní intenzitu rovnou . Označme takto vygenerované okamžiky příchodů událostí . Tyto časy považujeme za potenciální okamžiky událostí v nehomogenním procesu, přičemž o tom, jestli v čase skutečně nastala událost, rozhodneme na základě realizace alternativní náhodné veličiny s pravděpodobností úspěchu rovnou podílu . Celý algoritmus lze shrnout následovně: 

  • nalezneme
     
  •  vygenerujeme vzorek okamžiků příchodů událostí v homogenním Poissonově procesu s intenzitou
     
  •  pro každý okamžik , , vygenerujeme realizaci z alternativního rozdělení pravděpodobnosti s pravděpodobností úspěchu
     
  •  výsledný vzorek okamžiků příchodů událostí nehomogenního Poissonova procesu sestává právě z těch časů , pro něž jsou odpovídající realizace

Nyní ukážeme, jak trajektorii vygenerovat ve statistickém softwaru R [3] pro periodickou intenzitu (4). Využijeme přitom postup z kapitoly Statistické metody pro Poissonův proces. Po úvodním uložení hodnot do proměnných spočítáme maximální intenzitu Lambda a vygenerujeme okamžiky S homogenního Poissonova procesu s touto intenzitou. Spočítáme hodnoty proměnné intenzity lambda v těchto časech S a následně pomocí podílu pravděpodobnosti p. Ty využijeme jako parametry funkce rbinom pro generování náhodného výběru Z z alternativního rozdělení (které je speciálním případem binomického s parametrem size=1). Jako okamžiky událostí v nehomogenním Poissonově procesu Sn vybereme pomocí indexování S[ ] a funkce which právě ty z časů S, pro něž je odpovídající Z rovno jedné (Z == 1). Dále již analogicky jako u homogenního Poissonova procesu doplníme hodnoty čítacího procesu N (od nuly) a vytvoříme schodovitou funkci X odpovídající trajektorii procesu. Graf takto vygenerované trajektorie nehomogenního Poissonova procesu je na Obr.1.

K <- 150
a <- 1.5

b <- 1.5
omega <- 2 * pi / 24
phi <- -pi / 2
Lambda <- a + b
T <- rexp (K, rate=Lambda)
S <- cumsum (T)
lambda <- a + b * sin (phi + omega * S)
p <- lambda / Lambda
Z <- rbinom (K, size=1, prob=p)
Sn <- S[ which (Z == 1) ]
N <- seq (0, length (Sn), by=1)
X <- stepfun (Sn, N)
plot (X, verticals=FALSE, pch=20, col=2, lwd=1.5,
xlab="t", ylab="N(t)", main="")

 
vytvořil Institut biostatistiky a analýz Masarykovy univerzity | | zpětné odkazy | validní XHTML 1.0 Strict