数学のお勉強

数学 Ruby

某プログラム勉強会で、お題が出たので解いてみる。

[問1]6面体のサイコロを3回振る。同じ目が一度も出ない確率を計算せよ。

6P3/63 = 5/9

簡単な確率問題。
3回の試行で同じ目が出ない組み合わせ(6P3)を、3回の試行の全ての可能性(63)で割る。

[問2]128bitの乱数を1兆個発生させる。同じ値が一度も出ない確率のおおよその値を求めよ。（必要なら、数学公式集などを利用する）

[問1]と同じやり方で解こうとすると死ねる。

乱数のとりえる値の数をn、試行回数をkとすると
nPk/nk
ここでは、n=2128,k=1012
まともに計算したら、どのくらいかかるのだろう？

で、回答
0.9999999999999985306320614736099829776071878571538233981958611321639300600254765538643653
方法は後述。

[問3]128bitの乱数を一日に10万回発生させる。同じ値が一度も出ない確率が99.9%を下回るのは約何日後か。

128bitの乱数のうち、既に発生した数の割合が0.1%を超える日をn日目とする。
100000 * n / 0.001 > 2128
n > 2128 / 108

n > 3402823669209384634633746074317

その前に人類が死滅しているのは確実ですな。

[問4]「N bit の乱数を何個発生させたら同じ値が一度も出ない確率が p % を下回るか」のおおよその値を計算するプログラムを書け。また、Nは32〜512 辺り、p は50%〜99.9999% 辺りを想定している。

わーい、やっとプログラムだ。まずk回の試行で同じ値が一度も出ない確率（n=2N）をp(k)とする。
p(k) = nPk/nk
ここで、分子と分母を別に計算するとDoubleの精度ではあっというまに指数部がオーバーフローしてしまう。
BigDecimal*1を使えばいいけれど遅い。

ここで
p(2) = (n - 1) / n
p(k) = p(k-1) * (n-(k-1)) / n
なので、kを増やしながら (n-(k-1)) / n をかけていって、指定された確率以下になるkを求めてみる。

def try_count1(bits, rate)
  n = 2 ** bits
  r = 1.0
  k = 2
  while(k < n)
    r *= (n-(k-1)).to_f/n
    break if r < rate
    k += 1
  end
  k
end

予想されたことではあるけれど、使い物になるのはbit数が40くらいまで。

ここで、k回目の試行で重複する確率を考えてみる。

2回目の試行で重複する確率は(1/n)
3回目の試行で重複する確率は(2/n)
…
k回目の試行で重複する確率は((k-1)/n)

k回目までで重複する確率は
(1/n) + (2/n) + .. (k-1)/n = k*(k-1)/(n*2)（等差数列の和）

はい、おかしいですね。
k回目の試行で重複する確率が((k-1)/n)になるのは、それまでの試行で重複が無かった場合だけなので厳密にはこの計算じゃ駄目なのだけれど、kに対してnが十分大きければ、この値で近似してもいいじゃん（ほんとか？）。

require "bigdecimal"

def dup_rate(n, k)
  r = BigDecimal.new(k.to_s)
  r *= (k - 1)
  r /= (n*2)
end

def try_count2(bits, rate)
  rate /= 100
  rate = 1.0 - rate
  n = 2 ** bits
  min_try = 1
  max_try = n
  while(max_try > min_try + 1)
    try_cnt = (min_try + max_try) / 2
    if dup_rate(n, try_cnt) < rate
      min_try = try_cnt
    else
      max_try = try_cnt
    end
  end
  max_try
end

最初のプログラムと実行結果を比べてみる
N=40 p = 99.999%
try_count1: 4690
try_count2: 4690

N=40 p = 99%
try_count1: 148665
try_count2: 148292

N=40 p = 70%
try_count1: 885629
try_count2: 812224

N=40 p = 50%
try_count1: 1234605
try_count2: 1048577

まあ、使えるかな。
ということで、[問2]も

puts (1 - dup_rate(2**128,10**12))

で解いたのでした。

*1:scala(java)のBigDecimalは桁数が増えるとどんどん遅くなるけれど、RubyのBigDecimalは仮数部の精度を指定できるのでそんなに遅くならない。仮数部の精度はそんなに要らない用途には向いているかも。