本記事では、M/M/c/cモデルとM/M/cモデル、そしてアーラン分布について紹介します。M/M/c/cモデルにはアーランB公式、M/M/cモデルにはアーランC公式という有名な公式が存在します。アーラン分布はこれらの証明に直接使われるわけではありませんが、指数分布に従う複数の独立な確率変数の和として自然に現れるため、M/M/cモデルのように、複数のサービスステージを経るような待ち行列システムにおいて、待ち時間や滞在時間の分布を分析する際に本質的な役割を果たします。

1 M/M/cモデルとM/M/c/cモデル

M/M/cモデルは、ランダムにやってくる客を\(c\)個の窓口で応対する状況で、待ち行列（キュー）を許容する場合のモデルです。以前の記事「M/M/1行列モデルの公式の導出」で紹介したM/M/1モデルで窓口の数を一般化したものに当たります。仮定を整理すると次の通りです。

• 単位時間当たりの来客数：ポアソン分布\(\mathrm{Po}(\lambda)\)に従う

• サービス時間：指数分布\(\mathrm{Ex}(\mu)\)に従う

• 窓口の数：\(c\)

• 待機列：無限に並べる

このモデルは、銀行の窓口、病院の受付など、客が待つことが可能な場面によく用いられます。また、処理する主体は必ずしも実際の「窓口」である必要はなく、単位時間あたりの到着数がポアソン分布に従うようなリクエストを複数のサーバーで処理する場合などにも適用できます。この場合はリクエストが「客」、サーバーが「窓口」になります。

M/M/cモデルは、確率過程としては、システム内に受け入れられている（待機中と処理されている途中の両方を含む）客の数を状態とする連続時間マルコフ連鎖\(\{X(t)\}_{t\in [0,\infty)}\)として記述されます。状態空間は\(S=\{0,1,2,3,\ldots\}\)で、状態\(i\)から状態\(j\)への遷移率\(q_{i,j}\)は \[q_{i,i+1}=\lambda,\quad q_{i,i-1}=
\begin{cases}
i\mu & (1 \leq i < c) \\
c\mu & (i \geq c)
\end{cases},\quad
q_{i, i} = – (q_{i, i+1} + q_{i, i-1})\] となります。遷移率の非対角成分は単位時間当たりに状態\(i\)から状態\(j\)へ移る確率を表し、対角成分は列の和が0になるように補正しています。M/M/cモデルは、微小時間で各状態から番号が1つ増えた状態と減った状態へ移る確率がそれぞれ決まっている「出生・死亡過程」（birth–death process）と呼ばれる確率過程の一種です。

一方、M/M/c/cモデルは、ランダムにやってくる客を\(c\)個の窓口で応対するという点はM/M/cモデルと同じですが、待ち行列（キュー）を想定しないモデルです。つまり、到着時に窓口が\(c\)個全て埋まっていた場合、客は待ち行列に並ばずにサービスを受けるのを諦めて帰ります。「M/M/c/c」の2つ目の\(c\)というのは、システムが同時に受け入れることができる客の最大数を表しています。仮定は次の通りです。

• 単位時間当たりの来客数：ポアソン分布\(\mathrm{Po}(\lambda)\)に従う

• サービス時間：指数分布\(\mathrm{Ex}(\mu)\)に従う

• 窓口の数：\(c\)

• 待機列：なし（待てない）

M/M/c/cモデルを連続時間マルコフ連鎖と見る場合、状態空間は受付中の客の数としてあり得るものがなす有限集合\(S=\{0,1,2,\ldots,c\}\)になります。遷移率は \[q_{i, i+1} =
\begin{cases}
\lambda (i < c) \\
0 & (i = c)
\end{cases},\] \[q_{i, i-1} = i\mu \quad (1 \leq i \leq c),\] \[q_{i, i} = – (q_{i, i+1} + q_{i, i-1}).\]

2 定常分布

状態空間上の確率分布は、状態\(i\)にある確率を\(p_i\)とすると、\(p=(p_0,p_1,p_2,\ldots)\)のように行ベクトルで表されます。分布の時間変化は、推移速度行列\(Q\)（状態\(i\)から状態\(j\)への遷移率を\((i,j)\)-成分にもつ行列で、無限次であることも許す）を用いて、微分方程式 \[\frac{d}{dt}p(t)=p(t)Q\] により記述されます。

証明は省略しますが、M/M/cモデル（\(\lambda<c\mu\)の場合に限る）やM/M/c/cモデルでは、マルコフ連鎖の性質から、分布（各状態とその状態にある確率との対応）は十分時間が経過したとする極限では定常分布に収束します。

\(\lambda<c\mu\)のとき、M/M/cモデルの定常分布\(\pi=(\pi_0,\pi_1,\pi_2,\ldots)\)は次のように求められます。\(\pi Q=0\)より、 \[0=-\lambda\pi_0+\mu \pi_1,\] \[0=\lambda \pi_0-(\lambda+\mu)\pi_1+2\mu\pi_2,\] \[0=\lambda \pi_1-(\lambda+2\mu)\pi_2+3\mu\pi_3,\] \[\vdots\] \[0=\lambda \pi_{c-2}-(\lambda+(c-1)\mu)\pi_{c-1}+c\mu\pi_c,\] \[0=\lambda \pi_{c-1}-(\lambda+c\mu)\pi_{c}+c\mu\pi_{c+1},\] \[0=\lambda \pi_{c}-(\lambda+c\mu)\pi_{c+1}+c\mu\pi_{c+2},\] \[\vdots\] ここから、以下の漸化式が従います。 \[0=-\lambda \pi_{i-1}+i\mu\pi_i\quad (i=1,2,\ldots,c),\] \[0=-\lambda \pi_{i-1}+c\mu\pi_i\quad (i=c+1,c+2,\ldots).\] この漸化式を解いて \[\pi_n=\pi_0\prod_{i=1}^n\frac{\lambda}{i\mu}=\frac{1}{n!}\left(\frac{\lambda}{\mu}\right)^n \pi_0\quad (n\leq c),\] \[\pi_n=\pi_0\prod_{i=1}^c\frac{\lambda}{i\mu}\cdot \prod_{i=c+1}^n\frac{\lambda}{c\mu}=\frac{1}{c!c^{n-c}}\left(\frac{\lambda}{\mu}\right)^n \pi_0\quad (n\geq c).\] \(\lambda/(c\mu)<1\)より、これらの和は収束して、 \[\sum_{n=0}^{\infty}\pi_{n}=\pi_0\cdot \left(\sum_{n=0}^{c-1}\frac{1}{n!}\left(\frac{\lambda}{\mu}\right)^n+\frac{(\lambda/\mu)^c}{c!}\cdot\frac{1}{1-\frac{\lambda}{c\mu}}\right)\] 確率の合計は1でなければならないから、 \[\pi_0=\frac{1}{\sum_{n=0}^{c-1}\frac{1}{n!}\left(\frac{\lambda}{\mu}\right)^n+\frac{(\lambda/\mu)^c}{c!}\cdot\frac{c\mu}{c\mu-\lambda}}\] これを漸化式の解に代入することで定常分布が求まりました。

ちなみに\(\displaystyle \frac{\lambda}{c\mu}\)はM/M/cモデルの利用率と呼ばれ、客が到着する頻度と処理される頻度の比を表します。これが1以上になると列がどんどん長くなって定常分布が存在しなくなるので、M/M/cモデルは基本的に\(\lambda<c\mu\)という制約の下で考えます。

M/M/c/cモデルは有限状態なので定常分布が存在し、同様に次のように求められます。\(c\)人より多い状態が存在せず、場合分けが発生しないため、より簡単になります。 \[\pi_i = \frac{1}{i!} \left( \frac{\lambda}{\mu} \right)^i \pi_0 \quad (0 \leq i \leq c),\] \[\pi_0 = \left( \sum_{i=0}^{c} \frac{1}{i!} \left( \frac{\lambda}{\mu} \right)^i \right)^{-1}\]

3 アーランB公式

アーランB公式は、呼損系において資源の利用要求（「呼」と言う）が発生した際に、系がその要求を受け入れられずにブロックする確率を計算する公式です。呼損系とはM/M/c/cモデルで想定されているシステムとほぼ同じで、基本的には複数の回線が用意されている電話が想定されます。\(c\)本の回線が用意されているとき、M/M/c/cモデルの仮定の下で、電話が通話中でつながらない確率（呼損率）は、 \[\pi_c=\frac{\frac{1}{c!} \left( \frac{\lambda}{\mu} \right)^c}{\sum_{i=0}^{c} \frac{1}{i!} \left( \frac{\lambda}{\mu} \right)^i }\] になります。ここで、\(A\coloneqq\lambda/\mu\)は一定時間に回線が使用される割合を表し、呼量と呼ばれます。呼量はアーランという単位で表されます。1アーラン、つまり\(A=1\)とは、平均的に1つの回線（窓口、資源）が継続的に利用されるくらいの利用度合いを意味し、2つの資源が50％ずつ利用される場合も1アーランになります。ポアソン分布の母数の意味と指数分布の期待値が何だったかを思い出すと、\(\lambda\)は単位時間に何回利用要求（呼）が発生するかを表し、\(1/\mu\)は1回の要求が処理される時間の平均（平均サービス時間、平均利用時間）を表していたのでした。よって、 \[\text{呼量}=\frac{\text{呼数}\times\text{平均利用時間}}{\text{対象時間}}\] が成り立ちます。例えば、1時間に10回の通話があり、平均4分の長さだった場合、呼量は\(10\times 4\div 60=2/3\)アーランになります。

ちなみに、この公式はIPAのネットワークスペシャリスト試験によく出題されます。

4 アーランC公式

アーランC公式は、\(\displaystyle \frac{\lambda}{c\mu}<1\)が満たされる場合のM/M/cモデルにおいて、到着時にすぐに処理されずに待たさせる確率\(P_{\mathrm{w}}\)を計算する公式です。待たされるのは窓口が\(c\)全て埋まっているときかつそのときに限るので、確率は\(\pi_n\) (\(n\geq c\))の合計であり、前述の定常分布の式から \[P_{\mathrm{w}}=\sum_{n=c}^{\infty}\pi_n=\frac{\frac{(\lambda/\mu)^c}{c!}\cdot\frac{c}{c-\lambda/\mu}}{\sum_{n=0}^{c-1}\frac{1}{n!}\left(\frac{\lambda}{\mu}\right)^n+\frac{(\lambda/\mu)^c}{c!}\cdot\frac{c}{c-\lambda/\mu}}\] となります。

5 M/M/cモデルにおける他の量を求める公式

M/M/cモデルにおいて、サービスを受けている最中の客を含まない待ち行列の長さの期待値は \[\begin{align*}
L &= \sum_{n=c}^{\infty} (n – c) \pi_n=\frac{(\lambda/\mu)^{c+1}}{c!c}\frac{\pi_0}{(1-\frac{\lambda}{c\mu})^2}
\end{align*}\] なお、M/M/cモデルでは、最初から窓口の数だけ列を作るのではなく、一列に並んで、処理されるときにどれかの窓口に割り当てられることを想定しています。

サービスを受けている最中の客を含む系内の客の総数の期待値は、定常分布を求める途中で出てきた漸化式を上手く使うことで次のように求めることができます： \[\begin{align*}
\sum_{n=0}^{\infty}n\pi_n = \sum_{n=0}^{c-1}n\pi_n +\sum_{n=c}^{\infty}((n-c)+c)\pi_n \\
&= \sum_{n=0}^{c-1}n\pi_n +\sum_{n=c}^{\infty}c\pi_n +L\\
= \sum_{n=1}^{c-1}\frac{\lambda}{\mu}\pi_{n-1} +\sum_{n=c}^{\infty}\frac{\lambda}{\mu}\pi_{n-1} +L\\
&= \frac{\lambda}{\mu}\sum_{n=0}^{\infty}\pi_n +L=\frac{\lambda}{\mu}+L.\\
\end{align*}\] 待ち時間、すなわちサービスが開始されるまでに要する時間の期待値は \[\begin{align*}
W = \sum_{n=c}^{\infty}\frac{1}{c\mu}(n-c+1)\pi_n \\
&= \sum_{n=c}^{\infty}(n-c+1)\frac{\pi_{n+1}}{\lambda}=\frac{L}{\lambda}.
\end{align*}\] これはリトルの法則と呼ばれます。この式で計算できるのは若干非自明かもしれませんが、系内の客の数が\(c\)未満のときは待ち時間0で、全ての窓口が埋まっているときは処理中の客がいずれか一人処理されるまでの平均時間が\(1/(c\mu)\)（\(c\)個の\(\mathrm{Ex}(\mu)\)に従う独立確率変数のmin）であることから分かります。固定された待ち人数においては、待ち時間は指数分布に従う確率変数の和であるため、「アーラン分布」と呼ばれる分布に従います。\(W\)はこれを\(\pi\)で重み付けたものの期待値になります。

6 アーラン分布

\(n\in\mathbb{N}\)と\(\mu>0\)に対し、確率密度関数が \[f(t)=\frac{\mu^n}{(n-1)!}t^{n-1}e^{-\mu t}\chi_{[0,\infty)}(t)\] で与えられる確率分布をアーラン分布\(\Gamma(n,\mu)\)と言います。ガンマ分布の形状母数が正の整数の場合に当たります。

この証明方法はいくつか考えられます。
証明1：

確率密度関数が \[f_n(t)=\frac{\mu^n}{(n-1)!}t^{n-1}e^{-\mu t}\quad (t\geq 0)\] であることを\(n\)に関する帰納法で証明する。\(n=1\)のときは指数分布そのものである。\(n=k\)のときに成り立つと仮定して、\(n=k+1\)のときの確率密度関数を畳み込みにより計算する。 \[\begin{align*}
f_{k+1}(t)&=\int_0^{t}f_{k} (s)f_1(t-s)ds \\
= \int_0^{t}\frac{\mu^k}{(k-1)!}s^{k-1}e^{-\mu s}\cdot\mu e^{-\mu(t-s)}ds \\
&= \frac{\mu^{k+1}}{(k-1)!}e^{-\mu t}\int_0^{t}s^{k-1}ds \\
= \frac{\mu^{k+1}}{k!}t^{k}e^{-\mu t}
\end{align*}\]
証明2：

\(\{T_1+\cdots +T_n\leq t\}\)は「単位時間に平均\(\mu\)回起きる偶発的な事象について、\(n\)回目の事象が発生した時刻が\(t\)以前である」という事象であり、これはすなわち「単位時間に平均\(\mu\)回起きる偶発的な事象が、時刻\(t\)までに少なくとも\(n\)回起きる」ということを意味する。よって、パラメータ\(\mu t\)のポアソン分布に従う離散確率変数\(N(t)\)を用いて \[\begin{eqnarray*}
P(T_1+\cdots +T_n\leq t)&= P(N(t)\geq n) \\
&= \sum_{k=n}^{\infty}\frac{(\mu t)^k e^{-\mu t}}{k!}\\
&= 1-\sum_{k=0}^{n-1}\frac{(\mu t)^k e^{-\mu t}}{k!}.
\end{eqnarray*}\] これを\(t\)で微分することで確率密度関数\(f_n\)が得られる。この読み替えの妥当性については「指数分布とポアソン分布」の記事を参照してください。
証明3：

指数分布の特性関数は \[\varphi_{T_i}(t)=(1-it/\mu)^{-1}\] よって、 \[\varphi_{T_1+\cdots+T_n}(t)=\prod_{i=1}^{n}\varphi_{T_i}(t)=(1-it/\mu)^{-n}\] これをフーリエ逆変換するとアーラン分布の確率密度関数が得られる。（アーラン分布の確率密度関数をフーリエ変換してこれに一致することを示しても良い）
証明4：

特性関数の代わりにモーメント母関数（フーリエ変換の代わりにラプラス変換）を使っても同様である。