Blog

2011.01.15

乱択アルゴリズム紹介(行列乗算の検査&多項式等価性の検査)

Yuichi Yoshida

Engineer

吉田です。今回は乱数を用いたアルゴリズム(Randomized Algorithms、乱択アルゴリズム)を紹介したいと思います。
理論の世界では乱数を使ったアルゴリズムは既に当たり前のものになっているのですが、実際の応用で使われている所は残念ながら余り見たことが無いです。多分それは宣伝が足りないのだろうと思ったので、今回少し書いてみることにしました。実は他の場所で話すことになっていることの下準備も兼ねているのですが。これから書くことがそのまま実用に耐えるとは思っていませんが、それで乱択アルゴリズムに関する感覚を蓄えれば他の形で応用出来るんじゃないかと考えています。
最初に参考文献を挙げておくと、
Randomized Algorithms by Rajeev Motwani and Prabhakar Raghavan
Probability and Computing: Randomized Algorithms and Probabilistic Analysis by Michael Mitzenmacher and Eli Upfal
がお薦めです。
さてアルゴリズムに乱数を導入することのメリットの一つは、多くの問題で計算時間を小さくすることが出来ることです。例えばグラフの最小カットは決定性アルゴリズム(乱数を使わないアルゴリズム)では、僕の記憶では\(O(nm \log n^2/m)\)が最速のはずですが、乱択アルゴリズムで良ければ\(O(m \log^3 n)\)で実現出来ます。ここで\(n\)は頂点数、\(m\)は枝数です。
もう一つのメリットは、乱数を用いるとしばしばアルゴリズムが非常に簡単になることです。例えば一般のグラフが完全マッチングを持つかは多項式時間で判定出来るのですが、その実装は非常に重たいです。それに対して乱数を使って良いのであれば比較的簡単に実装することが出来ます。
乱択アルゴリズムは、当たり前ですが、乱数によってその挙動が変わります。その変わり方によって大きく二種類に分類されています。
ラスベガスアルゴリズム: 常に正しい解を出すが、乱数によって計算時間が変わる。
モンテカルロアルゴリズム: 常に同じ計算時間だが、乱数によっては間違った解を出すかもしれない。
ある解が与えられた時にそれが正しいかチェックするアルゴリズムが有れば、モンテカルロアルゴリズムをラスベガスアルゴリズムに変換するのは簡単ですが、逆の変換は一般には簡単ではないです。今回見る例は全てモンテカルロアルゴリズムになっています。
さて今回最初に考えるのは行列乗算です。
(行列乗算)
入力: 行列\(A, B \in \mathbb{F}^{n\times n}\)
出力: 行列\(C = AB\)
つまり二つの行列を受け取って、その積を出力しなさいという問題です。さて\(C\)を求める計算時間は何になるでしょうか?
単純には、
\(C_{ij} = \sum_{k}A_{ik}B_{kj}\)
なので、各\(C_{ij}\)を求めるのに\(O(n)\)時間、\(C\)全体を求めるのは\(O(n^3)\)時間になります。実はそれよりかなり速い、Coppersmith-Winogradによる\(O(n^{2.376})\)時間アルゴリズムも知られています。自明な下限は出力の数から\(\Omega(n^2)\)になり、多くの人が\(O(n^2\mathrm{polylog} n)\)時間で行列乗算は出来ると考えているのですが、まだ達成できていません。
少し問題設定を変えて次のようにしてみましょう。
(行列乗算の検査)
入力: 行列\(A, B, C \in \mathbb{F}^{n \times n}\)
出力: \(AB = C\)ならYes、そうでなければNo。
つまり行列\(A,B\)に対して行列\(C\)が\(A\)と\(B\)の積になっているかを検査する問題です。乱数を用いないと、結局\(AB\)を計算して、それを\(C\)と比べる以外の方法は知られていません。つまり\(O(n^{2.376})\)時間が最速ということになります。それに対して、乱数を用いるとこれを\(O(n^{2})\)時間に出来ることが知られていますので、それを紹介したいと思います。
アルゴリズムは非常に単純です。まず\(r \in \{0,1\}^{n}\)をランダムに選び、次に\(x = ABr\)と\(y = Cr\)を計算します。最後に\(x = y\)ならYes、そうでなければNoと出力します。
最初に計算時間を考えましょう。注意してほしいのは行列とベクトルの乗算は\(O(n^2)\)時間で計算出来ることです。なので\(x = A(Br)\)という順に計算すれば、\(x\)は\(O(n^2)\)時間で計算でき、\(y\)についても同様です。結局\(O(n^2)\)時間となり、決定性アルゴリズムの最速\(O(n^{2.376})\)よりも随分と速くすることが出来ました。
次にアルゴリズムの正しさを見ましょう。
(定理) 上記のアルゴリズムは、\(AB=C\)の時は確率\(1\)でYesと出力し、\(AB\neq C\)の時は確率\(1/2\)でNoと出力する。
(証明) \(AB=C\)の時: 任意の\(r\)に対して\(x = y\)になるので、確率\(1\)でYesを出力する。
\(AB \neq C\)の時: \(D = AB-C\)と置くと、ある\(s,t\)が存在して\(D_{st} \neq 0\)である。\(x – y = Dr\)だが、\(Dr\)の\(s\)行目が\(0\)になる必要十分条件は
\(r_t = \frac{\sum_{k \neq t}D_{sk}r_k}{D_{st}}\)
である。これが成り立つ確率は高々\(1/2\)であるので、少なくとも確率\(1/2\)で\(x – y \neq 0\)となり、Noと出力する。
\(AB \neq C\)の時にNoと出力する確率が\(1/2\)なのが不満ですが、これは上のアルゴリズムを繰り返し実行することで、小さくすることが出来ます。
(定理) 行列乗算の検査に対して、\(AB=C\)の時は確率\(1\)でYesと出力し、\(AB\neq C\)の時は確率\((1/2)^k\)でNoと出力する\(O(kn^2)\)時間アルゴリズムが存在する。
\(k = 100\)とでもしてあげれば、\((1/2)^k \approx 10^{-30}\)となりもう十分でしょう。
上で紹介したアルゴリズムはFreivalsのアルゴリズムと呼ばれている非常に古典的なものです。Freivalsの手法はとても単純なので、他にも色々な問題に応用されてきました。その一例として多項式等価性について見たいと思います。最初にいくつか言葉を定義しておきます。
\(P(x_1,\ldots,x_n) \in \mathbb{F}[x_1,\ldots,x_n]\)を多変数多項式とします。\(P\)の項の次数は項中の変数の次数の総和で定義されます。例えば\(x_1^2x_2x_4^3\)という項であれば、次数は\(6\)です。\(P\)の全次数は項の次数の最大値のことを表します。例えば\(P(x_1,x_2,x_3,x_4)=x_1+x_1^2x_3+x_1^2x_2x_4^3\)であれば、全次数は\(6\)です。
(多項式等価性)
入力: 多変数多項式 \(P_1(x_1,\ldots,x_n),P_2(x_1,\ldots,x_n) \in \mathbb{F}[x_1,\ldots,x_n]\), \(P_1,P_2\)の次数は高々\(d\).
出力: \(P_1(x_1,\ldots,x_n) = P_2(x_1,\ldots,x_n)?\)
\(P_1\)と\(P_2\)が明示的に与えられているのであれば、単に各係数を比較するだけで終わります。なので、具体的な\(x_1,\ldots,x_n\)に対して\(P_1(x_1,\ldots,x_n)\)と\(P_2(x_1,\ldots,x_n)\)を計算するのは簡単だが、その係数を求めるのは大変という状況を考えてみましょう。こういう状況は実際にあって、例えば\(x_1,\ldots,x_n\)を使って記述された行列の行列式は、具体的に\(x_1,\ldots,x_n\)の値が決まれば多項式時間で計算出来ますが、\(x_1,\ldots,x_n\)を記号として置いたままだと\(O(n!)\)時間かかる方法しか知られていません。このような行列が現れる場面としては、グラフに完全マッチングが有るかの判定などが挙げられるのですが、それはまた別の機会に話をしたいと思います。
さて先ほどのFreivalsの手法を真似ると次のようなアルゴリズムが思いつきます。\(Q(x_1,\ldots,x_n) = P_1(x_1,\ldots,x_n) – P_2(x_1,\ldots,x_n)\)としましょう。まず有限集合\(\mathbb{S} \subseteq \mathbb{F}\)を適当に固定して、\(r_1,\ldots,r_n \in \mathbb{S}\)をランダムに選びます。\(\mathbb{S}\)を固定する理由は、\(\mathbb{F}\)が無限集合の時(例えば実数全体)、\(\mathbb{F}\)からランダムに選ぶには無限ビットの乱数が必要になってしまうからです。こういう所だけは妙に厳密になるのが理論の世界なのでしょう。その後は、もし\(Q_1(x_1,\ldots,x_n) = 0\)であればYesと出力し、そうでなければNoと出力します。
先ほどの行列乗算の検査では\(\mathbb{S} = \{0,1\}\)として選んでいましたが、今回はどのようにすれば良いでしょうか。\(n = 1\)であれば、\(\mathbb{S} \geq d+1\)と取れば良いことが分かります。実際\(Q\)の次数は高々\(d\)なので、恒等的に零でない\(Q\)は高々\(d\)個しか解を持ちません。よって\(\Pr[Q(r_1) = 0 \mid Q(x_1) \not \equiv 0] \leq d / \mathbb{S} < 1\)となります。確率が\(1\)未満であれば、行列乗算の検査の時と同じように、このアルゴリズムを繰り返し適用することで任意の精度が得られます。
さて上のアルゴリズムがが任意の\(n\)でも成り立つと示したのが以下の定理です。
(Schwartz-Zippelの定理) \(Q(x_1,\ldots,x_n) \in \mathbb{F}[x_1,\ldots,x_n]\)を全次数\(d\)の多変数多項式とする。\(\mathbb{S} \subseteq \mathbb{F}\)とし、\(r_1,\ldots,r_n \in \mathbb{S}\)をランダムに選ぶ。この時
\(\Pr[Q(r_1,\ldots,r_n) = 0 \mid Q(x_1,\ldots,x_n) \not \equiv 0] \leq \frac{d}{|\mathbb{S}|}\)
が成り立つ。
(証明) \(n\)に関する帰納法を用いる。\(n = 1\)の時は既に示したので、任意の\(n’ < n\)で定理が成り立つと仮定する。
\(Q(x_1,\ldots,x_n)\)から\(x_1\)を外に出すと次のように書ける。
\(Q(x_1,\ldots,x_n)=\sum_{i=0}^k x_1^iQ_i(x_2,\ldots,x_n)\)
ここで\(k \leq d\)は\(x_1\)の最大次数である。また\(k = 0\)であれば帰納法の仮定から定理が成り立つので、\(k > 0\)としてよい。\(k\)の選び方から\(x_1^k, Q_k(x_2,\ldots,x_k) \not \equiv 0\)であることに注意。\(Q_k\)の全次数は高々\(d – k\)であるので、帰納法の仮定から\(Q_k(r_2,\ldots,r_n) = 0\)となる確率は高々\((d-k)/|\mathbb{S}|\)である。
\(Q_k(r_2,\ldots,r_n) \neq 0\) であったとして、次の一変数多項式を考える。
\(q(x_1) = Q(x_1,r_2,\ldots,r_n) = \sum_{i=0}^kx_1^iQ_i(r_2,\ldots,r_n)\)
多項式\(q(x_1)\)は次数\(k\)であり、\(x_1^k\)の係数が\(Q_k(r_2,\ldots,r_n)\)であることから、恒等的に零ではない。
よって帰納法の仮定から\(q(r_1)=Q(r_1,\ldots,r_n)\)が\(0\)になる確率は高々\(k / |\mathbb{S}|\)である。
ここまでで以下の二つを示した。
\(\Pr[Q_k(r_2,\ldots,r_n) = 0] \leq \frac{d-k}{|\mathbb{S}|}\)
\(\Pr[Q(r_1,\ldots,r_n) = 0 \mid Q_k(r_2,\ldots,r_n) \neq 0] \leq \frac{k}{|\mathbb{S}|}\)
ここで確率事象\(X_1,X_2\)に対して、
\(\Pr[X_1] = \Pr[X_1 \mid \overline{X_2}]\Pr[\overline{X_2}] + \Pr[X_1 \mid X_2]\Pr[X_2] \leq \Pr[X_1 \mid \overline{X_2}] + \Pr[X_2]\)
が成り立つことに注意すると、
\(\Pr[Q(r_1,\ldots,r_n) = 0] \leq \Pr[Q(r_1,\ldots,r_n) = 0 | Q_k(r_2,\ldots,r_n) \neq 0] + \Pr[Q_k(r_2,\ldots,r_n)] = \frac{d}{|\mathbb{S}|}\)
  • Twitter
  • Facebook