📂位相データ分析

ベルツルアルゴリズム: 最小内包ディスク問題の解法

定義

最小包含円

$n > d$ としよう。$d$次元のユークリッド空間で与えられた有限な集合$P = \left\{ p_{k} \right\}_{k=1}^{n} \subset \mathbb{R}^{d}$に対して、以下のような最適化問題を最小包含円問題^{smallest Enclosing Disk Problem}という。 $$\begin{matrix} \text{Minimize} & r \ge 0 \\ \text{subject to} & \left\| c - p_{k} \right\|_{2} \le r \end{matrix} \\ c \in \mathbb{R}^{d} , k = 1, \cdots , n$$

ヒント

これは正確には定理と呼ぶほどのものではないが、この問題を扱う際に知っておくと良い事実である。

• $P$がアフィン独立である場合、円の境界には$P$の点が最少で$2$から最大で$d+1$個まで存在する。つまり、点が重なっていたり、3点以上が一直線上にあるような場合を除き、正確に$2 \le m \le d+1$個の点で最小包含円が一意に定まる。
  • 例えば、$d = 2$次元平面では、円は正確に3つの点で一意に定まる。
• 一般に$n > d+1$個の点が与えられた場合、アルゴリズム^algorithmではなく明示的な公式^{explicit formula}を見つけることは不可能とされている。境界上に最小包含円を一意に定める点をサポート^supportと呼ぶが、単に点を持っているだけでは誰がサポートかを知る方法がないためである。
  • このように境界上の点をサポートと呼ぶのは、幾何学の問題で一般的である。サポートベクターマシンでは、サポートが境界上の点を意味する。
• したがって、この問題を解くアルゴリズムを開発するということは、境界上のサポートを見つけることを保証するか、それを迅速に見つける方法に関する研究と言っても過言ではない。ただし、現在使用されているアルゴリズムでは、その根本的なアイデアはほぼウェルツルのものに基づいており、その改良や変種をまとめて単にウェルツルアルゴリズム^{welzl Algorithm}と呼んでいる¹。少なくともこの投稿で紹介されている後続の研究はすべてその系統に属している。

解法

ウェルツルアルゴリズム ²

ウェルツルアルゴリズム^{welzl Algorithm}は、最小包含円問題を解く再帰的な^recursive解法である。基本的には、点を一つずつ追加したり削除したりしながらサポートを見つけ、そのようにして得られた円が与えられたすべての点を包含しているかどうかを繰り返し確認する。

点が$n \le d+1$個の場合にそれらを包含する円を正確に見つけることは比較的簡単であるため、そのような関数が存在すると仮定する。実際の実装では、これが思ったほど単純ではないが、最小包含円問題の核心ではない。

擬似コード ³

$(c,r)$ = welzl$\left( P, S \right)$

• Input: 与えられた点の集合$P = \left\{ p_{k} \right\}_{k=1}^{n} \subset \mathbb{R}^{d}$とサポーターの候補$S \subset P$を受け取る。ヒントで述べたように、$\left| S \right| \le d+1$である。
• Output: $P$のすべての点を包含する最も小さい円の中心$c$と半径$r$のタプル$(c,r)$を得る。中心が$c$で半径が$r$の閉じた球を$D = B \left[ c,r \right]$と表す。

$(c,r)$ = trivial$\left( S \right)$

• Input: $\left| S \right| \le d+1$の点の集合$S \subset P$を受け取る。
• Output: $S$のすべての点を包含する最も小さい円の中心$c$と半径$r$のタプル$(c,r)$を得る。welzlに比べて単純であるという仮定に従い、trivialの擬似コードは別途記述しない。

---

function welzl$\left( P, S \right)$
　　$S := \emptyset$
　　if $P = \emptyset$ or $\left| S \right| = d+1$ then
　　　　return trivial$\left( S \right)$
　　else
　　　　choose $p \in P$
　　　　$D :=$ welzl$\left( P \setminus \left\{ p \right\}, S \right)$
　　　　if $p \in D$ then
　　　　　　return $D$
　　　　end if
　　end if
　　return welzl$\left( P \setminus \left\{ p \right\} , S \cup \left\{ p \right\} \right)$
end function

---

welzlは再帰関数として記述されるため、実際のプログラムではtrivialから計算が始まることになる。welzlの擬似コードで登場するchooseは、choose$x \in X$のように書かれ、集合$X$から要素$x$を一様ランダムに選ぶキーワード^keywordである。

最後にwelzl$\left( P \setminus \left\{ p \right\} , S \cup \left\{ p \right\} \right)$をリターンすること自体が、$P$にある点を一つずつ取り除いて$S$に入れてみてサポートを見つけるプロセスであることが分かれば、このアルゴリズムを理解したことになる。

説明

最小包含円の制約条件を文字通りに解釈してみよう。$\left\| c - p_{k} \right\|_{2} \le r$を満たす$c$と$r$を見つけるということは、少なくとも与えられた点をすべて包むことができる中心^centre$c$と半径^radius$r$を見つけることを意味する。しかし、ユークリッド空間は非常に広いため、$c$がどこであっても$r$を無制限に大きくすれば、この制約条件は必ず満たせる。当然のことながら、私たちの関心事は、それをしながら$r$を最小化すること、つまり与えられた点を包む最も小さい球を見つけることである。

単純化の問題点

包含^enclosingは、数学全般で広く使われている表現ではないが、発音そのままにインクロージングと言うにはあまりに長くて感じがこないため、任意に翻訳した。そもそもインクロージング^enclosing自体がバウンディング^boundingを圧倒しているわけではない。また、円盤^diskも内部を含む閉じた球の意味で使われたが、実際の英語表現では単にボール^ballやスフィア^sphereもよく使われる。単語や翻訳にこだわらないようにしよう。

歴史

早くもライムント・ザイデル^{raimund Seidel}によって線形計画法に基づく解法が知られていた。

1991年にエモ・ヴェルツル^{emo Welzl}が彼の論文で再帰的な^recursiveアルゴリズムを提案し、いわゆるSOTA^{state Of The Art}を記録した。2022年現在でも、一般的な最小包含ディスク問題ではこのヴェルツルのアルゴリズムが最も優れているとされており、以降の多くの研究はこれを改良する形で行われてきた。

1999年にはベルント・ゲルトナー^{bernd Gärtner}が、ヴェルツルのアルゴリズムに従いつつも、二次計画法^{quadratic Programming}の応用を導入してこれを改善した。⁴ 彼のコードはC++で書かれており、チューリッヒ連邦工科大学のウェブサイト⁵で実際の実装を見ることができる。

2003年にはカスパー・フィッシャー^{kaspar Fischer}が、ゲルトナーとの研究で線形計画法のシンプレックス法で登場するブランドのルールを導入し、高速なコードを作成した。⁶ 2013年にはトーマス・ラーソン^{thomas Larsson}が、近似的なものではなく、速度と堅牢性^robustnessを備えた方法を提案した。⁷

ここまで紹介された研究を参照すると、ヴェルツル-ゲルトナー-フィッシャーと続く大きな流れを確認することができる。

応用

ヴェルツルアルゴリズムの代表的な応用は、チェックコンプレックスの構築である。