球状CNNについて

### 更新履歴
-2019年5月31日：初版

論文["Spherical CNNs"](https://arxiv.org/pdf/1801.10130.pdf%EF%BC%89)を紹介ことになったので，ついでに解説を行う．
本稿は筆者の解釈であるので，誤りや追記などがあればコメントいただきたい．

# モチベーション
3次元オブジェクト認識などで頻繁に利用される3D-CNNは，回転運動に対して耐性が弱いことが知られている．もともとの球状CNNのモチベーションは，三次元オブジェクトを球状で入力することでその問題を改善しようというものらしい．

論文中では，**立面$\boldsymbol{R}^3$上の移動は平行移動群$T(3)$で表され，これらは群同型$\boldsymbol{R}^3\cong T(3)$であるのに対して，球面$S^2$はそれ上の移動である回転移動群$SO(3)$とは非同型である**[^1]と述べられている．つまり，回転移動を$SO(3)$上に作用する関数として考える必要がある．

ただし，球状CNNには2つの問題がある．
ひとつめは，球面を構成するピクセルボックスが離散的な並進対称性を有するために，球に対する完全に対照的なボックスが存在しないこと．ふたつめに，$SO(3)$の相互相関[^2]の単純な実装が$O(n^6)$のオーダーであることである．本論文のモチベーションは，一般化フーリエ変換（$\text{GFT}$）を使用してこれら二つの問題に対処したことである．先に概要を示した図を掲載しておく．

&emsp;&emsp;&emsp;&emsp;&emsp;&emsp;&emsp;&emsp;&emsp;&emsp;&emsp;&emsp;&emsp;&emsp;&emsp;<img src="https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/409783/18be6c96-f5a8-c911-af4a-ff73b3d9dc9c.png" width=40%>
結局これを実現したい、だという理解だ．


[^1]: 正確には，$m$次元数空間上の合同写像は行列式1のアフィン変換で表される（運動群）．論文中では，「移動」という単語をそれぞれのモデルが耐性を持って欲しいもの（CNNには平行移動，球状CNNには回転移動）といった意図で使用している．（と思った．間違ってたらゴメン）

[^2]: 実際，CNNの畳み込み層では数学的によく知られる畳み込み演算ではなく，相互相関というそれによく似た操作を行っている．

# 球面上の相互相関と回転群
一般の平面$\boldsymbol{Z}^2$上の相互相関（平面相関）と同様に，$S^2$と$SO(3)$の相互相関（球面相関）を説明する．

>平面相関の畳み込み層に置ける出力特徴マップは，入力特徴マップとフィルタの内積と計算され，$\boldsymbol{x}$だけ平行移動される．

同様にして，球面相関は次のように考えることができる．

>回転$R\in SO(3)$で評価された出力特徴マップは，入力特徴マップとのフィルタとの間の内積と計算され，$R$だけ回転される．

本章では，必要な概念を調べて，球状CNNの正確な数学モデルを提示する．

### *m*次元球面（*m*-dimensional Sphere）
$\boldsymbol{R}^{m+1}$の部分集合$S^m$を

```math
S^m=\{(x_1, x_2, ...,x_{m+1})\mid x_1^2+x_2^2+...+x_n^2=1\}
```
と定義し，**$m$次元球面**と呼ぶ．$S^m$は，$\boldsymbol{R}^{m+1}$の原点からの距離が$1$であるような点$\boldsymbol{x}=(x_1, x_2,...,x_{m+1})$の集合である．ただし，ここで変形を$1$としたのは便宜上である．$S^m$に，$\boldsymbol{R}^{m+1}$の自然な位相$\mathcal{O}^{(m+1)}$から導かれた相対位相を入れ，それ自身ひとつの位相空間と考える．**$S^m$は$m$次元$C^\infty$球多様体である．**

### 球面信号（Spherical Signals）
球面関数とフィルタを以下のように連続関数としてモデル化する．

```math
f:S^2\rightarrow\boldsymbol{R}^K
```
ただし，$K$はチャネル数を表す．

### 直交群（Orthogonal Group）
集合$O(m)$を次の式で定義する．

```math
O(m)=\{A\mid A:m次正方行列，A^\mathsf{T}A=I_m\}
```

$O(m)$に属する行列$A$を**$m$次直行行列**といい，$O(m)$を**$m$次直交群**という．$A$の行ベクトルを$\boldsymbol{a}\_1,\boldsymbol{a}\_2,...,\boldsymbol{a}\_m$と書くと，$A^\mathsf{T}A=I\_m$の条件は，

```math
\|\boldsymbol{a}_i\|^2=1\hspace{15pt}(1\leq i\leq m)\\
\boldsymbol{a}_i\cdot\boldsymbol{a}_j=0\hspace{30pt}(i\neq j)
```

と同値である．したがって，$(\boldsymbol{a}\_1,\boldsymbol{a}\_2,...,\boldsymbol{a}\_m)$は$\boldsymbol{R}^m$の直行$m$枠であり，対応$A\mapsto(\boldsymbol{a}\_1,\boldsymbol{a}\_2,...,\boldsymbol{a}\_m)$により$O(m)$はシュフィーフェル多様体$V_{m,m}$と同一視され，

```math
\text{dim}V_{m,m}=m^2-\frac{m(m+1)}{2}=\frac{m(m-1)}{2}
```

となり，**$O(m)$は$\frac{m(m-1)}{2}$次元コンパクト$C^\infty$級多様体である．**

### 特殊直行群，回転群(Rotation Group)
n次元の回転の集合$SO(n)$は

```math
SO(n)=\{R\mid R\in O(n),\ \text{det}(R)=1\}
```

### 球面信号の回転（Rotation of Spherical Signals）
回転不変性を持たせるために点$x\in S^2$の回転だけでなく，球面相関の回転を定義する．関数$f$に対して，回転演算子$L_Rf$を導入する．

```math
[L_Rf](x):=f(R^{-1}x)
```

論文中では，$L\_{RR'}=L\_RL_\{R'}$を得ると記述されている．<details><summary>証明</summary><div>

```math
\begin{align}
[L_{RR'}f](x)&=f((RR')^{-1}x)\\
&=f(R'^{-1}(R^{-1}x))\\
&=[L_{R'}f](R^{-1}x)\\
&=[L_RL_{R'}f](x)\\
\Rightarrow &L_{RR'}=L_RL_{R'}\hspace{15pt}\square
\end{align}
```
</div></details>

### 内積（Inner Products）
$S^2$上の内積を次のように定義する．

```math
\langle\psi,f\rangle:=\int_{S^2}\sum_k \psi_k(x)f_k(x)dx
```

積分測度$dx$は球上の標準的な回転不変積分測度を表し，これは球面座標で$d\alpha \sin{\beta}d\beta/4\pi$と表現することが出来る．<details><summary>証明</summary><div>

$SO(3)$において$ZYZ$系オイラー角表現を用いる．

```math
R=R(\alpha, \beta, \gamma)=Z(\alpha)Y(\beta)Z(\gamma)
```
ここで，$\alpha\in [0, 2\pi],\ \beta\in[0, \pi],\ \gamma\in [0,2\pi]$である．よってハール測度は

```math
\begin{align}
&dR=\frac{d\alpha}{2\pi}\frac{d\beta\sin{\beta}}{2}\frac{d\gamma}{2\pi}\\
\Rightarrow&\int_{SO(3)}dR=\int_0^{2\pi}\frac{d\alpha}{2\pi}\int_0^\pi\frac{d\beta\sin{\beta}}{2}\int_0^{2\pi}\frac{d\gamma}{2\pi}=1
\end{align}
```

となる．ハール測度は以下の性質を持つため不変測度とも呼ばれる．

```math
\int_{SO(3)}f(R'R)dR=\int_{SO(3)}f(R)dR
```

さて，球体$S^2$について考えてみよう．球体のパラメータ化は次の通り．

```math
x(\alpha, \beta)=Z(\alpha)Y(\beta)n\hspace{15pt}\text{where}\ x\in S^2, n\text{ is North Pole}
```
ここで，Z軸に関する回転群を$H:=SO(2)$とすると，$S^2=SO(3)/SO(2)$であることが確認できる．このとき，$SO(3)$の元を$\bar{x}=R(\alpha,\beta,0)\in SO(3)$とすると，

```math
\bar{x}H=\{R(\alpha,\beta,\gamma)\mid \gamma\in[0,2\pi]\}
```

は$SO(3)$における$H$の剰余類となる．$S^2$と$SO(2)$に関してハール測度を計算すると，

```math
\begin{align}
dx&=\frac{d\alpha}{2\pi}\frac{d\beta\sin{\beta}}{2}\\
dh&=\frac{d\gamma}{2\pi}
\end{align}
```

ここで，$SO(3)$上の$\gamma$不変関数，すなわち$S^2$上の関数$f:S^2\rightarrow\mathbb{C}$を考え，$\bar{f}(\alpha, \beta,\gamma)=f(\alpha, \beta)$を設ける．

```math
\begin{align}
\int_{SO(3)}\bar{f}(R)dR&=\int_0^{2\pi}\frac{d\alpha}{2\pi}\int_0^\pi\frac{d\beta\sin{\beta}}{2}\int_0^{2\pi}\frac{d\gamma}{2\pi}\ \bar{f}(R)\\
&=\int_0^{2\pi}\frac{d\alpha}{2\pi}\int_0^\pi\frac{\sin{\beta}d\beta}{2}\ f(\alpha,\beta)\int_0^{2\pi}d\gamma\\
&=\int_0^{2\pi}\frac{d\alpha}{2\pi}\int_0^\pi\frac{\sin{\beta}d\beta}{2}\ f(\alpha, \beta)\\
&=\int_{S^2}f(x)dx\hspace{15pt}\square
\end{align}
```

これにより，$d\alpha \sin{\beta}d\beta/4\pi$が得られるのがわかる．さらに，$S^2$上の関数を$SO(3)$上の$\gamma$不変関数と見なし，その$SO(3)$フーリエ変換を取ることによって，**$SO(3)$上のフーリエ変換から$S^2$上のフーリエ変換を定義することができる．**
</div></details>測度の不変性は任意の回転$R\in SO(3)$に対して，

```math
\int_{S^2}f(Rx)dx=\int_{S^2}f(x)dx
```
であることが保証される．このことを利用すれば，

```math
\begin{align}
\langle L_R\psi,f\rangle&=\int_{S^2}\sum_k\psi_k(R^{-1}x)\ f_k(x)dx\\
&=\int_{S^2}\sum_k\psi_k(x)\ f_k(Rx)dx=\langle\psi,L_{R^{-1}}f \rangle
\end{align}
```
となり，$L_R$がユニタリ行列であることがわかる．

### 球面相関（Spherical Correlation）
上記のことを踏まえ，2つの球面信号$\psi,\ f$の相関を次のように定義する．

```math
[\psi\star f](R):=\langle L_R\psi,f\rangle=\int_{S^2}\sum_k\psi_k(R^{-1}x)\ f_k(x)dx
```
前述のように，球面相関の出力は$SO(3)$上の関数となる．

### SO(3)信号の回転（Rotation of SO(3) Signals）
球面信号の回転（Rotation of Spherical Signals）で定義した回転演算子$L_Rf$を一般化する．
$f:SO(3)\rightarrow\boldsymbol{R}^k$，$R, Q\in SO(3)$のとき

```math
[L_Rf](Q)=f(R^{-1}Q)
```
これは単なる回転の合成を表現する行列乗算である．

### 回転群相関（Rotation Group Correlation）
回転群$f,\psi$上の二つの信号の相互相関を次のように定義する．
$f,\psi:SO(3)\rightarrow\boldsymbol{R}^K$として，

```math
[\psi\star f](R)=\langle L_R\psi,f \rangle=\int_{SO(3)}\sum_k \psi_k(R^{-1}Q)\ f_k(Q)dQ
```

積分測度$dQ$は$ZYZ$系オイラー角で表される$S0(3)$の不変測度である．

### 同値性（Equivariance）
上記のように，回転群相関を回転演算子$L_Rf$を導入することによって定義した．
ここでは，あるネットワークの層$\Phi$が行う演算$T_R$について同値$\Phi\circ L_R=T_R\circ\Phi$を示す．

```math
[\psi\star [L_Qf]](R)=\langle L_R \psi, L_Qf\rangle=\langle L_{Q^{-1}R}\psi, f\rangle
=[\psi\star f](Q^{-1}R)=[L_Q[\psi\star f]](R)
```
故に，あらゆる種類の演算（畳み込み，相互相関）において等価性が保証される．この導出は回転群相関だけでなく球面相関にも有効である．しかし，これは$f$と$\psi$が連続関数のときに成立するのであって，離散関数のときは厳密には成立しない． 

# 一般化フーリエ変換による高速球面相関
以下の畳み込み定理に高速フーリエ変換を適用することで，畳み込みおよび相互相関を効率的に計算する手法がある．

```math
\widehat{f*\psi}=\hat{f}\cdot\hat{\psi}
```
<details><summary>証明</summary><div>

```math
\begin{align}
\widehat{f*g}&=\int_{-\infty}^\infty dx\ e^{-2\pi ix\xi}\int_{-\infty}^\infty f(x)g(y-x)dy\\
&=\int_{-\infty}^\infty dx\ e^{-2\pi i(x+y)\xi}\int_{-\infty}^\infty f(x)g(y)dy\\
&=\int_{-\infty}^\infty f(x)e^{-2\pi ix\xi}dx\cdot\int_{-\infty}^\infty g(y)e^{-2\pi iy\xi}dy\\
&=\hat{f}\cdot\hat{g}\hspace{15pt}\square
\end{align}
```
</div></details>

球と回転群の関数については，これに似た変換である一般化フーリエ変換（Generalized-FT）を利用した高速一般化フーリエ変換（GFFT）を用いる．参考は[こちら](https://qiita.com/bokuikkun/items/e6f29707f32d303aed3d).

論文中では，概念的に，GFTは[既約ユニタリ表現](https://qiita.com/bokuikkun/items/6dbc34d194e28c83cb26)の行列要素と呼ばれる一連の直交基底関数への関数の線形射影に他ならない．と記述されている．

例えば，円$S^1$と線$\boldsymbol{R}$の場合，複素指数$e^{i\theta}$がそうである．これはユニタリ作用素であると同時に，次の直交条件を満足する．

```math
\int_{-\pi}^{\pi}e^{im\theta}\cdot\overline{e^{in\theta}}d\theta=2\pi\delta_{m,n}
```
これと同様に$SO(3)$の場合は，$l\geq 0, -l\leq m,n\leq l$の次数を持つウィグナーの$D$関数$D_{mn}^l(R)$がある．$S^2$の場合は，$l\geq 0, -l\leq m\leq l$の次数を持つ[球面調和関数](https://ja.wikipedia.org/wiki/%E7%90%83%E9%9D%A2%E8%AA%BF%E5%92%8C%E9%96%A2%E6%95%B0)$Y_m^l(x)$である[^3].

多様体（$S^2$または$SO(3)$）を$X$で表現し，対応する基底関数を$U^l$（ベクトル$Y^l$または行列$D^l$）で表すと，関数$f$のGFTを次のように記述できる．

```math
\hat{f}^l:=\int_X f(x)\overline{U^l(x)}dx
```
この積分は後述のGFFTアルゴリズムにより高速に計算できる．逆$SO(3)$フーリエ変換は次の通り．

```math
f(R)=\sum_{l=0}^b(2l+1)\sum_{m=-l}^l\sum_{n=-l}^l
\hat{f}_{mn}^lD_{mn}^l(R)
```
これは$S^2$についても同様である．最大周波数$b$は帯域幅と呼ばれ，空間グリッドの分解能に関与する．
ウィグナーの$D$関数は次の等式を満足する．

```math
\begin{align}
&D^l(RR')=D^l(R)D^l(R')\\
&D^l(R^{-1})=D^l(R)^*
\end{align}
```
これにより先述の$SO(3)$についての畳み込み定理の成立が容易に導出される．$S^2$の場合の畳み込み定理の導出には，

```math
\begin{align}
&Y(Rx)=D(R)Y(x)\\
&Y_m^l=\left.D_{m0}^l\right|_{S^2}
\end{align}
```
を用いれば先ほどと同様に導ける．
これは，**2つの球面信号の$S^2$相関の$SO(3)$フーリエ変換が，それぞれの信号の$S^2$フーリエ変換のテンソル積をとることにより計算できることを示す．**

![test.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/409783/b1c676e1-7e98-da31-d545-b44945cddb5d.png)

上図はア○パ○マン(球体)を相互相関する概略図である．
スペクトル内の球面相関信号$f$とフィルタ$\psi$はそれぞれ高速フーリエ変換，離散フーリエ変換が施され[^4]， ブロック毎にテンソル積を計算し，入力チャネルの合計を取る．そして最後に逆フーリエ変換$SO(3)-\text{IFFT}$される．

[^3]: 厳密には，$S^2$は群でないため，既約表現は存在しないが群$SO(3)/SO(2)$の商であり，等式$Y\_m^l=D\_{m0}^l|\_{S^2}$が成り立つ．

[^4]: フィルター$\psi$は局所的に定義されているため，$S^2-\text{FFT}$よりも行列乗算である$S^2-\text{DFT}$の方が高速．

### G-FFTとスペクトルG-Convolutionの実装
この節の詳細についても[こちら](https://qiita.com/bokuikkun/items/e6f29707f32d303aed3d "こちら")を参考にすると良い．
$SO(3)-\text{FFT}$の入力は，$SO(3)$上の空間信号$f$で離散グリッド上でサンプリングされ，軸が$ZYZ$系オイラー角表現$\alpha,\beta,\gamma$に対応した3D配列に格納される．（$c.f.$ 上図）

$2B\times 2B\times 2B$グリッド$\alpha_{j_1}=\cfrac{2\pi j_1}{2B},\ \beta_k=\cfrac{\pi(2k+1)}{4B},\ \gamma_{j_2}=\cfrac{2\pi j_2}{2B}$上でサンプリングされた，関数$f\in SO(3)$の帯域幅$B$での離散$SO(3)$フーリエ変換（$\text{DSOFT}$）は，

```math
\hat{f}_{MM'}^l=\frac{\pi}{(2B)^2}
\sum_{k=0}^{2B-1}w_B(k)\tilde{d}_{MM'}^l(\beta_k)
\sum_{j_2=0}^{2B-1}e^{iM'\gamma_{j_2}}
\sum_{j_1=0}^{2B-1}e^{iM'\alpha_{j_1}}\ f(\alpha_{j_1},\beta_{k},\gamma_{j_2})
```
ただし，$l=0,1,...,B-1,\ -l\leq M,M'\leq l,\ w_B(j)$は帯域幅$B$の変換に用いる直交重み関数である．
この式からわかるように，$SO(3)-\text{FFT}$ははじめに，$\alpha$軸と$\gamma$軸に対して標準2D並進$\text{FFT}$を計算する．この結果は式中では$M,M'$に寄与する．次に，$\text{FFT}$配列の$\beta$軸と，$L^2$正規化されたウィグナーの$d$関数$\tilde{d}_{mn}^l(\beta)$の線形収縮である．

# まとめ
実験の章に関しては，記事にする必要がないので理論部分だけで書き終えることにする．