下記の記事を読んで、子供の頃に「天皇が交代するたびに祝日が増えれば毎日が祝日になるなぁ(不謹慎」と思っていたのを思い出しました
せっかくなので回数の期待値(平均)を計算してみた
2月29日(うるう日)を考えない場合
\(E(x)\)を残り\(x\)個の日付が埋まるまでの回数の期待値とします。
このときもちろん\(E(0)=0\)です
残り\(x\)個のとき、残りの日付が\(1\)個減る確率は\(\frac{x}{365}\)となります
この場合期待値は確率の逆数として考えてよいので、\(E(x) = E(x - 1) + \frac{365}{x}=\sum_{k=1}^{x}\frac{365}{k}\)となる
計算すると\(E(365)\approx 2364\)回となる。
以下の表によると日本の男性の平均寿命は\(79\)年なので(この先も男性とは限らないですが)、あと\(19\)万年弱ぐらい経過すれば毎日が天皇誕生日になることがわかる
追記
年数について重大な指摘がされていました
しかしよく見ると前提がおかしいな。79歳で死んで0歳のひ孫が即位ということはなさそう。生殖で考えると平均30歳ぐらいとして、8万年ぐらいか。
\(E(x)\)のグラフ(横軸が\(x\))
\(x\)が大きい時にはあまり回数を必要としない(縦の減少幅が小さい)が、\(x\)が小さくなると多くの回数が必要になる
ちなみにこういうタイプの問題はクーポンコレクター問題という名前で知られています。
ソースコード
from pylab import * s = 0 result = [s] for k in xrange(1, 365): s += 365.0 / k result.append(s) xlabel('#day') xticks(xrange(1, 366, 50), xrange(365, 0, -50)) plot(xrange(1, 366), result[::-1]) show()
2月29日(うるう日)を考える場合
うるう日を加えると計算が複雑になります(少し合っているか不安
\(E(x, y)\)を残り\(x\)個の日付(うるう日以外)と残り\(y\)個の日付(うるう日)が埋まるまでの回数の期待値とします。
このとき\(E(0, 0)=0\)です
まずうるう日以外とうるう日になる確率を考えます。
うるう日は\(400\)年に\(97\)回出現するので、うるう日になる確率は\(p_l=\frac{97}{365\times400+97}\)
うるう日以外の日付が残り\(x\)個あるとき、\(x\)個のうちのいずれかが選ばれる確率は\(p_n(x)=\frac{400x}{365\times400+97}\)
うるう日を考えない場合、\(x>0\)について次の式が成り立つ
$$E(x, 0) = p_n(x)E(x - 1, 0)+(1-p_n(x))E(x, 0)+1$$
これを式変形すると上でやった計算と同様の式になる
$$E(x, 0) = E(x - 1, 0) + \frac{1}{p_n(x)}$$
次にうるう日を考えた場合、まず\(x=0\)について
$$E(0, 1) = p_lE(0, 0)+(1-p_l)E(0, 1)+1$$
これを式変形すると\(E(0, 1) = \frac{1}{p_l}\)
次に\(x>0\)について
$$E(x, 1) = p_lE(x, 0)+p_n(x)E(x - 1, 1)+(1-p_l-p_n(x))E(x, 1)+1$$
式変形すると
$$E(x, 1)=\frac{p_lE(x, 0) + p_n(x)E(x - 1, 1)+ 1}{p_l+p_n(x)}$$
以上より動的計画法を使えば簡単に計算できて、期待値\(E(365, 1) \approx 2693\)回となる
ソースコード
pl = 97.0 / (365 * 400 + 97) pn = 400.0 / (365 * 400 + 97) noleap = 0 leap = 1 / pl for x in xrange(1, 366): noleap += 1 / pn / x leap = (pl * noleap + pn * x * leap + 1) / (pl + pn * x) print leap