あみだくじの確率を計算してみた

あみだくじの確率には,かなりの偏りがある。

あみだくじの確率を計算しようとしたら,思ったよりも面白い数学の話になったので紹介します。

問題設定

縦線が m+1m+1 本,横棒が nn 本であるようなあみだくじを考えます。

あみだくじ

nn 本の横棒は「ランダムに」引かれています。ここで言うランダムとは,各々の横棒が mm 箇所のどこに引かれる確率も 1m\dfrac{1}{m} であるという意味です。実現されうるあみだくじのパターンは mnm^n 通りあります。

このとき,aa 番からスタートして bb 番に到達する確率 Pm,n(a,b)P_{m,n}(a,b) がどうなるのかという問題を考えます(左端を 00 番,右端を mm 番とします)。

実験と考察

実際に確率を計算する部分はけっこう大変な数学になるので後回しにして,まずは m=5m=5(縦線が6本)の場合の計算結果を紹介します。

n=10n=10 の場合〜

下の表は縦がスタート地点の場所,横がゴール地点の場所,数値が確率を表しています。例えば二行目の三列目が 0.2030.203 なので,左から二番目からスタートすると左から三番目にゴールする確率が 20.320.3 %くらいということです。

(0.3740.2970.1870.0930.0360.0130.2970.2640.2030.1310.0700.0360.1870.2030.2070.1800.1310.0930.0930.1310.1800.2070.2030.1870.0360.0700.1310.2030.2640.2970.0130.0360.0930.1870.2970.374)\begin{pmatrix} 0.374 & 0.297 & 0.187 & 0.093 & 0.036 & 0.013\\ 0.297& 0.264 & 0.203 & 0.131 & 0.070 & 0.036\\ 0.187 & 0.203 & 0.207 & 0.180 & 0.131 & 0.093\\ 0.093 & 0.131 & 0.180 & 0.207 & 0.203 & 0.187\\ 0.036 & 0.070 & 0.131 & 0.203 & 0.264 & 0.297\\ 0.013 & 0.036 & 0.093 & 0.187 & 0.297 & 0.374 \end{pmatrix}

  • 端からスタートすると同じ場所に到着する確率が 3737 %以上もあり,反対側の端に到着する確率は 11 %程度しかありません。確率がかなり偏っています!
  • 対角成分が大きい,つまりスタート地点と同じ場所にゴールする確率が高い傾向があります。

n=50n=50 の場合〜

(0.1870.1810.1720.1610.1520.1470.1810.1770.1710.1630.1560.1520.1720.1710.1680.1650.1630.1610.1610.1630.1650.1680.1710.1720.1520.1560.1630.1710.1770.1810.1470.1520.1610.1720.1810.187)\begin{pmatrix} 0.187 & 0.181 & 0.172 & 0.161 & 0.152 & 0.147\\ 0.181 & 0.177 & 0.171 & 0.163 & 0.156 & 0.152\\ 0.172 & 0.171 & 0.168 & 0.165 & 0.163 & 0.161\\ 0.161 & 0.163 & 0.165 & 0.168 & 0.171 & 0.172\\ 0.152 & 0.156 & 0.163 & 0.171 & 0.177 & 0.181\\ 0.147 & 0.152 & 0.161 & 0.172 & 0.181 & 0.187 \end{pmatrix}

偏りがだいぶなくなりましたが,横棒を 5050 本も引いたのにまだ完全に均等にはなっていません! 66 人であみだくじをやるときに横棒を 5050 本も引く人はあまりいないと思います。

確率の計算

ここからはどうやって確率を計算するのか説明します。行列積の定義を知っている必要があります。

確率の計算

まず n=1n=1,つまり横棒が一本の場合を考えてみる。

丁寧に場合分けして考えてみると,

Pm,1(0,0)=Pm,1(m,m)=11mP_{m,1}(0,0)=P_{m,1}(m,m)=1-\frac{1}{m}

Pm,1(a,a)=12m(a=1,2,,m1)P_{m,1}(a,a)=1-\frac{2}{m}\:(a=1,2,\cdots,m-1)

Pm,1(a,a+1)=1m(a=0,1,,m1)P_{m,1}(a,a+1)=\frac{1}{m}\:(a=0,1,\cdots, m-1)

Pm,1(a,a1)=1m(a=1,2,,m)P_{m,1}(a,a-1)=\frac{1}{m}\:(a=1,2,\cdots, m)

上記以外のとき Pm,1(a,b)=0P_{m,1}(a,b)=0 となる。これを行列で表すと,

(11m1m1m12m1m1m12m1m1m11m)\begin{pmatrix}1-\frac{1}{m}&\frac{1}{m}& & & \\\frac{1}{m}&1-\frac{2}{m}&\frac{1}{m}& &\\&\ddots&\ddots&\ddots&\\&&\frac{1}{m}&1-\frac{2}{m}&\frac{1}{m}\\&&&\frac{1}{m}&1-\frac{1}{m}\end{pmatrix}

となる。この行列を PmP_m とおく。

また,Pm,n(a,b)=c=0mPm,n1(a,c)Pm,1(c,b)P_{m,n}(a,b)=\displaystyle\sum_{c=0}^mP_{m,n-1}(a,c)P_{m,1}(c,b) であることに注意すると,行列 PmP_mnn 乗の各成分が Pm,n(a,b)P_{m,n}(a,b) となることが分かる。

注:推移確率行列が PmP_m であるようなマルコフ連鎖です。→マルコフ連鎖の基本とコルモゴロフ方程式

これで,与えられた m,n,a,bm,n,a,b に対して Pm,n(a,b)P_{m,n}(a,b) を求めることができます。コンピュータに行列のべき乗を計算させればよいのです!

例えば,さきほどの実験結果で紹介した表は P510P_5^{10}P550P_5^{50} を計算したものです。

解析的に確率を求める

行列 PmP_mnn 乗を解析的に計算するのは無理そうだと思っていたのですが,知人が気合いで(固有値,固有ベクトルを求めて)計算してくれました。

Pm,n(a,b)=1m+1+2m+1k=1m{14msin2(kπ2(m+1))}nFk(a)Fk(b)P_{m,n}(a,b)\\=\dfrac{1}{m+1}+\dfrac{2}{m+1}\displaystyle\sum_{k=1}^m\left\{1-\frac{4}{m}\sin^2\left(\frac{k\pi}{2(m+1)}\right)\right\}^nF_k(a)F_k(b)

ただし,Fk(a)=cos{kπm+1(a+12)}F_k(a)=\cos\left\{\dfrac{k\pi}{m+1}\left(a+\dfrac{1}{2}\right)\right\}

三角関数が登場するのが非常に面白いですね。

m2m\geq 2 のとき)mm を固定して nn\to \infty とすると確率は 1m+1\dfrac{1}{m+1} となり,均等に混ざることも分かります。

計算には,固有値など線形代数の知識を使います:

計算の方針

PP の固有値を並べた対角行列を DD,固有ベクトルを並べた行列を XX とすると,PP は対称行列なので P=XDXP=XDX^{\top} となる。つまり,

Pn=XDnXP^n=XD^nX^{\top}

と書ける。よって,PP の固有値と固有ベクトルを求めればよい。

まず,(1,1,,1)(1,1,\dots,1) が固有値 11 に対応する固有ベクトルであることがわかる。また,cos\cos の和積公式を使うと,「cos\cos の等差数列」

(cosθ0,cos(θ0+θ1),cos(θ0+2θ1),,cos(θ0+mθ1))(\cos\theta_0,\cos(\theta_0+\theta_1),\cos(\theta_0+2\theta_1),\dots,\cos(\theta_0+m\theta_1))

が固有ベクトルになりそうだと分かる(難しい)。

境界条件(最初の成分と最後の成分)に注意して θ0,θ1\theta_0,\theta_1 を探すと,残りの固有ベクトルが mm 本見つかる。具体的には,Fk(a)F_k(a)kk 本目の固有ベクトルの第 aa 成分。

最後の計算をはじめ,この記事に全面的な協力をしてくださった知人に感謝。

Tag:難しめの数学雑学・ネタまとめ