ギルバートハウス › ポアソン過程を見る

ポアソン過程

標準ポアソン過程の標本路

 ランダム生成した標本路を表示しています。ページ読み込みのたびに再計算されるので、グラフの見た目が変化します。

コード(生成箇所)

const deltaT = 1.0 / 256.0;
let p1 = 0, p2 = 0, p3 = 0, p4 = 0;
poissonProcessPathsData.addRow([0.0, p1, p2, p3, p4]);
/* while t <= 5.0 */
for (let i = 1; i <= 1280; i++) {
    p1 += rbernoulli(deltaT);
    p2 += rbernoulli(deltaT);
    p3 += rbernoulli(deltaT);
    p4 += rbernoulli(deltaT);
    poissonProcessPathsData.addRow([i * deltaT, p1, p2, p3, p4]);
}

二項過程とポアソン過程

  • 二項過程は、ステップごとの増加がベルヌーイ分布にしたがう離散時間・離散状態の確率過程である。
    • $X_0 = 0$
    • $B_i \sim \text{Bernoulli}(p)$として、$X_i = X_{i - 1} + B_i$
  • 二項過程で
    • 1ステップの経過時間を$\Delta t$
    • ステップごとの増加の成功確率を$\Delta t\,\lambda$(1未満)
    とすると、時間あたりの平均増加数は$\Delta t$によらず$\lambda$(一定)に保たれる。数学的な厳密性を抜きにすれば、$\Delta t \to 0$で連続時間・離散状態の確率過程となり、これがポアソン過程である。$\lambda$はポアソン過程のパラメータで、強度または到着率と呼ばれる。$\lambda = 1$なら標準ポアソン過程である。
  • $\{X(t),\ t \in [0, +\infty)\}$が強度$\lambda$のポアソン過程にしたがうとき、
    • 任意の$t_1, t_2$(ただし、$0 \leq t_1 \lt t_2$)に対して、$X(t_2) - X(t_1)$は$\text{Poisson}((t_2 - t_1)\,\lambda)$にしたがう。〔定常増分性(ポアソン分布)〕
    • 任意の$t_1, \ldots, t_n$(ただし、$0 \leq t_1 \lt \cdots \lt t_n$)に対して、$X(t_1), X(t_2) - X(t_1), \ldots, X(t_n) - X(t_{n - 1})$は独立である。〔独立増分性〕

ポアソン乱数の生成

  • $\{X(t),\ t \in [0, +\infty)\}$が強度$\lambda$のポアソン過程にしたがうとき、到着時刻の間隔は$\text{Exp}(\lambda)$にしたがう。また、時刻$t$にかかわらず次の到着までにかかる時間も$\text{Exp}(\lambda)$にしたがう。1
    $i$番目の到着時刻
    $T_i = \min\{t \mid X(t) = i\}$
    到着時刻の間隔
    $D_i = T_i - T_{i - 1}$
  • 逆方向、つまり到着時刻の間隔が同一の$\text{Exp}(\lambda)$に独立にしたがうとき、$X(t)$を時刻$t \in [0, +\infty)$までの累積到着数とすると、$\{X(t)\}$は強度$\lambda$のポアソン過程にしたがう。1
    $i$番目の到着時刻
    $T_i = T_{i - 1} + D_i$(ただし、$T_0 = 0$)
    累積到着数
    $X(t) = \max\{i \mid T_i \leq t\}$
  • このことを利用して、以下の手順でポアソン乱数を得ることができる。
    1. 標準指数乱数$E_i$を順次生成する。
    2. $\frac{E_1 + \cdots + E_n}{\lambda}$が1以上になった時点で止める。
    3. $n - 1$を乱数の値とする。
    前節のポアソン過程とポアソン分布の関係で、$t_1 = 0,\ t_2 = 1$の場合にあたる。
  • Knuthのアルゴリズム2,3は上の手順を効率化したものである。
    1. 標準一様乱数$U_i$を順次生成する。
    2. $U_1 \cdots U_n$が$e^{-\lambda}$以下になった時点で止める。
      〈注意〉$U_i$の値が0をとる可能性がある場合は、$(1 - U_1) \cdots (1 - U_n)$のようにする。3
    3. $n - 1$を乱数の値とする。
    効率化前、Knuthのどちらの手順でも、繰り返しの期待回数は$\lambda + 1$となる。そのため、計算量は$\lambda$にほぼ比例する。
function rpoisson(lambda) {
    let L = Math.exp(-lambda);
    let n = 0;
    let product = 1.0;
    do {
        n += 1;
        product *= 1.0 - Math.random();
    } while (product > L);
    return n - 1;
}
  1. Kyle Siegrist. "Introduction – The Poisson Process". Random.
  2. Donald E. Knuth. The Art of Computer Programming, Volume 2: Seminumerical Algorithms. 3rd ed., Addison–Wesley, 1997.
  3. JIS Z 9031:2012. 乱数生成及びランダム化の手順.