次の式は 3 つの部分に分けて掛け合わせたものになります。
まず、
は Python で書くとこうですね。
次に、
は 順列組み合わせ の記号で、Python では次の関数で表せます。
python
1def comb(n, p):
2 return math.factorial(n) // (math.factorial(p) * math.factorial(n - p))
最後に、
は シグマ といって、例えば
は h(0) + h(1) + ... + h(p) なので、Python ではこんなふうに書けます。
python
1sum = 0
2for j in range(p + 1):
3 sum += h(j)
また、シグマの中身
はこうですね。
python
1pow(-1, p - j) * comb(p, j) * pow(j, i)
以上をまとめると、最初の式は次の関数 f で表すことができます。
(なお、引数 n はカードの種類、i は購入した枚数、p は何種類揃うかを意味します。)
python
1import math
2
3def comb(n, p):
4 return math.factorial(n) // (math.factorial(p) * math.factorial(n - p))
5
6def f(n, i, p):
7 sum = 0
8 for j in range(p + 1):
9 sum += pow(-1, p - j) * comb(p, j) * pow(j, i)
10 return comb(n, p) * sum / pow(n, i)
これを使って、5 種類のカードを i 枚購入して全種類をコンプリートする確率のグラフを描画してみます。
python
1import matplotlib.pyplot as plt
2
3x = range(0, 31)
4y = [f(5, i, 5) for i in x]
5plt.plot(x, y)
また、社員が i 人いる場合に年間 100 日以上パーティーが開かれる確率ですが、1 年は n = 365 日と考えると、これがカードの種類に相当し、社員 i 人というのはカードを i 枚購入した場合に相当するので、ちょうど p = 100 日開かれる可能性は f(365, i, 100) です。
で、100 日以上になる確率は、ちょうど 100 日の確率と、ちょうど 101 日の確率と、(中略) ちょうど i 日の確率をすべて足し合わせたものになるので、次のコードで描画できます。
python
1def g(n, i, p):
2 sum = 0
3 for j in range(p, min(i, n) + 1):
4 sum += f(n, i, j)
5 return sum
6
7x = range(100, 131)
8y = [g(365, i, 100) for i in x]
9plt.plot(x, y)