M/M/1行列モデルの公式の導出

今回は待ち行列理論のM/M/1モデルの公式の数学的背景と導出について解説します。これは典型的には一列に並んだ客を一つの窓口が対応したり、リクエストの列を一台のサーバーが処理したりする状況で使われる公式です。

応用情報技術者試験でもよく問われるものですが、公式は基本的に暗記するものとされ、数学的理論について出題されることはありません。情報系の書籍でもさらっと公式だけ載っているケースが多いです。しかし、導出過程を知っておいた方が何かと気持ちよく覚えることができるでしょう。

1 M/M/1モデルとは

待ち行列理論とは、ざっくり言うと、列に並ぶときに平均でどれくらい待たされるかを予測する理論です。同様の状況であれば処理の内容や処理対象が何でも同じですが、以下では分かりやすいように客が窓口の前に並ぶケースを考えます。待ち行列理論では、待ち時間に影響する要素として主に以下の3つを想定します（応用として、待ち行列の容量の上限など4つ以上の要素を考えることもあります）。

単位時間当たりに何人の客が来るか
サービス時間：一人に対応するのにどれくらいの時間がかかるか
窓口（サービスチャネル）の数：一つの列に対していくつの窓口があるか

待ち行列理論にはいくつかの種類があり、種類を識別するために上記のそれぞれの要素を記号で表現してスラッシュで区切った表記（ケンドールの記号）がよく用いられます。中でもM/M/1モデルは代表的なもので、単位時間当たりの来客数が（一定ではなく）ポアソン分布に従い、サービス時間がそれとは独立な指数分布に従い、窓口の数が1である状況のモデルです。「M」というのはMarkovianの略で、マルコフ過程であることを表しています。

2 ポアソン分布と指数分布

母数\(\lambda\)のポアソン分布\(\mathrm{Po}(\lambda)\)は、大雑把に言うと「単位時間に平均では\(\lambda\)回発生するようなランダムな出来事が単位時間に何回発生するか」（\(X\)回とする）が従う確率分布で、確率質量関数は \[P(X=k)=\frac{\lambda^k}{k!}e^{-\lambda}\] で表されます。ポアソン分布の期待値（平均）は\(\lambda\)であることが容易に分かります。ここで言う「ランダムな出来事」というのは単に発生回数が確率変数であるという意味ではなく、隕石の落下回数や電話の本数などのような定常的で一定の起こりやすさをもつ無記憶で偶発的な出来事を意味します。ポアソン分布がそのような出来事の発生回数のモデル化として妥当であることの数学的定式化の一つは次の命題です。

Proposition 1. \(\{N(t)\}_{t\geq 0}\)を非負整数全体を状態空間とする確率過程とする。任意の\(t,s>0\)に対して\(N(t+s)-N(s)\)と\(N(s)\)が独立であるとする。このとき、\(N(t+s)-N(s)\sim \mathrm{Po}(\lambda t)\) (\(\forall t,s>0\))であることと、以下の2条件が成立することは同値である。

\(P(N(t+h)-N(t)=1)=\lambda h+O(h^2)\) as \(h\to +0\) (\(\forall t>0\)),
\(P(N(t+h)-N(t)\geq 2)=O(h^2)\) as \(h\to +0\) (\(\forall t>0\)).

1つ目の条件は生起率が時間によらず一定値\(\lambda\)であること、2つ目の条件は微小時間では2回以上発生する確率が無視できる、すなわち偶発的で多重発生が起きにくいということを表しています。証明は別の記事「指数分布とポアソン分布」で扱います。

一方、母数\(\mu\)の指数分布\(\mathrm{Ex}(\mu)\)は、「平均\(1/\mu\)の時間間隔で（つまり、平均的に単位時間当たり\(\mu\)回）起こるようなランダムな出来事」の発生する時間間隔が従う確率分布で、確率密度関数は \[f(x)=\mu e^{-\mu x}\] で表されます。発生回数がポアソン過程に従うような出来事の発生時間間隔は指数分布に従います。数学的にまとめると次の命題になります。

Proposition 2. \(\{N(t)\}_{t\geq 0}\)を生起率\(\lambda>0\)のポアソン過程とする。 \[T_j=\inf\{t\geq 0\mid N(t)=j\}\] に対し、\(T_1,T_2-T_1,T_3-T_2,T_4-T_3,\ldots\)は独立同分布で\(\mathrm{Ex}(\lambda)\)に従う。

逆に、発生間隔が独立に指数分布に従う場合、一定時間に発生する回数はポアソン分布に従います。これらの証明も「指数分布とポアソン分布」で扱います。

つまり、前述の意味で偶発的な出来事について、単位時間当たりの発生回数はポアソン分布に従い、発生時間間隔は指数分布に従うということです。

3 公式

M/M/1モデルを考えます。平均到着率（単位時間当たりの来客数の期待値）を\(\lambda\)、平均サービス率（サービス時間の期待値の逆数）を\(\mu\)とおきます。\(1/\lambda\)は平均到着間隔、\(1/\mu\)は平均サービス時間と呼ばれます。さらに、平均利用率と呼ばれる量を \[\rho\coloneqq \frac{\lambda}{\mu}\] で定義します。\(\rho\geq 1\)の場合、すなわち窓口の処理が客がやってくるスピードよりも平均的に遅い場合は、列がどんどん長くなり、定常状態では待ち時間が無限大に発散してしまうので、以下では\(\rho<1\)とします。このとき、平均待ち時間（待ち時間の期待値）は \[\frac{\rho}{1-\rho}\cdot \frac{1}{\mu}\] となります。

4 導出

以下、上記の公式を証明します。マルコフ過程のエルゴード性についての議論は厳密な証明をするためには必要ですが、どのようにして上記のような式が出てくるのかを理解する上では必要ないので、難しければ読み飛ばしてください。

Proof. M/M/1モデルの仮定から、単位時間に増える客の人数は平均\(\lambda\)のポアソン分布、単位時間に窓口が処理する人数は平均\(\mu\)のポアソン分布に従う。よって、微小時間\(\Delta t\)の間に客が新たに1人並ぶ確率は\(\lambda \Delta t+O(\Delta t^2)\)、客が1人処理される確率は\(\mu \Delta t+O(\Delta t^2)\)となる。従って、時刻\(t\)に並んでいる人数がちょうど\(n\)人である確率を\(p_n(t)\)とおくと、\(n\geq 1\)のとき、 \[\begin{align}
p_n(t+\Delta t) &= p_{n-1}(t)\cdot \lambda \Delta t+p_n (t)\cdot (1-\lambda \Delta t+O(\Delta t^2))(1-\mu \Delta t+O(\Delta t^2)) \\
&&\quad +p_{n+1}(t)\cdot \mu \Delta t+O(\Delta t^2) \\
&= p_{n-1}(t)\cdot \lambda \Delta t+p_n (t)\cdot (1-\lambda \Delta t-\mu \Delta t)+p_{n+1}(t)\cdot \mu \Delta t+O(\Delta t^2),\\
p_0 (t+\Delta t) &=p_0(t)\cdot(1-\lambda\Delta t)+p_1(t)\cdot \mu\Delta t+O(\Delta t^2).
\end{align}\] \(\{p_n(t)\}\)の定常過程\(\pi_n\)においては時間微分が0になる（実は微分を持ち出さなくても、時間に依存しない漸化式の解として求まる）から、 \[\lambda \pi_{n-1}-(\lambda+\mu)\pi_n+\mu \pi_{n+1}=0,\] \[\lambda\pi_0=\mu\pi_1\] すなわち \[\pi_{n+1}-(\rho+1)\pi_n+\rho \pi_{n-1}=0, \quad \pi_1=\rho \pi_0.\] これを隣接3項間漸化式として解く。特性方程式の解は\(1,\rho\)なので、\(n\)に依存しないある定数\(C_1,C_2\)が存在して、 \[\pi_n=C_1 +C_2 \rho^n\] となる。全確率の条件 \[\sum_{n=0}^\infty \pi_n=1\] より、 \[C_1=0,\quad C_2=1-\rho.\] よって、 \[\pi_n=(1-\rho)\rho^n.\] すなわち、定常分布\(\{\pi_n\}\)において、待ち行列の長さ\(n\)は母数\(1-\rho\)の幾何分布に従う（これは必要条件から求めたものだが、これが実際に定常分布になっていることは直接計算で確認できる）。また、\(\{p_n(t)\}\)はエルゴード的な連続時間マルコフ連鎖だから、その定常分布は極限分布でもある。よって、平均待ち時間は、母数\(1-\rho\)の幾何分布の期待値（待ち行列の長さの平均）に平均サービス時間（待ち行列中の1人当たりの処理時間）を掛けたものを計算して、 \[\sum_{n=1}^\infty n \pi_n\cdot \frac{1}{\mu}=\frac{\rho}{1-\rho}\cdot \frac{1}{\mu}\] と求まる。なお、ここで、サービス時間と待ち行列の長さが独立であることを用いた。 ◻

Remark 1. 証明の途中で出てきた漸化式を差分商の形に変形し、微分可能性を仮定して\(\Delta t\to 0\)とすると、微分漸化式 \[\frac{dp_n}{dt}=\lambda p_{n-1}(t)-(\lambda+\mu)p_n(t)+\mu p_{n+1}(t),\] \[\frac{dp_0}{dt}=-\lambda p_0(t)+\mu p_1(t)\] が得られ、ここから更に極限分布への漸近の仕方に関する情報を引き出すこともできます。

Remark 2. \(\rho\)に「利用率」という名前が付いているのは、定常分布が幾何分布であることから分かるように、定常状態で客が訪れた際に窓口が利用されて塞がっている（既に他の人が一人以上並んでいる）確率が\(\rho\)であることに由来しています。