あざらしとペンギンの問題

主に漫画、数値計算、幾何計算、TCS、一鰭旅、水族館、鰭脚類のことを書きます。

誕生日のパラドックスと素因数分解・その1

皆さんは「誕生日のパラドックス」と呼ばれるものを耳にしたことがあるでしょうか?それは具体的には次のようなものです。

人が23人集まれば、その中には50%の確率で同じ誕生日を持つ人の組が存在する。

もっと言えば、

人が57人集まれば、その中には99%の確率で同じ誕生日を持つ人の組が存在する。

有名な話なのでご存知の方は少なくないと思われますが、初めてこの話を聞いたとすればさぞ驚くことでしょう。誕生日のバリエーションは1年365日(あるいは366日)もあるのに、それよりもずっと少ない人数が集まっただけで、十分な確率でその中に重複が出てしまうのです。

この意外性をもって、上の言明は「誕生日のパラドックス」と呼ばれています。「パラドックス」とは言っても、論理的な問題があるわけではなくて、あくまで直感に反するということです。実際、上の数値を計算で出すのは、ちょっとした近似を使えばあっさりとできます。そのため、この話は確率計算の入門的話題としてよく引き合いに出されます。

さて、お話だけならこれで終わりなのですが、実はこの誕生日のパラドックス情報科学において様々な応用があります。代表的なのものでは、ハッシュテーブルというデータ構造*1におけるテーブルの大きさを決めるのに、誕生日のパラドックスの結果を利用したりします。

今回は、誕生日のパラドックスを流用して、大きな整数を素因数分解するアルゴリズムを紹介します。このアルゴリズムはポラード(John M. Pollard)のρ(ロー)アルゴリズムと呼ばれています。この \rhoアルゴリズムの「軌跡」に由来していることは後に明らかになるでしょう。

誕生日のパラドックス

まずは出発点の誕生日のパラドックスについて、問題をより一般的に定義しましょう。そのために、1年の日を n 個の整数に、集団の人数を m に置き換えます。そうすると、考えるべき問題は次のように言うことができます。

問題: 整数 1, 2, \dots, n の中から独立に等しくランダムに m 個を選ぶとき、同じ数が重複して選ばれる確率 p(n, m) はどの程度か?あるいは、その確率が 1 - \delta を超えるには m をどの程度にすれば十分か?

ここで、「独立に」とは前の選択の結果に依存せずに次の選択を行うということです。また、「等しくランダムに」とは、1, 2, \dots, n のどの1つも確率 1/n で選ばれるということです。

確率の計算

確率の式を出せばそれを達成する m は出せるので、まずは確率の方から計算しましょう。

確率の計算には、余事象を使います。すなわち、重複が発生しない確率 1 - p(n, m) を先に計算します。

具体的に m = 1, 2, 3, \dots と増やしながら計算してみましょう。

m = 1 の場合は、明らかに重複は発生しないので、

1 - p(n, 1) = 1

です。

m = 2 の場合、1個目を選ぶ時点では上と同じですが、2個目を選ぶときに重複が発生する可能性があります。その確率は、等しくランダムに選ぶことから 1/n です。よって、2個目を選ぶときに重複が発生しない確率は 1 - 1/n となります。さて、2個を選んだ時点で重複が発生しないとは、1個目で重複が発生せず、かつ、2個目でも重複が発生しないことです。1個目と2個目の選択は独立に行われているので、確率の積法則が成り立ちます。よって、求める確率 1 - p(n, 2)1 - p(n, 1)1 - 1/n の積となります。すなわち、

1 - p(n, 2) = 1 - \frac{1}{n}

です。

m = 3 の場合、3個目を選ぶときに重複が発生する確率を考えると、既に2個の異なる数が選ばれているので、2/n となります。よって、3個目で重複が発生しない確率は 1 - 2/n です。ここでもやはり確率の積法則より、求める確率 1 - p(n, 3)1 - p(n, 2)1 - 2/n の積となります。すなわち、

1 - p(n, 3) = \left( 1 - \frac{1}{n} \right) \left( 1 - \frac{2}{n} \right)

ここまで来れば、m = 4 のとき

1 - p(n, 4) = \left( 1 - \frac{1}{n} \right) \left( 1 - \frac{2}{n} \right) \left( 1 - \frac{3}{n} \right)

であることも解るでしょう。

この計算を帰納的に繰り返せば、一般の m については、

1 - p(n, m) = \left( 1 - \frac{1}{n} \right) \left( 1 - \frac{2}{n} \right) \dots \left( 1 - \frac{m - 1}{n} \right)
\displaystyle = \prod_{i = 1}^{m - 1} \left( 1 - \frac{i}{n} \right)

となります。(\displaystyle \prod_{i = 1}^{m - 1}i = 1, 2, \dots, m - 1 で積をとる記号です。)

上と同じ式をもっと直接的に出すこともできます。それには、高校入試などでよくお目にかかる順列を数える式を利用します。

まず、n 個のものの中から重複を許して m 個のものを選んだ列の総数は n^m 通りあります。その中で重複がない列は _n \mathrm{P}_m = n(n-1) \dots (n - m + 1) 通りあります。

1 - p(n, m) は、重複を許した列すべての中から等しくランダムに一つ選んだとき、重複のない列にヒットする確率に等しいので、

\displaystyle 1 - p(n, m) = \frac{_n \mathrm{P}_m}{n^m} = \left( 1 - \frac{1}{n} \right) \left( 1 - \frac{2}{n} \right) \dots \left( 1 - \frac{m - 1}{n} \right)
\displaystyle = \prod_{i = 1}^{m - 1} \left( 1 - \frac{i}{n} \right)

と出ます。これは上で導いた式と同じです。

以上から、p(n, m) の具体的な数値を計算することができるようになりました。めでたしめでたし……ということはないのですね。この式からでは n, m を変数として動かしたときの p(n, m) の大まかな振る舞いが見えてきません。n, m が大きい場合などは数値計算もままならないので、p(n, m) の振る舞いを見るためにはもっと工夫が要ります。

近似計算法

そこで考えるのが近似です。簡単のため m - 1 を偶数とします。まず、上の式の最初と最後、最初から2番目と最後から2番目、最初から3番目と最後から3番目…とペアにして、次のように書き換えます。

\displaystyle 1 - p(n, m) = \prod_{i = 1}^{(m - 1)/2} \left( 1 - \frac{i}{n} \right) \left( 1 - \frac{m - i}{n} \right)
=\displaystyle \prod_{i = 1}^{(m - 1)/2} \left( 1 - \frac{m}{n} + \frac{i(m - i)}{n^2} \right)

ここで、n \gg m として思い切って \frac{i(m - i)}{n^2} の項を無視してしまいましょう。そうすると、

\displaystyle 1 - p(n, m) \approx \prod_{i = 1}^{(m - 1)/2} \left( 1 - \frac{m}{n} \right) = \left( 1 - \frac{m}{n} \right)^{(m - 1)/2}

となります。さらなる近似として、

\left( 1 - \frac{m}{n} \right)^n \approx e^{-m} \; (n \gg m)

を利用すると、

\left( 1 - \frac{m}{n} \right)^{(m - 1)/2} = \left( 1 - \frac{m}{n} \right)^{n (m - 1)/2n}
\approx (e^{-m})^{(m - 1)/2n} = \exp \left( -\frac{m(m - 1)}{2n} \right)

となります。よって、

p(n, m) \approx 1 - \exp \left( -\frac{m(m - 1)}{2n} \right)

という式が得られました。

m - 1 が奇数のときは式が煩雑になるため、上の近似式をそのまま流用することにします。

試しに、n=365 として m を動かして、元の式で計算した正しい確率と近似式で計算した値とを並べると、近似式でもかなり真の確率に近い値が出ることが分かります。表によれば、最も大きくずれているところでも相対誤差は1.5%程度に収まっています。

m真の確率近似式の値相対誤差 m真の確率近似式の値相対誤差
10.00.0- 310.73045463370.72028179760.01392671854
20.0027397260270.0027359764030.001368612854 320.75334752790.74305779680.01365867768
30.0082041658850.0081854929890.002276026097 330.77497185420.76462501190.01335124914
40.016355912470.016303983680.003174924291 340.79531686460.7849718440.01300742014
50.02713557370.027025359430.004061615485 350.81438323890.80409726370.01263038663
60.040462483650.040262904050.004932460439 360.83218210640.82200990870.01222352364
70.05623570310.0559104420.005783889534 370.84873400820.83872713490.01179035272
80.074335292350.07384375620.006612419658 380.86406782110.85427403850.01133450675
90.094623833890.093922229290.007414671033 390.87821966440.86868246880.01085969259
100.11694817770.11599067810.008187383999 400.89123180980.88199004590.01036965221
110.14114137830.13988134770.008927435717 410.90315161150.89423920040.009868122872
120.16702478880.165416030.009631856773 420.91403047160.90547624590.009358797065
130.19441027520.19240826780.01029784764 430.92392285570.91575049650.008845283048
140.2231025120.2206656090.01092279495 440.93288536860.92511343840.008331066659
150.25290131980.24999187030.01150428749 450.94097589950.93361796170.007819475224
160.28360400530.28018937560.01204013184 460.94825284340.94131765940.007313644301
170.31500766530.31106113350.0125283676 470.95477440280.94826619460.006816487936
180.34691141790.34241291970.01296728194 480.96059797290.95451674120.006330673044
190.3791185260.37405523760.01335542351 490.96577960930.96012149440.0058585984
200.41143838360.40580512750.01369161538 500.97037357960.96513125410.005402378632
210.44368833520.43748780530.01397496695 510.97443199330.96959507520.004963833423
220.47569530770.46893811080.01420488446 520.97800450930.97355998540.004544482033
230.50729723430.50000175220.01438108006 530.98113811350.97707076320.004145543066
240.53834425790.53053633940.01450357906 540.98387696280.98016977410.003767939274
250.5686997040.56041219950.01457272512 550.98626228880.9828968590.003412307062
260.59824082010.58951297520.0145891832 560.98833235490.98528926950.003079010233
270.62685928230.61773600990.0145539399 570.99012245930.98738164450.002768157424
280.65446147230.64499252680.01446830101 580.99166497940.98920602450.002479622609
290.68096853750.67120761210.01433388604 590.99298944840.99079189530.002213067984
300.70631624270.69632001770.01415261945 600.99412266090.99216625870.001967968572

それでは、n が大きいときの p(n, m) の振る舞いについて見てみましょう。まず、m = 0 の付近では小さい値で推移します。その後増加のスピードが速くなり、m \approx 1.177 \sqrt{n} で50%を通過します。そして、m \gg \sqrt{n} のときはほとんど1に近い値となります。

n=365 に対して 1 - \exp \left( -\frac{m(m - 1)}{2n} \right)m を横軸にとってプロットした図を示します。
f:id:azapen6:20130501073732p:plain

次に、p(n, m) = 1 - \delta となる m の値を見積もってみましょう。これは、少し計算すれば

m = \frac{1 + \sqrt{1 + 8n \ln (1/\delta)}}{2} \approx \sqrt{2n \ln (1/\delta)} + 0.5

と出ます。(もう一つの解は負になるので除きました。)例えば、誕生日が重複する確率が90%以上になるためには、\delta = 0.1 として

m = \sqrt{2 * 365 * \ln 10} + 0.5 \approx 41.5

から42人程度がいれば十分ということになります。(正しい値は41人です。)

今回はここまで

ちょっと話が長くなりそうなので、ここで切ります。次回はポラードのρアルゴリズムによる素因数分解について説明します。まぁ、言ってしまえばアホみたいに簡単なアルゴリズムなのですが、どうして上手くいくかがよく解っていないという意外に手強い面もあります。

それでは次回に。

*1:データを目的に応じて効率的に出し入れするための枠組みのこと。