📂位相データ分析

ゾモロディアンのアルゴリズムの導出

概要

ゾモロディアン^zomorodianとカールソン^Carlssonの論文’Computing Persistent Homology’で紹介されたアルゴリズムの導出過程を説明する¹。抽象単体的複体から作られたフィルタード複体を受け取って$\mathcal{P}$-インターバルを返し、コンピュータで扱いにくいパーシステンスモジュールの構築を省略して行列のリダクションで持続性ホモロジーを計算する。

導出

Part 0. 事前調査

アルゴリズムの本格的な導出に先立ち、上のような図で描写されるパーシステンス複体が数式的にどのような形なのかからまず見てみよう。この過程をしっかり固めておかないと、論文を読むのが非常につらいだろう。

• まず下段にある数を$\deg$といい、これが$0$から$5$まで増加し、次のようにフィルタード複体を成す。 $$ \left\{ a,b \right\} = K_{0} \subset K_{1} \subset K_{2} \subset K_{3} \subset K_{4} \subset \left( K_{4} \cup \left\{ acd \right\} \right) = K_{5} $$
• $\deg$と関係なく$K$は$2$-単体であり、ホモロジーを考えるという文脈で次のような鎖複体を成す。 $$ \mathsf{C}_{2} \overset{\partial_{2}}{\longrightarrow} \mathsf{C}_{1} \overset{\partial_{1}}{\longrightarrow} \mathsf{C}_{0} $$

アルゴリズムの目標は、このような$\partial_{2}$と$\partial_{1}$がデータに対して与える代数的位相幾何学的情報、たとえばベッチ数$\beta_{k}$のようなものがどの$\deg$で生じ、いつ$\deg$で消えるのかを次のように計算することである。 $$ \begin{align*} L_{0} =& \left\{ \left( 0, \infty \right) , \left( 0,1 \right) , \left( 1,1 \right) , \left( 1,2 \right) \right\} \\ L_{1} =& \left\{ \left( 2,5 \right) , \left( 3,4 \right) \right\} \end{align*} $$

$L_{0}$は$\beta_{0}$に該当する情報、すなわちコンポーネントがいつ生じいつ消えるかを示す$\mathcal{P}$-インターバルからなり、$L_{1}$は$\beta_{1}$に該当する情報、すなわち空間で「穴」というべきものがいつ生じいつ消えるかの$\mathcal{P}$-インターバルからなる。

---

Part 1. $\partial_{1}$

論文で著者らはそのための計算がすべての体で可能だと主張するが、簡単に次数付き加群である$\mathbb{Z}_{2} [t]$-加群でどのような計算が行われるのかを見てみよう。さて$\mathsf{C}_{k}$の斉次^homogeneous基底を$\left\{ e_{j} \right\}$、$\mathsf{C}_{k-1}$の斉次基底を$\left\{ \hat{e}_{i} \right\}$と表記しよう。ここで斉次であるとは、$\mathsf{C}_{k}$を次数付き加群として見たとき項が一つしかないという意味に受け取っても、すなわち$t^{2} + t$のようなものではなく$t^{4}$のように単項式の形だとみなしてもかまわないということである。

$$ \deg M_{k} (i,j) = \deg e_{j} - \deg \hat{e}_{i} $$ ホモロジー代数にある程度慣れているなら、これから$\partial_{k}$に対応する境界行列$M_{k}$を上の表と方程式に合わせて構成し、そのスミス標準形$\tilde{M}_{1}$を求めに行くのだという感じが来るだろう。まず$k=1$の場合を考えてみると、先ほど行列の基底が斉次だと言ったので、次のように唯一の$M_{1}$を得ることができる。

このように基底を用いて行列を構成することは、$\partial_{k}$の役割の一つが$t^{n}$を掛けること(群作用を取ることによって次数付き加群で次数が上がったこと)を逆に行うことだと見ればよい。感覚をつかむためにちょうど三回だけ直接計算してみよう。 $$ \begin{align*} \deg M_{1} (2,5) =& \deg ac - \deg c = 3 - 1 = 2 = \deg t^{2} \\ \deg M_{1} (4,5) =& \deg ac - \deg a = 3 - 0 = 3 = \deg t^{3} \\ \deg M_{1} (2,2) =& \deg bc - \deg c = 1 - 1 = 0 = \deg t^{0} = \deg 1 \end{align*} $$

先ほど述べたように、これから、これの階段形、特に列-階段形を作れば次のようになる。

学部のときに学んだ線形代数を思い返してみると、各列で最も上にあって$0$でない、図のように四角で囲んだ部分のようなものをピボットと呼んだ。ここで次の二つの補助定理を紹介する。

• (1): 列-階段形の対角成分はスミス標準形の対角成分と同じである。
• (2): $\tilde{M}_{k}$の$i$行のピボットが$\tilde{M}_{k} (i,j) = t^{n}$であればホモロジー群$H_{k-1}$の$\sum^{\deg \hat{e}_{i}} F[t] / t^{n}$に該当するものであり、その他は$H_{k-1}$の$\sum^{\deg \hat{e}_{i}} F[t]$に該当するものである。これは$L_{k-1}$が$\left( \deg \hat{e}_{i} , \deg \hat{e}_{i} + n \right)$と$\left( \deg \hat{e}_{i} , \infty \right)$からなることと同値である。

言い換えると、

• 補助定理 (1)によれば、持続性ホモロジーを計算するときは行演算が必要なく、列演算だけあればよい。
• 補助定理 (2)によれば、$L_{k-1}$は$\left( \deg \hat{e}_{i} , \deg \hat{e}_{i} + n \right)$と$\left( \deg \hat{e}_{i} , \infty \right)$からなる。
  • 一行目のピボットが$t^{1}$で$\deg d = 1$なので$(1,1+1)$を得る。
  • 二行目のピボットが$t^{0}$で$\deg c = 1$なので$(1,1+0)$を得る。
  • 三行目のピボットが$t^{1}$で$\deg b = 0$なので$(0,0+1)$を得る。
  • 四行目のピボットがなく$\deg a = 0$なので$(0,\infty)$を得る。

これはアルゴリズムを導出する前に言及した$L_{0}$と正確に一致する。 $$ L_{0} = \left\{ \left( 0, \infty \right) , \left( 0,1 \right) , \left( 1,1 \right) , \left( 1,2 \right) \right\} $$

---

Part 2. $\partial_{2}$

$L_{1}$を求めるための$\partial_{2}$の行列形$M_{2}$は上のとおりである。ところがここで次の補助定理で計算を減らして楽に進むことができる。

• (3): $\mathsf{C}_{k+1}$の標準基底と$\mathsf{Z}_{k}$に対する$\partial_{k+1}$を表現するためには、$\tilde{M}_{k}$に対応する行を$M_{k+1}$からそのまま除去してもよい。

言い回しが少し難しそうに見えるが、今の我々の具体的な状況では、$\tilde{M}_{1}$の$1$-単体$ab,bc,cd,ad,ac$のうち$cd,bc,ab$のピボットだけが残ったので、これをそのまま$M_{2}$から削除してもよいということである。直感的に考えてみると、これらはすでに$k$次元で使われたので$k+1$では見る必要もない、という程度に受け取ってもかまわない。このように列-階段形$\tilde{M}_{2}$を直接求める過程を省略し、その三つの行を削除してみると、次のように下が切られた$\check{M}_{2}$を得る。

$$ \begin{align*} z_{2} =& ac - bc - ab \\ z_{1} =& ad - bc - cd - ab \end{align*} $$

再び補助定理 (2)によって計算してみよう。

• 一行目のピボットが$t^{1}$で $$ \deg z_{2} = \deg \left( ac - bc - ab \right) = \max \deg \left\{ ac , bc , ab \right\} = 3 $$ なので$(3,3+1)$を得る。
• 二行目のピボットが$t^{3}$で $$ \deg z_{1} = \deg \left( ad - bc - cd - ab \right) = \max \deg \left\{ ad , bc , cd , ab \right\} = 2 $$ なので$(2,2+3)$を得る。

これはアルゴリズムを導出する前に言及した$L_{1}$と正確に一致する。 $$ L_{1} = \left\{ \left( 2,5 \right) , \left( 3,4 \right) \right\} $$

このような過程を複体$K$の次元$\dim K$の分だけ繰り返せば、我々が望んでいたアルゴリズムを得る。行列の左右のサイズは$\partial_{k}$に従い、その成分は$\deg$に従って埋めると考えれば、少し混乱が減るだろう。

■

一方、補助定理 (1)で列演算だけで十分だということは、これまでの導出で見てきたように行列表現に固執する理由が特にないという意味でもある。また補助定理 (3)によって「過去にすでに計算が終わった」部分に対して大胆に行を捨てる式の効率的なプロシージャが含まれているが、そのためにはピボットでない列を「マーキング」する能力などが必要である。結果的に実際のアルゴリズムの疑似コード^{pseudo Code}は、行列をそのまま使うのではなく、もう少し高水準のデータ型、辞書あるいはデータフレームのようなもので説明することになる。これは直接体験してみると非常に戸惑い、難しい。

実装

• ゾモロディアンのアルゴリズム実装: 科学界に従事しているなら誰でも読みやすいJulia言語を通じて、論文の疑似コードをほぼ文学的に移した実装を紹介する。