2021年11月20日土曜日

確率変数の中央値と期待値の差

確率集中不等式はしばし実数値確率変数$Z$に対して$\Pr[|Z-\mathbf{E}Z|\geq t]$の上界を与えますが, Talagrandの不等式に代表されるようなisoperimetric inequalityを考えると期待値$\mathbf{E}Z$の代わりに中央値$\mathbf{M}Z$を考えることがよくあります. ここで確率変数$Z$の中央値(median)とは, $\Pr[Z\geq \mathbf{M}Z]\geq 1/2$ かつ $\Pr[Z\leq\mathbf{M}Z]\leq 1/2$を満たす実数$\mathbf{M}Z$のことです. Isoperimetric inequalityを使うと$\Pr[|Z-\mathbf{M}Z|\geq t]$の上界を与えることができます.

平均値と中央値は一般に異なりますが, 例えば$Z$の分散が小さいとバラつきが小さいので $\mathbf{E}Z$と$\mathbf{M}Z$の間はそれほどギャップがないように思われます. この直感は実際正しく, 以下の事実が成り立ちます.

定理

$\mathbf{E}Z^2<\infty$を満たす確率変数$Z$に対し,

$$|\mathbf{M}Z-\mathbf{E}Z|\leq \sqrt{\mathrm{Var}Z}$$


本記事ではこの事実を証明していきます.


補題 (Cantelli's inequality)

任意の$t>0$に対し,

$$\Pr\left[\frac{Z-\mathbf{E}Z}{\sqrt{\mathrm{Var}Z}} \geq t\right] \leq \frac{1}{1+t^2}.$$


(補題の証明)

平行移動とスケーリングによって一般性を失わず$\mathbf{E}Z=0$, $\mathrm{Var}Z=1$と仮定してよい. すなわち, 平均0で分散1の確率変数$Z$に対して$\Pr[Z\geq t]\leq 1/(1+t^2)$を示せば良い.

マルコフの不等式より, 任意の$\lambda>0$に対して

$$\Pr[Z\geq t] \leq \Pr[(Z+\lambda)^2\geq (t+\lambda)^2] \leq \frac{\mathbf{E}Z^2+\lambda^2}{(t+\lambda)^2}=\frac{1+\lambda^2}{(t+\lambda)^2}.$$

関数$\lambda\mapsto \frac{1+\lambda^2}{(t+\lambda)^2}$は$\lambda=1/t$の時に最小化され, このとき$1/(1+t^2)$をとる (証明終).


(定理の証明)

補題において$t=1$を代入すると, $\Pr[Z\geq \mathbf{E}Z+\sqrt{\mathrm{Var}Z}] \leq 1/2$を得る. 中央値の定義より

$$\Pr[Z\geq \mathbf{E}Z+\sqrt{\mathrm{Var}Z}] \leq 1/2 \leq \Pr[Z\geq \mathrm{M}Z].$$

よって$\mathrm{M}Z\geq\mathbf{E}Z+\sqrt{\mathrm{Var}Z}$を得る.

一方, 新たな確率変数$-Z$に対して$t=1$として改めて補題を適用すると, $-1$倍しても分散は変わらないため, $\Pr[Z\leq \mathbf{E}Z-\sqrt{\mathrm{Var}Z}]\leq 1/2$となり, 同様にして$\mathrm{M}Z\leq\mathbf{E}Z-\sqrt{\mathrm{Var}Z}$を得る (証明終).



2021年11月19日金曜日

モーメント母関数のlogについて

実数値をとる確率変数 $X$ と適当な値$t\in\mathbb{R}$に対して, 確率$\Pr[X\geq t]$を抑えるときにモーメント母関数 (MGF: Moment Generating Function) を考えたくなります. 確率変数$X$のモーメント母関数は$M_X(c)=\mathrm{E}[e^{cX}]$で定義されます. これを使うと, 任意の$\lambda\geq 0$に対して

$\Pr[X\geq t] = \Pr[\lambda X \geq \lambda t] = \Pr[e^{\lambda X} \geq e^{\lambda t}] \leq e^{-\lambda t} M_X(\lambda)$

が成り立つ (最後の不等式ではMarkovの不等式を適用した) わけです. 一般に浸透された言い回しではありませんが, モーメント母関数のlogをとったものをここではlogMGFと呼ぶこととします. すなわち, 確率変数$X$のlogMGFは

$$\Psi_X(\lambda) = \log M_X(\lambda) = \log \mathrm{E}[e^{\lambda X}]$$

です. これを使うと任意の$\lambda\geq 0$に対して

$$\Pr[X\geq t] \leq \exp(-\lambda t + \Psi_X(\lambda))$$

となるので, 最適な$\lambda$をチューニングしたくなります. ここで$\Psi_X^*(t) = \sup_{\lambda\geq 0}\{\lambda t - \Psi_X(\lambda)\}$ を考えてやると,

$$\Pr[X\geq t] \leq \exp(-\Psi^*_X(t))$$

を得ますが, $\Psi_X^*(t)$をよくみると関数$\Psi_X(\lambda)$のルジャンドル変換っぽいことが分かります ($\sup$で動く$\lambda$の範囲が$\lambda\geq 0$でなく$\lambda\in\mathbb{R}$ならばルジャンドル変換と一致します. なお, $t\geq \mathrm{E}[X]$ならば, $\Psi_X^*(t)=\sum_{\lambda\in\mathbb{R}}\{\lambda t - \Psi_X(\lambda)\}$ となり, ルジャンドル変換と一致します).

$\Psi_X(\lambda)$の有名な性質として凸性が挙げられます. つまり, $S=\{\lambda>0\colon \Psi_X(\lambda)<\infty\}$とすると$\Psi_X(\lambda)$は$S$上で凸関数になっています. これは, $x,y\in S$を$0<x<y$を満たす実数として$p\in(0,1)$を任意に選んで$\Psi_X(px+(1-p)y)$をHölderの不等式を使って上から抑えることによって証明することができます.

2021年10月11日月曜日

直径の1.5近似アルゴリズム

概要
$n$頂点$m$辺グラフの直径の1.5近似を$\tilde{O}(m\sqrt{n}+n^2)$時間で計算するAingworthら[1]のアルゴリズムを紹介します ($\tilde{O}(\cdot)$は$\mathrm{polylog}(n)$ factorを無視). 次にこのアルゴリズムの計算時間を改良したRoditty and Williams[2]による$\tilde{O}(m\sqrt{n})$時間乱択アルゴリズムを紹介します.


1. 導入

重みなし無向グラフ $G=(V,E)$ に対し二頂点$u,v$間の距離 $\mathrm{dist}(u,v)$ を最短経路長とします. このとき, $G$の直径 $\mathrm{diam}(G)$ は $\mathrm{diam}(G):=\max_{u,v\in V}\mathrm{dist}(u,v)$ で定義されます.

本記事のトピックは$n$頂点$m$辺グラフ$G$が与えられた時にその直径$\mathrm{diam}(G)$を計算する問題の計算量です. なお, 計算モデルは$O(\log n)$-word RAMを仮定します.

各頂点からBFSをすることによって$O(mn)$時間で計算することができます. 従ってこれよりも高速に計算できるかが焦点になります.

グラフの直径の計算は精緻な計算量(fine-grained complexity)においてかなり活発に議論されているトピックで, 直径3以下のグラフが与えられた時にその直径が$2$であるか$3$であるかを$O(m^{2-\epsilon})$時間で判別することはSETHと呼ばれる仮定の下では不可能であることが証明されています ($\epsilon>0$は任意の定数)[2]. すなわち, SETHを仮定すると直径の1.49近似を$O(m^{2-\epsilon})$時間で計算できないので, 例えば疎なグラフ ($m=O(n)$) では各点BFSによる$O(n^2)$時間アルゴリズムが最適であることを意味します. ここで, 直径の$c$近似 ($c>1$) というのは直径の$1/c$倍以上の値を出力することを意味します (例えば最大カットなどの文脈では$c^{-1}$近似の言い回しをよく使いますが, 直径近似研究の文脈では$c$近似と言うことが多く, 本記事もそれにならいます). ですので1.49近似アルゴリズムがあれば, 直径3のグラフに対してそのアルゴリズムを走らせると出力値は$\geq 3/1.49>2$となるため, 直径2と直径3のグラフを区別することができます.

それでは, SETHの下での近似困難性はタイトなのでしょうか? 具体的には, 疎なグラフに対し$O(n^{1.99})$時間で直径の1.5近似が計算できるかどうかが気になります. 

まず, 与えられるグラフ$G$が木である場合を考えてみます. このときはよく知られたdouble sweepと呼ばれる方法を使えば$O(m)$時間で直径の厳密値を計算できます. また, 一般グラフに対して適用すると, 直径の2近似が$O(m)$時間で計算できます.

本記事ではまず, Aingworthら[1]による$\tilde{O}(m\sqrt{n}+n^2)$時間で直径の"ほぼ"1.5近似を計算するアルゴリズムを紹介します. 具体的には彼らのアルゴリズムは $\lfloor 2/3 \mathrm{diam}(G) \rfloor \leq D \leq \mathrm{diam}(G)$を満たす値$D$を出力します. なお, 乱択を許しちょっと修正することによってAingworthらのアルゴリズムと同じ条件を満たす値$D$を$\tilde{O}(m\sqrt{n})$時間で出力する乱択アルゴリズムが知られています[2].


2. Aingworthらのアルゴリズム


Aingworthらによるアルゴリズム[1]を紹介します. double sweepを知っておくと感覚を若干ながら掴みやすいかもしれませんが, 知らなくても大丈夫です. また, 以下の説明ではBFSが出てきますが, これをDijkstra法に置き換えるとweighted graphに対して良い感じで動くと思うのですが, 私はあまり深く考えたことがないのでこの記事では割愛します. 以下, $O(m+n)$と書くのが面倒なので, $m\geq n$を仮定し, 単一始点BFSの時間計算量は$O(m)$と書いていきます.


アルゴリズム

1. 各頂点$v\in V$に対しBFSを行い$v$に近い順に頂点を$\sqrt{n}$個選び$S_v$とする (BFSを行い$\sqrt{n}$個の頂点に訪問したら途中で打ち切るようにする). なお, $v\in S_v$とする.

2. ステップ1で得られた集合族$(S_v)_{v\in V}$に対し, サイズ$O(\sqrt{n}\log n)$の hitting set $H\subseteq V$ を計算する (つまり, $H$は全ての$v\in V$に対し$H\cap S_v\neq\emptyset$を満たし, かつ$|H|=O(\sqrt{n}\log n)$を満たす集合. 存在性は後で示します).

3. $H$の各頂点$h\in H$に対し, $h$を始点としたBFSを行い, 全ての$h\in H,v\in V$に対し$\mathrm{dist}(h,v)$を計算する. その後, $D_1:= \max_{h\in H,v\in V}\mathrm{dist}(h,v)$とする.

4. $H$から最も遠い頂点を$u$とする. つまり$u$は$\min_{h\in H}\mathrm{dist}(h,v)$を最大にする$v$である.

5. ステップ4で得られた頂点$u$に対し, 各$x\in S_u$を始点としてそれぞれBFSを行い, 得られた距離の中での最大値を$D_2$とする. すなわち, $D_2:= \max_{x\in S_u,v\in V}\mathrm{dist}(x,v)$. 最後に$\max\{D_1,D_2\}$を出力.




計算時間:

- ステップ1では各頂点に対してBFSを行うので, 各$v\in V$に対し$O(|E(G[S_v])|)=O(n)$時間かかるので計$O(n^2)$時間です.

- ステップ2を考えます. まず存在性を議論します. 各頂点を独立に確率$q:=\frac{2\log n}{\sqrt{n}}$で選んで得られるランダムな頂点部分集合を$R$とします. このとき$\mathrm{E}[|R|]=2\sqrt{n}\log n$であって, Chernoff boundにより$|R|$は高確率 (具体的には確率$1-n^{-1}$で) $|R|=O(\sqrt{n}\log n)$となります. また, 固定した頂点$v\in V$に対して「$S_v$の全ての頂点が選ばれない」確率は

$(1-q)^{|S_v|}\leq \exp(-q|S_v|)=\exp(-2\log n)=n^{-2}$

なので, 確率$1-n^{-2}$で$R$は$S_v$の頂点を一つ以上含みます. 従って, union boundにより確率$1-1/n$で$R$は hitting set になっています. 以上より, $R$は高確率でサイズ$O(\sqrt{n}\log n)$の hitting set になってます. もし条件を満たすhitting setが存在しないならばこの確率は$0$にならなければいけないので, 存在性を示せたことになります. また, ランダムにとってくればほぼ確実に条件を満たすhitting setになっているので, 乱択を許せば$O(n)$時間でステップ2を計算できます. また, greedyにやれば ($H=\emptyset$からスタートして$S_{v_1}$の頂点を一つ選んで$H$に追加し, 暫定の$H$がhitしていないような$S_{v_2}$を選んで...を繰り返す) 条件を満たすようにできます.

-ステップ3は$O(|H|m)=O(m\sqrt{n}\log n)$時間で計算できます.

-ステップ4はステップ3の情報を使えば計算できます.

-ステップ5は単一始点BFSを$\sqrt{n}$回行うので$O(m\sqrt{n})$時間で計算できます.


以上より, 全体の時間計算量は$O(m\sqrt{n}\log n+n^2)=\tilde{O}(m\sqrt{n}+n^2)$です.


近似率の解析

アルゴリズムの出力を$D=\max\{D_1,D_2\}$とします. 単純のため, $\mathrm{diam}(G)$は3の倍数とします. 直径を達成する二頂点を$a,b\in V$とします (つまり$\mathrm{dist}(a,b)=\mathrm{diam}(G)$).

二つのケースを考えます:

(1). ある $h\in H$, $v\in V$ が存在して $\mathrm{dist}(h,v)\geq \frac{2}{3}\mathrm{diam}$ がなりたつならば, $D_1$が直径の1.5近似になっているので大丈夫です.

(2). 全ての$h\in H$, $v\in V$に対して $\mathrm{dist}(h,v)<\frac{2}{3}\mathrm{diam}$が成り立つとします. ステップ4で計算した頂点$u$に対し, $\ell:=\min_{h\in H}\mathrm{dist}(u,h)$とすると, $u$の定義より任意の$v\in V$は距離$\leq\ell$で$H$内のいずれかの頂点にたどり着くことができます. ここで$\mathrm{dist}(a,h)\leq \ell$を満たす頂点$h\in H$をとります. すると三角不等式より

$\mathrm{diam}\leq\mathrm{dist}(a,h)+\mathrm{dist}(h,b)<\ell+\frac{2}{3}\mathrm{diam}$,

すなわち $\ell>\frac{\mathrm{diam}}{3}$ を得ます (下図). また, $H$ は $(S_v)_{v\in V}$の hitting set であるため, 特に $H \cap S_u \neq \emptyset$ です. 従って (*)$u$ から距離$\leq\frac{\mathrm{diam}}{3}$ なる頂点は全て$S_u$に属しているはずです (行間補足: $S_u$ は $u$から近い順に頂点を取っていて, 距離$\ell>\frac{\mathrm{diam}}{3}$ だけ離れている頂点 ($H$との共通部分) も取っているのでそれより近い頂点は全て$S_u$に入っていなければならない).




ここで, $D_2 < \frac{2}{3}\mathrm{diam}$と仮定します. このとき, $D_2$の定義より, 全ての頂点は$u$から距離$<\frac{2}{3}\mathrm{diam}$となっています. 次に, $\mathrm{diam}=\mathrm{dist}(a,b)\leq \mathrm{dist}(a,u)+\mathrm{dist}(u,b) < \mathrm{dist}(a,u)+\frac{2}{3}\mathrm{diam}$ より, まとめると $\frac{1}{3}\mathrm{diam}<\mathrm{dist}(u,a) < \frac{2}{3}\mathrm{diam}$ を得ます.

最短$ua$経路上の頂点であって, $u$からの距離がちょうど$\mathrm{diam}/3$となる頂点を$p$とします. $\mathrm{dist}(a,u)=\mathrm{dist}(a,p)+\mathrm{dist}(p,u)<\frac{2}{3}\mathrm{diam}$となるため, $\mathrm{dist}(a,p)<\frac{1}{3}\mathrm{diam}$です. さらに, $p$の取り方から主張(*)より, $p\in S_u$です.

三角不等式より,

$$\mathrm{diam}=\mathrm{dist}(a,b)\leq \mathrm{dist}(a,p)+\mathrm{dist}(p,b)<\mathrm{dist}(p,b)+\frac{1}{3}\mathrm{diam} $$


よって, $D_2\geq \mathrm{dist}(p,b) > \frac{2}{3}\mathrm{diam}$ となってしまい仮定に矛盾します. よって$D_2\geq \frac{2}{3}\mathrm{diam}$ となります (証明終).


ちょっとした改善


最後に, Roditty と Williams [2] がによるちょっとした改善を紹介します. 基本は上記のAingworthらのアルゴリズムと同じですが計算時間を$\tilde{O}(m\sqrt{n}+n^2)$ から $\tilde{O}(m\sqrt{n})$にしています. つまりグラフが疎であれば高速になります.

アイデアは単純で, ステップ1を行わず, ステップ2で得られた hitting set $H$ を, 全頂点を独立に確率$2\log n/\sqrt{n}$で選んで得られるランダムな頂点部分集合とします. そしてステップ5では$u$に対して$S_u$を計算するだけです. Aingworthらのアルゴリズムではステップ1の計算で$O(n^2)$時間がかかっていましたが, そもそも$H$が $(S_v)_{v\in V}$ の hitting setであれば良いだけなので, 計算時間の解析でも紹介したようにランダムにとるだけで良いので省くことができます.


感想

賢すぎて草


参考文献

[1] D. Aingworth, C. Chekuri, P. Indyk, and R. Motwani. Fast estimation of diameter and shortest paths (without matrix multiplication). SIAM J. Comput., 28(4):1167–1181, 1999 (Preliminary version is in SODA96).

[2] L. Roditty and V. V. Williams, Fast Approximation Algorithms for the Diameter and Radius of Sparse Graphs, STOC13.


2021年9月17日金曜日

3SUMのLinear Decision Tree Complexity

Linear Decision Treeと呼ばれる計算モデルについて解説します.

1. 決定木 (Decision Tree)


「ソートアルゴリズムは$\Omega(n\log n)$時間かかる」という話をよく聞くと思いますが, これは「ソートアルゴリズムは2要素の比較を$\Omega(n\log n)$回しなければならない」ことを意味します. 具体的には各ノードに「$a_i\leq a_j$?」と書かれた二分木があり, YESならば右の子ノード, NOならば左の子ノードに辿っていき葉に至ったらそこで終了しそこに書かれているラベル (例えば「$a_3\leq a_1 \leq a_2 \leq a_4$」のような要素の順位が書かれている)を出力するという計算モデルを考えています. このように二分木にYES/NOの質問が書かれており, その答えに対応する子ノードを辿っていき葉に至ったら終了するという計算モデルを決定木 (decision tree) といいます.


例えば上記の決定木では, どんな入力が与えられても質問に3回答えれば葉に至ります. このように. 決定木の効率性は主にその深さによって測ります (もちろん全体のノード数が小さい方が良いですが). ここで, 深さ$D$の決定木の葉の個数は$2^D$くらいです. つまり深さ$D$の決定木は高々$2^D$通りの出力しか区別することができないので, 出力しうるラベルの種類が$N$であるような問題を考えるとその問題を解く決定木の深さは$\lceil\log_2 N\rceil$必要であるという情報論的な下界が得られます. 例えば$n$要素のソートでは$n!= 2^{\Theta(n\log_2 n)}$通りのラベルを区別しなければならないため深さは$D=\Omega(n\log n)$となります. ある問題に対してその問題を解く決定木の深さの最小値を決定木複雑性 (decision tree complexity) といいます. ソートの決定木複雑性は$\Omega(n\log n)$ですが, 実際にマージソートなどは高々$O(n\log n)$回の数値比較でソートをするので, この下界は定数倍を除いてタイトとなります.

決定木の深さはちょうどクエリに回答する回数に対応するため, 計算時間に対応していると解釈できるわけです.

2. 線形決定木


一般に決定木のノードに書かれるラベルはどんな質問でも構いません. 極端な例で言うとノードのラベルにSATのインスタンスが書いてあってそれがYESインスタンスならば左の子に, そうでなければ右の子に辿るようなものを考えてもよいです. しかしながらそれでもなおソートの決定木複雑性は$\Omega(n\log n)$になります.

決定木のノードラベルの質問を線形不等式「$c_1a_1+c_2a_2+\dots+c_na_n\geq 0?$」の形に制限したものを線形決定木 (linear decision tree) といいます. ここで$a_1,\dots,a_n$は入力, $c_1,\dots,c_n$は適当な定数とします. また, 全てのノードのラベルにおける係数に含まれる非ゼロ要素が高々$k$個であるとき, 特に$k$-線形決定木 ($k$-LDT) といいます. 例えばソートの二要素の比較は$a_i-a_j\geq 0$と書き換えることができるので2-LDTとなります.

幾何的な解釈としては, 入力全体の空間を超平面で分割 (線形不等式クエリに対するYES/NOの応答に対応)していき「ラベル1を出力すべき入力の全体」「ラベル2を出力すべき入力の全体」...を区別していくことになります.



3. 決定木とアルゴリズムの違い: uniform vs nonuniform


一つの決定木は固定サイズの問題しか解くことができません. 従って決定木を使って問題を解く場合は決定木の族$(\mathcal{T}_n)_{n\in\mathbb{N}}$を考えます. 問題Pのサイズ$n$の任意の入力$x$に対し, $x$の答えを$\mathcal{T}_n$を使って計算できる時に$(\mathcal{T}_n)$は問題Pを解くといいます. つまり, 与えられる問題のサイズに応じて異なる決定木を使って良いということになります. 

一方でアルゴリズムを使って問題を解く時は入力サイズがどのような場合であっても必ず同一のアルゴリズムを使わなければなりません. このように, 決定木とアルゴリズム(正確にはチューリングマシン)の間には問題サイズによって異なるものを使って良いかどうかという重要な違いがあります. 決定木のように問題サイズごとに異なるものを使う計算モデルをnonuniform, チューリングマシンのようにどのような問題サイズであっても同じものを使う計算モデルをuniformといいます.

UniformとNonuniformな計算モデルのギャップは未解決な部分が多いですが, 少なくともギャップが存在することが分かっています. 具体的にはnonuniformな計算モデルを使うと停止性判定問題を解くことができます (以前の記事参照). ちなみにnonuniformな計算モデルの計算量として最もよく研究されているのは(私の主観ですが)回路計算量だと思います.

今回の記事では線形計算機とアルゴリズムのギャップを, 3SUMという問題に焦点を当てて見ていきます.

4. 3SUMのLDT複雑性


3SUMという問題は

$n$個の数字$A=\{a_1,\dots,a_n\}\in\mathbb{R}$が与えられた時に$a_i+a_j+a_k=0$となる$i,j,k$が存在するか?

を判定する問題です. $a_i+a_j$を格納したサイズ$n^2$の配列Xを用意してこれをソートした後, 各$a_k$に対して$X$が$-a_k$を含むかどうかを二分探索でチェックすれば$O(n^2\log n)$時間で解くことができます. また, 前回の記事で紹介したように$O(n^2)$時間アルゴリズムも存在します. しかしながら任意の定数$\epsilon>0$に対し$O(n^{2-\epsilon})$時間アルゴリズムが存在するかどうかは未解決で, 精緻な計算量 (fine-grained complexity)の分野では$n^{2-o(1)}$時間かかると予想されています. とりあえず詳しくない人は「3SUMを$n^{1.99999}$時間で解くアルゴリズムは知られていない」と思っておけば大丈夫です.

本記事では3SUMのLDT複雑性に関する以下の驚くべき結果の証明(概略)を説明します.

定理1 ([1])
3SUMを解く深さ$O(n^{1.5}\log n)$の4-LDTが存在する.

実は3SUMを解く任意の3-LDTは深さが$\Omega(n^2)$あることが知られていたため[2], 4-LDTにすると一気に深さが減らせるというのはLDT界隈の人にとっては驚くべき事実だったことは想像に難くありません. 何よりも「3SUMは$O(n^{1.5}\log n)$回の要素の比較で解ける」ということになるわけです.

しかし定理1の証明はかなり単純かつ非常に賢いアイデアに基づいています. 定理1の4-LDTは以下のステップからなります.

ステップ1. 入力$A=\{a_1,\dots,a_n\}$を昇順にソートする.
ステップ2. $A$の各要素を連続した$\sqrt{n}$個ずつに分割してそれらを$A_1,\dots,A_{\sqrt{n}}$とする.
ステップ3. $L:=\cup_{i\in[\sqrt{n}]}\{x-y\colon x,y\in A_i\}$をソートする.
ステップ4. 以下のしゃくとり法っぽい計算を行う. ただし文中に出てくる$A_{st}$は$A_{st}:= \{x+y\colon x\in A_s,y\in A_t\}$とする.

for a in A:
    $s \leftarrow 0$
    $t \leftarrow \sqrt{n}$
    while $s\leq t$:
        if $-a\in A_{st}$ then:                            (1)
            return True
        elif $-a > \max(A_s)+\min(A_t)$:            (2)
            $s\leftarrow s+1$
        else:                                                     (3)
            $t\leftarrow t-1$
return False

正当性 (イメージ).


(1)において$-a\in A_{st}$となるならば, ある$b\in A_s,c\in A_t$に対して$-a=b+c$, すなわち$a+b+c=0$となるので, 確かにYESインスタンスとなります. 逆に, 入力がYESインスタンスであるとします. このときwhile文以下のループ処理は以下のような幾何的な解釈ができます:

入力$A=\{a_1,\dots,a_n\}$に対し, 点集合$P=\{(a_i,a_j):a_i,a_j\in A\}$を$xy$平面にプロットする. 各$a\in A$に対して直線$x+y=-a$上に乗っているような$P$の点の有無を判定していきます. ここで, $\mathbb{R}^2$を図のよう長方形で区切っていき, 各グリッドが$P$の点を$n$個ずつ含むようにします (実際には図のようなグリッドにはならないとは思いますが, 説明の簡単のためグリッドで区切っています). すると, while文の最初の時点では$s=0,t=\sqrt{n}$で点$(b,c)\in A_s\times A_t$に着目するわけですが, これは最も左上にある長方形に着目することに対応します. (1)では着目している短形内に直線上の点があるかどうかを判定しています. (2)と(3)では短形内の右下にある点に着目し, この点と直線の位置関係によって右隣か下の短形どちらに進むべきかを判断しています.
3SUMのLDTのしゃくとり法のイメージ (右の吹き出しの$\sqrt{n}$はタイポで正しくは$n$です)


深さの解析.


- ステップ1はマージソートをそのまま実装すれば良く, 深さ$O(n\log n)$の2-LDTで実装できます.

- ステップ2は分割しただけです.

- ステップ3はステップ1と同様にソートするだけです. $L$の要素数は$n^{1.5}$なので, 深さは$O(n^{1.5}\log n)$となります. また, $L$の各要素は$a_i-a_j$の形で書けるため, 4-LDTとなります.

- ステップ4では, while文が$s\leq t$である限り続けますが1回の反復で$t-s$は1ずつ減少するため, 高々$\sqrt{n}$回の反復で終了します. 各反復では(1)の処理がボトルネックになります. もし各$A_{st}$があらかじめソートしてあるならば二分探索ができるので(1)は$O(\log n)$時間で行うことができます. しかしながら全ての$A_{st}$をソートしようとすると, $0\leq s\leq t\leq \sqrt{n}$に対し各$A_{st}$は要素数$n$を持つため, $O(n^2\log n)$時間または同程度の深さのLDTが必要になってしまいます.

ここで, ステップ3で$L$をソートしているという事実に着目しましょう. $L$の全ての要素の大小を比較しているため, 各$i,j$と$a_i,a'_i\in A_i$, $a_j,a'_j\in A_j$に対し
$$a_i-a'_i \leq a_j-a'_j$$
またはその逆方向の不等式どちらが成り立つかが分かっています. これを移項すると
$$a_i+a'_j \leq a'_i+a_j$$
となりますが, この両辺はともに$A_{ij}$の要素です. これら4つの要素は任意にとって良かったので, 
$L$の大小関係の情報が確定すれば各$A_{ij}$の全ての2要素の大小関係の情報は必ず確定します.
決定木でステップ1,2,3を経たノードにいるならば, $L$の要素の大小の順位の情報が確定しているということになるわけなので, 必然的に全ての$A_{ij}$の大小関係の情報が確定しており, 実質的に全ての$A_{ij}$がソートされているとみなすことができます. そこでこの情報を使ってステップ4で二分探索を行うことができるので, 全体の深さは$O(n^{1.5}\log n)$となります.

5. 何が肝なのか?


4章で紹介したLDTでは, ステップ3で計算した$L$の大小関係を使ってステップ4(1)で行う二分探索のために必要な$A_{st}$の要素の大小関係の情報を得られる, という部分が本質的に最も効いています. あくまでも「得られる」というわけであって「簡単に計算できる」というわけではないのがミソです。例えばソートアルゴリズムをそのままLDTに翻訳した場合は「$a_i-a_j\geq 0?$」というノードに対してインデックス$i,j$を簡単に計算して求められていたのですが, 今回考えたLDTではこのようなインデックスは存在性が保証されているだけであって簡単に計算できるとは限らないという点がポイントです。このように存在性と計算可能性のギャップがチューリングマシンとnonuniformモデルのギャップに繋がってくるのは面白い観点だと思いました。


参考文献

[1] Allan Grønlund and Seth Pettie, Threesomes, Degenerates, and Love Triangles, JACM, 65(4), pp.1--25, 2018 (preliminary version is in FOCS14).
[2] Jeff Erickson, Bounds for linear satisfiability problems, Chicago J. Theor. Comput. Sci.,1999(8), 1999 (preliminary version is in SODA95).

2021年8月27日金曜日

小記事: 3SUMのO(n^2)時間アルゴリズム

 3SUMと呼ばれる次の判定問題を考えます:

入力 : $n$個の整数値 $A=(a_1,\dots,a_n)$ ただし$a_i\in[-n^c,n^c]$.

出力 : YES iff 相異なる三つ組 $(i,j,k)$ が存在して $a_i+a_j+a_k=0$


この問題は


for i in range(n):

    for j in range(n):

        for k in range(n):

            if a[i]+a[j]+a[k]==0: return YES

return NO


というコードによって$O(n^3)$時間で解けますが, 少し工夫して次のようにすれば $O(n^2\log n)$ 時間で解けます.


L = []

for i in range(n):

    for j in range(n):

        L.append(a[i]+a[j])

L.sort()

for k in range(n):

    if -a[k] is in L:     // 二分探索をして-a[k]がLにあるかどうかを判定

        return True

return False


しかし尺取り法に基づいた次のアルゴリズムを用いると$O(n^2)$時間で解くことができます (簡単のため, 入力配列$A$の中身は全て相異なる要素であるとします).

A.sort()

for i in range(n):

    l = 0

    u = n-1

    while True:

        if u==l: return False

        if A[i]+A[u]+A[l]==0: return True

        elif A[i]+A[u]+A[l]>0: u-=1

        elif A[i]+A[u]+A[l]<0: l+=1

            

3SUMという問題は計算幾何学や精緻な計算量理論において非常に重要な問題の一つとされていて, 任意の定数$\epsilon>0$に対して$O(n^{2-\epsilon})$時間では解けないと予想されています.

一方で$\mathrm{polylog}(n)$倍の改善はなされています.


2021年7月27日火曜日

講義を完走した感想

人生で初めて講義(教える側)をしたのですが、とても濃い経験になったので備忘録として残そうかと思います。大学で講義を「する」という経験は一般的に見ればなかなかない経験であまり想像がつかないと思いますので、どういう感じなのかを説明します。講義資料の準備で意識したことや悩んだことなども綴っていくので、これから講義をする機会のある方には参考になる部分もあるんじゃないかと思います。既に幾つも講義を受け持っているプロフェッショナルな方々には自身の経験を思い出して温かい目で見守っていただけたら幸いです。コロナ禍という特殊な環境ではあったものの初の対面講義で奮闘した一人の助教の記録です。


0. 背景

自分が担当した講義はプログラミングの発展的な内容を扱うもので、学部2年生向けの講義です。まだ学部2年生の前期なのでまだあまり専門分野を学んでいない状態です。とはいえプログラミングの基礎的な部分はある程度学んでいてPythonは少しは分かるという状態です。また、必修科目というわけでもないので、受講者はある程度プログラミングに興味を持っているという想定です。また、ORなどに興味を持ちそうな理系向けの講義になっています。


1. 講義内容と形式の決定

自分が最初に苦労したのはこの部分です。特に何をすべきかが厳密に定まっておらず、内容はほとんどを任されていましたが、昨年度までの資料を見る限りでは基礎的なアルゴリズムをやっていたためその方向性は守ろうと思いました。他の先生に話を聞くと「ダイクストラ法はやってほしい」みたいな声も聞いたので、この辺りはちゃんと取り込もうと思いました(もちろんダイクストラは元々やるつもりでしたが)。昨年度までは本当に幅広く(平面幾何、乱択アルゴリズム、組合せ最適化など)さまざまなアルゴリズムを軽く紹介するという感じだったので、組合せ最適化の内容を主に扱い、最後に自分の趣味(ランダムグラフ、ランダムウォークなど)を少し取り入れた感じにしようと考えていました。

私の所属する大学はクォーター制度を導入していたため、講義は全7回です (6月頭〜7月末)。週1回講義と演習の時間がそれぞれ100分ずつあります(なので計200分あります!)。また、2020年度はコロナ禍の影響でオンライン講義になっていたようですが、私が赴任した2021年度は演習ということもあり対面授業になるとのことでした。

私がD3だった2020年度はコロナ禍が始まった年ということもあり仕方のないことですが対面講義と比べるとオンライン講義はあまり楽しいものではありませんでした。なので対面で講義できるというのは嬉しく思いました。特に受講生の顔で反応を見ながらインタラクティブに教えたりするのは教育の楽しい部分だと個人的に思います。

とはいえ緊急事態宣言が発出されている状況なのでいつオンライン講義に変更になっても対応できるように、昨年度と同様に google colab (jupyter notebookをブラウザ上で動かすgoogleのサービス) のドキュメントを講義資料とすることにしました。特にgoogle colabだと学生のPCにpythonを動かすためのソフトを入れる必要もなく, googleアカウントさえあれば誰でも動かせるのもあって導入が非常に簡単です (今の時代だとほとんどの人がYouTubeを観たりするのでgoogleアカウントは持っています). ipynbファイルを講義資料として配布し、受講生は自分のPCを持参してもらえれば普通に講義を行うことができます。

2. 講義資料の準備 (第1回〜第3回)

4月に着任して早々に講義をどうするかが頭をよぎったので、かなり前もって準備をすることにしました。今思えばこの選択はかなり良かったです(何事も余裕を持って準備しておくことは非常に大切です)。

講義資料を作る前に自分の中ではグラフアルゴリズムをメインに据えるということを意識しました。グラフ理論やグラフアルゴリズムは自分が非常に得意とする分野であり、何度も実装したことがあるからです。そこで、当初は次の二つの事項を講義で扱おうと考えました:

1. Pythonを用いて基礎的なアルゴリズムを理解

  キーワード: 再帰, DFS, BFS, 貪欲法, 最短経路問題 (Dijkstra, Warshall—Floyd), DP.

2. 応用を知る: ネットワーク解析

  キーワード: ランダムウォーク, ランダムグラフ (BAモデル),  ネットワーク中心性 (PageRank, 媒介中心性) など.

ひとまず第1回の講義はPythonの復習を行い、その後は基礎的なアルゴリズムをひたすら紹介していき、最後の第7回講義で自分の趣味を交えた講義を行うという計画を立てました。

この方向で講義を組み立てて実際に講義資料を作り始めたのが4月の中旬 (講義開始2ヶ月前) でした。初回の講義資料が完成したので見返しました。初回はPythonの文法の復習だけやって「1から1000までの7の倍数の和を求めるコード」とかを1時間やるものだったので、純粋に「これはツマランな」と心配になりました。そこで少しでも退屈を紛らわそうと再帰関数の紹介を行いました。再帰を使ってフィボナッチ数列や階乗の計算のプログラムを書き、最後に非自明な例として分割数の漸化式を示しつつ分割数を計算するコードも紹介しました。

第1回の課題は「与えられた自然数が素数かどうか判定せよ」みたいな問題を出していましたが、なんとなく味気ないのでコラッツ予想を絡めた問題も出題してみました。こんな感じです。


この調子で第2回講義資料も取り掛かり、計算量とグラフの定義を行いました。グラフの隣接行列や隣接リストを紹介し、「グラフの次数列を計算せよ」「三角形(長さ3の閉路)の個数を求める関数を書け」といった課題を出していました。

課題を考える時に意識したのはデキる人でも退屈しない内容にしたいということです。講義そのものは流石に集中が散ったりするのでしょうがない部分もあるかもしれませんが、せめて競技プログラミングのように、問題を解く部分に楽しみを見出して楽しんでもらいたいと思っていました。第2回講義以降は毎週の課題の中に一つだけ学部2年生向けの講義にしてはありえないくらレベルの高い難問を設けることにしました (難問はチャレンジ課題枠としてあり、解けなくても成績に影響は出ないことを明言した上で出題しました)。

ちなみに第2回講義で出題した難問は以下の問題です:

要するに長さ4の閉路の判定を$O(n^2)$時間で行えという問題で、この記事にも紹介したように鳩の巣原理を使うと計算量が抑えられるという問題です。これを計算量の概念とグラフの定義を学んだばかりの学部2年生に出すのは正気の沙汰とは思えません。こんなん知らない状態から解けるやつおらんやん...

第1回と第2回の講義資料と課題の準備が一通り終わったので、他の先生や学生TAに難易度や分量を含めて確認してもらい、大丈夫そうとのことだったのでここでようやく一息つきました。ここまでの講義の内容は必須事項だったのであまり迷いはありませんでしたが、以降の構成に頭を悩ませることになりました。

ひとまず考えた計画としては残りの5回で

・可視化ライブラリ(matplotlibとnetworkx)の紹介、幅優先探索、深さ優先探索
・貪欲法、クラスカル法
・動的計画法
・Dijkstra、Warshall--Floyd、Bellman--Ford
・フロー (Ford--Fulkerson) 

をやるか、または最後のフローはやらずに代わりに自分の趣味(ランダムグラフなど)を紹介するか、どちらかにしようと決めました。

いずれにしてもグラフアルゴリズムをやることは決まっていたので、せっかくなのでアルゴリズムの可視化を行うことを思いつきました。ひとまずgoogle colab上でのアニメーションの方法を学び、試しにDFSやBFSのアニメーションを書いて遊びつつ講義資料の準備を進めました。

ちなみにこんな感じのアニメーションを作りました。受講生はgoogle colabで開いている講義資料にあるPythonコードを実行するとこのようなアニメーションが出てくるというようになっています。アルゴリズム自体を可視化しているのでグラフは自由に変えられます。また、ボタンを押した時に一コマだけ進ませるというようなことも可能です。


ちなみにこの動画のソースコードは本記事の5章で公開しています。

この調子でクラスカル法やダイクストラ法もアニメーションプログラムを作っていこうと決めていくうちに、とうとう初回講義の日がやってきました。

3. 対面講義初回


結論からいうと初回講義は機器の不調なども重なって非常に大変でした。そもそも自分の担当は講義100分+演習100分だったのですが、実は当時の自分は勘違いをしていて講義100分の中に演習が組み込まれているものと思っていました(なので実際の講義は想定の倍あったということになります)。

元々70分程度の想定の下で準備していた講義を緊張もあって少し早口で進めていき、めちゃ早く終わっちゃうと内心焦っていたら、突然PCの調子が悪くなり学生TAのPCを急遽借りて講義を再開して案の定めちゃくちゃ早く講義部分が終わって演習時間がものすごく長くなるなど、思い返しても大変な1日になってしまいました。結局初回は課題が簡単だったのでほとんどの学生がすぐに帰る感じになっていました。

今にして思えばもっとやりようはあっただろうとも思えますが、2ヶ月弱前から準備していたしやれることはやったのであまり気にせずとりあえず講義を進めていくことになります。

第2回講義は使用PCも変えたので第1回と比べるとかなりスムーズに終えられました。第2回課題には第1回とは違い難問も用意していましたが、これは予想以上に反応が良かったです。計算量を気にしてプログラムを書くということをあまり経験してこなかった学生がほとんどだと思うので四角形の高速な判定にどこまで取り組んでくれるか心配したのですが、10人弱の学生が18時くらいまで残って議論しながら取り組んでくれました。これには教員側からしてもとても嬉しかったです。


4. 講義の改善点の洗い出し


まだ2回ではあるもののここまで講義をやってきていくつか変更点が見えました。まず、そもそも講義時間を勘違いしたまま準備していたので単純問題として内容は倍にできる。どうするか?

また、学生の質問を受け付けていくうちに色々改善した方がいいかもなぁというものも多く出てきました。例えば

  • google colabのインデント幅の設定が人によって違う
  • 自身が書いたコードがちゃんと正しく動いてるかを確認するのが面倒(学生が自分でテストケースを用意して確認しなければならないが、グラフアルゴリズムだと自分で隣接リストなどを用意しなければならない)
  • printfデバッグみたいなことを教えた方がいいのか?
などです。printfデバッグは質問してきた学生に「こうしてprintをするとどこで詰まったか分かるよ」と直接教えることにしました。

インデントなど、全体に連絡すべき事項はちゃんと講義冒頭にスライドを使って説明するようにしました。

課題プログラムの自己チェックを可能にするために、Checkerというものを導入しました。これはつまるところ競技プログラミングのサンプルケースとその答えを事前にこちらで用意しておき、学生が書いたコードがこれらのサンプルケースで正しい答えを出すかどうかをチェックするプログラムを講義資料に入れておいたのです。こんな感じです。





*)完全に問題がAtCoder感があるかもしれませんが、こういう文章題は少なくて、普通の問題はグラフ理論の用語を使った問題文になっています。

また、課題の採点を楽にするためにテストケースを用意しました。基本的にはテストケースで正解していて大丈夫そうなコードを書いていたらその課題を満点にするという感じにしました。これのおかげで採点は劇的に楽になりました。

もちろん、サンプルケース, Checker, テストケースによる自動採点の準備はかなり大変でこれを講義資料の準備と並行して行うのは非常に時間がかかり大変でした。しかし、特にCheckerの導入は受講生にはかなり高評価だったようで、さらに「Checkerを走らせたらエラーが出た/Wrong answerだった」という類の質問がきてくれたことは良かったです(Checker導入以前はあまり質問をしてくれる学生がいなくて寂しかったです)。やはり、自分の書いたプログラムがちゃんと動いているかの確認が簡単にできるということは準備の大変さに目を瞑ればメリットしかありません。


5. 以降の講義資料の準備と講義


ここまで来ると講義自体も少しずつ慣れてきて、講義内容も固まってくるので特に悩むこともなくなってきます。以前立てた計画通り

・貪欲法、クラスカル法
・Dijkstra、Warshall--Floyd、Bellman--Ford
・動的計画法
・フロー (Ford--Fulkerson) 

の順番で講義資料を準備していきました。結局最終回はフローと最大二部マッチングをやることにしました。

ちなみに講義資料のために作ったアルゴリズムのアニメーションは

https://colab.research.google.com/drive/1ed0fAC6rUMcScJYxv0Kc0Ox-wDBf2FSR

で閲覧することができます(ダウンロードすれば誰でも編集とかができるので、講義などで使ってみたいという方ご自由にどうぞ。ソースコードの下から二行目のanim.save(...)の部分をコメントアウトすると動画をmp4ファイルで保存できます)。

ちなみに細かい講義内容は

第4回:

  • 貪欲法 (コイン支払い問題, 区間スケジューリング問題)
  • 最小全域木問題
  • クラスカル法


第5回

  • 最短経路問題
  • ダイクストラ
  • ベルマンフォード
  • ワーシャルフロイド


第6回

  • 動的計画法
  • 最長共通部分列問題
  • ナップザック問題
  • 重み付き区間スケジューリング問題
  • Held—Karp のTSP


第7回

  • フロー
    • 残余グラフ
    • Ford—Fulkerson
    • Edmonds—Karp
  • 二部マッチング
    • 割り当て問題

という感じになりました。難問枠としてmeet-in-the-middleなどを出題しました。

反省点は幾つもあります。
  • 第5回以外は意外と早く講義が終わってしまったというのがあります。アルゴリズムを広く浅く紹介して課題で問題を解きながら理解を深めてもらうというスタンスだったので、内容としてはちょうど良いのかもしれませんが講義自体は説明は割と早めに終わってしまいました。この辺りは説明をもっとゆっくりやればちょうど良くなりそうな気配があるので、経験を積むと感覚が掴めてくるのかもしれません。
  • 初回の講義についてはPythonの復習にたくさん時間を割くよりも、エラトステネスの篩、ユークリッドの互助法、二分探索などの軽いアルゴリズムの紹介などが出来た気がしています。復習に時間を割くよりもコードで出てくる度に「リストはこうやって要素を追加するよ」と言えば良さそうです。
  • 課題の提出についても注意喚起が足りなかった部分があります。課題を解いたのに保存で失敗して白紙の状態で提出してしまったという学生もいました。こういった不運な事故を防ぐ努力(注意喚起)はもうちょっとすべきだったと思います(流石に初めての講義でそこまで気が回らなかったのもしょうがない気がしますが。。。)

このように、色々と改善点はまだまだあるものの、初めてやったにしては上出来なんじゃないかと開き直ることにしています。

なお、講義資料の準備段階では他の先生の講義動画などが非常に参考になりました。説明の順序やその問題を考える動機などが自然に導入されており、流石だなぁと思いました(もはやこの動画が講義資料でいいじゃんと思っちゃうくらい)

6. 感想


元々自分は人に教えるのが好きだったので、とても楽しく完走することができました(もちろん準備はかなり大変でしたが。。。)来年度以降はゼロからアニメーションコードを開発したりまではしなくて良いのである程度苦労は減ることを期待していますが、今年度できなかった内容も取り入れたいし、実は少し楽しみにしています。例えば乱択アルゴリズムを取り入れて難問枠としてcolor codingを入れたりとかを妄想しています笑

一方でコロナ禍での対面授業ということもあり、感染拡大の防止にものすごく気をつけました。対面講義をするというのもあって、4月の緊急事態宣言以降は今もなお一度も人と一緒に外食をしておらず、そもそも食事はほとんど全てをテイクアウトか自炊で済ましています(講義が終わった今後も同様の生活を続けていくつもりです)

受講生自身からも対面授業が嬉しいというコメントがきている一方で彼らも感染防止にはかなり気を遣ってくれていて、もちろん完全に会話をゼロにするには流石に無理ですが、それでも常にマスクをつけて暑い中換気のために窓を開けるのに協力してくれたり、入室時には必ず手の消毒をしてくれていました。

私が担当した2年生は昨年度はコロナ禍でほとんど大学に来れなかったということもあり、対面授業で盛り上がって大声で話し出したらどうしようと危惧していて(実際そうなったとしたら流石に注意せざるを得ないが彼らの心情を考えると非常に心苦しい)悩んでいた時期もありましたが、基本的には近くの2~3人組で相談し、意欲のある学生は難問の相談を行う感じになってくれており、当初危惧していたようなこともなかったため非常に安堵しています。

私の周りや受講生で感染者が出たという話は今のところなく、おそらく2週間後くらいまでは完全に安心できないですが、ひとまず講義中にこのような問題が起きなかったことが一番ホッとしています。



2021年4月2日金曜日

次数分離: スパースなグラフで役立つアルゴリズムのアイデア

単純無向グラフ$G=(V,E)$上で長さ$k$の閉路を探したり直径を求めたくなったりする時があります。そして時には$G$はスパースになっていることを利用したくなります。今回の記事では、スパース性を利用した単純ながら面白く賢いトリック (次数分割: high-degree low-degree technique)を紹介し、これに基づいてグラフに含まれる三角形 (長さ3の閉路) の個数の数え上げを行う非自明なアルゴリズムを紹介します。次数分離という和名は私が勝手に名付けたものであり世間で浸透しているものではないので注意してください (もし世間で浸透している名前をご存知の方がいらっしゃいましたら教えていただけると嬉しいです...!!)。このアイデアはAlon, Yuster, Zwick [1] によるものです。

$|E|=m$, $|V|=n$とします。一般には$m=O(n^2)$になりえますが、今回は$G$がスパースであることを仮定してるので$m=O(n)$を想定します。ただ次数分離は一般の$m$に対しても適用でき、例えば$m=O(n^{1.5})$とかでも似たような解析で計算時間を見積もることができます。

次数分割の基本アイデア


頂点$v$の次数を$\deg(v)$とします。$\alpha\in[0,1]$をパラメータとし、
$$H=\{v\in V: \deg(v)>n^\alpha\},$$
$$L=V\setminus H$$
とします ($H$は high degree、$L$は low degree のイメージです)。すると、$2m=\sum_v \deg(v)$ なので
$$2m \geq \sum_{v\in H} \deg(v) > |H|n^\alpha$$
より$|H|<2m/n^\alpha$となります。特に$m=O(n)$ならば$|H|=O(n^{1-\alpha})$となり、$n$と比べて小さくなります。これらの頂点部分集合は$O(m)$時間で得られます (計算モデルは$O(\log n)$bit-Word-RAM モデルを仮定)。


応用. 三角形(長さ3の閉路$C_3$)の数え上げ


入力としてグラフ$G=(V,E)$が与えられ、$G$に含まれる長さ$3$の閉路$C_3$と同型な部分グラフを数え上げる問題を解いてみましょう。以下の二つの素朴なアルゴリズムを考えます。

・アルゴリズム1:

$G$の隣接行列の三乗を計算し対角成分の総和を見ると、長さ3の閉じたウォークの個数となり、これは三角形の個数の$3!=6$倍になっています。従って$O(n^{\omega})$時間で数え上げられます。 ($\omega<2.373$は高速行列積の指数)。

・アルゴリズム2:

各頂点に対して隣接点を列挙して隣接点のペアで隣接しているものがあるかどうかを見れば$O(\sum_{v\in V}\deg(v)^2)$時間で数え上げが解けます。$(2m)^2=\left(\sum_{v\in V}\deg(v)\right)^2\geq \sum_{v\in V}\deg(v)^2$より、このアルゴリズム2は$O(m^2)$時間なので、例えば$m=O(n)$なら$O(n^2)$となりアルゴリズム1以上に高速になります ($\omega=2$なら同等のオーダーになる)。

ここで、冒頭に説明したアイデアを考えます。適当なパラメータ$\alpha\in[0,1]$に対して上で定義した頂点部分集合$H,L$を得ます。入力グラフ$G$に含まれる三角形$C_3$は以下の二つのケースに分類されます:

ケース1. 三頂点全てが$H$に含まれる場合
ケース2. 一つ以上の頂点が$L$に含まれる場合

ケース1は誘導部分グラフ$G[H]$内に含まれる三角形の個数を数え上げれば良く, $|H|=O(n^{1-\alpha})$より素朴なアルゴリズム1を用いれば$O(n^{\omega(1-\alpha)})$時間で数え上げることができます。

ケース2は$L$に含まれる全ての頂点に対し素朴なアルゴリズム2を用いることで列挙できるので、数え上げることができます。このとき、$k$頂点が$L$に含まれるような三角形 $(k=2,3)$ は$k$回列挙されるため適当な値で割ったりする必要があるので気をつけてください。列挙にかかる時間は、$L$に属する頂点の次数が高々$n^\alpha$であるので $\sum_{v\in L}\deg(v)^2 \leq n^{1+2\alpha}$ となります。

以上より、$m=O(n)$の場合は任意の$\alpha\in[0,1]$に対して三角形の数え上げを
$$O(n^{\omega(1-\alpha)}+n^{1+2\alpha})$$
時間で解くことができます。$\alpha=\frac{\omega-1}{\omega+2}$とおくとこの計算時間は$O(n^{3-\frac{6}{\omega+2}})$時間になります。現状では$\omega<2.373$なのでこの計算時間は$O(n^{1.628})$となり, $O(n^2)$よりも高速になっています。もし将来$\omega=2$になることがあれば$O(n^{1.5})$となります。

上の解析では$m=O(n)$を仮定していましたが一般の$m$に対しても同様の解析で$O(m^{3-\frac{6}{\omega+2}})$時間で計算できることを示せます ($|H|=O(m^{1-\alpha})$となるので)。


まとめ


スパースなグラフに対して、ある閾値を設定しそれより高い次数を持つ頂点集合とそれ以外の頂点集合に分割するテクニックを紹介し、その応用として三角形の数え上げ問題を取り上げました。本質的には、このような分離を行い、低次数の頂点の探索は近傍を見るコストが低いことを利用し、高次数の頂点の探索はまとめて処理できる賢いアルゴリズムを利用する、というだけのことです。

三角形と同様の解析を用いれば、例えば入力グラフが四角形(長さ4の閉路)を部分グラフとして持つかどうかの判定問題を$O(m^{4/3})$時間で組合せ的 (高速行列積という飛び道具を使わず初等的に)解くことができます。$H$と$L$に分解したあと、$H$内部の四角形を鳩の巣原理に基づくアルゴリズム (昔の記事参照) で判定し、$L$の頂点を含む四角形は近傍探索を駆使して判定します。結構実装も簡単でシンプルかつ汎用性がありそうなアイデアなので、競技プログラミングで出ても良さそうな気が個人的にはします(というか既にある??)


参考文献


[1] N. Alon, R. Yuster, and U. Zwick, Finding and counting given length cycles, Algorithmica, 1997.

アルゴリズムの理論研究の最高峰とは?

競プロで道具として用いられる様々な賢いアルゴリズムやデータ構造の多くは, 非常に賢い研究者たちによって発見されており, ほとんどは理論計算機科学の論文として国際学会の会議録や学術雑誌の形として出版されています. ここではアルゴリズム系の研究論文がよく出てくる様々な国際会議を紹介し...