[グラフニューラルネットワーク] グラフ信号処理における空間領域での畳み込みと周波数領域での要素積について

2022年08月03日

グラフ信号処理における空間領域の畳み込みと周波数領域での要素積について

グラフ信号処理において、空間領域の畳み込みと周波数領域での要素積は本質的に同じことを意味する。逆に、空間領域での要素積（乗算）は、周波数領域での畳み込みに対応する。

この定理を式で表すと以下の通り：

空間領域の畳み込み → 周波数領域での要素積

$$ \mathcal{F}(f * g) = \mathcal{F}(f) \cdot \mathcal{F}(g) $$

$ f * g $ は空間領域での畳み込み
$ \mathcal{F}(f) $ と $ \mathcal{F}(g) $ はそれぞれのフーリエ変換
$\cdot$ は周波数領域での要素積（成分ごとの乗算）

周波数領域での要素積 → 空間領域での畳み込み

$$ \mathcal{F}^{-1}(\mathcal{F}(f) \cdot \mathcal{F}(g)) = f * g $$

直感的な解釈

空間領域での畳み込みは、信号（画像や音など）を滑らかにしたり、特定の特徴を抽出する操作
周波数領域での要素積は、信号の周波数成分に対してフィルタリングを行う操作

この記事では、グラフニューラルネットワークの基本的な概念を理解するために、グラフ信号処理における空間領域の畳み込みと周波数領域での要素積について理論とPythonによる実装について説明する。

ソースコード

github

jupyter notebook形式のファイルはこちら

google colaboratory

google colaboratory で実行する場合はこちら

実行環境

OSはmacOSである。LinuxやUnixのコマンドとはオプションが異なりますので注意していただきたい。

!sw_vers

ProductName:		macOS
ProductVersion:		15.2
BuildVersion:		24C101

!python -V

Python 3.9.17

基本的なライブラリをインポートし watermark を利用してそのバージョンを確認しておきます。ついでに乱数のseedの設定をします。

%matplotlib inline
%config InlineBackend.figure_format = 'svg'

import random

import scipy
import numpy as np

import matplotlib
import matplotlib.pyplot as plt

seed = 123
random_state = 123

random.seed(seed)
np.random.seed(seed)


from watermark import watermark

print(watermark(python=True, watermark=True, iversions=True, globals_=globals()))

Python implementation: CPython
Python version       : 3.9.17
IPython version      : 8.17.2

matplotlib: 3.8.1
scipy     : 1.11.2
numpy     : 1.25.2

Watermark: 2.4.3

1. グラフ信号処理の基礎

1.1 グラフと信号の定義

グラフ $G$ は、頂点集合 $V$ と辺集合 $E$ からなる。頂点数を $n = |V|$ とすると、辺の接続関係は隣接行列 $\mathbf{A}$ によって表現される。ここで、$\mathbf{A}$ は $n \times n$ の行列となる。辺が存在する場合は1、存在しない場合は0、などの定義を用いることが多い。 (重み付きグラフの場合は非0の重みを記入する。)

各頂点 $v_i$ 上に実数またはベクトル値の信号 $x_i$ が定義されるとき、これを「グラフ信号」と呼ぶ。頂点ごとにスカラーの信号を考える場合は、$\mathbf{x}$ を長さ $n$ のベクトルとみなせる。すなわち、 $$ \mathbf{x} = \begin{bmatrix} x_1 \\ x_2 \\ \vdots \\ x_n \end{bmatrix} $$ という形で表現できる。

このようなグラフ信号に対し、畳み込みの操作を定義したい。しかし、通常のCNNで行う畳み込みとは異なり、グラフ構造をもつ場合の畳み込みは容易ではない。そこで、空間領域と周波数領域の二つのアプローチが生まれている。

1.2 空間領域の畳み込み: 近傍集約

空間領域(Spatial Domain)の畳み込みは、各頂点の近傍頂点との情報を直接やりとりする。例えば、GNNの代表例であるGCN (Graph Convolutional Network) では、 $$ \mathbf{X}^{(l+1)} = \sigma\left( \tilde{\mathbf{D}}^{-\frac12} \tilde{\mathbf{A}} \tilde{\mathbf{D}}^{-\frac12} \mathbf{X}^{(l)} \mathbf{W}^{(l)} \right) $$ のような更新式が使われる。ここで、$\tilde{\mathbf{A}} = \mathbf{A} + \mathbf{I}$ は自己ループを追加した隣接行列、$\tilde{\mathbf{D}}$ は $\tilde{\mathbf{A}}$ の対角要素を用いた次数行列である。 $\mathbf{X}^{(l)}$ は第 $l$ 層目の頂点特徴ベクトルの集合、$\mathbf{W}^{(l)}$ は学習パラメータ行列、$\sigma$ は活性化関数である。

要するに、この手法では隣接行列を使って近傍頂点との平均や加重和を取る形となる。これを「空間領域の畳み込み」と呼ぶことが多いが、実際には近傍集約(Neighbor Aggregation)に近い。頂点 $v_i$ の特徴ベクトルは、その周囲の頂点集合からの情報を統合する。

この方法は直感的であり、実装もしやすい。多くの実用的なGNNライブラリでも、この空間領域の手法が取り入れられている。企業でSNSデータを解析する際などは、ユーザの近傍ユーザからの特徴を集約する手法となる。

1.3 周波数領域における要素積: グラフフーリエ変換

一方、周波数領域(Frequency Domain)での要素積を考える際は、グラフ上のフーリエ変換を用いる。通常の画像や音声でも、空間領域での畳み込みは周波数領域での要素積に対応する。グラフでも同様に、グラフラプラシアン $\mathbf{L}$ の固有ベクトルを用いてフーリエ変換を行う考え方が生まれた。ラプラシアンを $$ \mathbf{L} = \mathbf{D} - \mathbf{A} $$ と定義した場合、ここで $\mathbf{D}$ は次数行列である。 $\mathbf{L}$ は対称行列となる。従って、その固有ベクトルは直交基底を形成する。これを用いて、グラフ信号 $\mathbf{x}$ は以下のように変換できる。 $$ \hat{\mathbf{x}} = \mathbf{U}^T \mathbf{x} $$ ただし、$\mathbf{U}$ は $\mathbf{L}$ の固有ベクトルを列ベクトルに並べた行列である。

この $\hat{\mathbf{x}}$ は、周波数領域における表現である。通常のフーリエ変換において、畳み込みは周波数領域での積になる。よって、グラフ信号同士の畳み込みは、周波数領域での要素積として実装可能である。ただし、グラフごとに固有ベクトルが異なるため、大規模なグラフでは固有分解が多くの時間とメモリが必要となる。これが周波数領域で考える際のデメリットの一つとなる。

2. 空間領域と周波数領域の対応

2.1 畳み込みの基本性質

通常の1次元または2次元の離散フーリエ変換(DFT)では、次の性質が成立する。

時間(または空間)領域での畳み込み $$ \mathbf{x} * \mathbf{y} $$
周波数領域での要素積 $$ \hat{\mathbf{x}} \odot \hat{\mathbf{y}} $$ ここで、$\odot$ は要素積である。グラフ上でも、フーリエ変換をグラフラプラシアンの固有ベクトルを用いて定義すれば、同様の対応が成り立つ。

2.2 GCNにおける近似手法

周波数領域の畳み込みを直接行うには、固有分解を計算する必要がある。しかし、頂点数が膨大な場合、この計算は非常にコストが高い。そこで、多くの手法ではラプラシアンを多項式近似する。代表的な例としては、Chebyshev多項式を用いた近似が知られている。 Chebyshev多項式を用いれば、高次の畳み込み演算を繰り返すことができるため、周波数領域の計算を回避することができる。

その結果、実装上は、空間領域での隣接行列操作に近い形で計算できる。また、GCNの単純な更新式は、Chebyshev多項式の次数を1程度に落としたものとみなされる。このように、空間領域の手法も、もともとは周波数領域での畳み込みを近似した解釈が可能である。

2.3 メリットとデメリットのまとめ

空間領域手法
- メリット: 直感的に理解しやすい。
- メリット: 実装が容易で、大規模グラフにも適用しやすい。
周波数領域手法
- メリット: 理論的には畳み込みの定義が厳密。
- デメリット: 固有分解や変換が大規模では時間やメモリを必要とする。
- デメリット: 近似や省略なしでは実用に難がある。

3. 具体例: Pythonによる実装

ここでは簡単なサンプルグラフを用いて、空間領域と周波数領域の実装イメージを示す。プログラムはPythonを用いる。

3.1 サンプルグラフの作成

import numpy as np

# 頂点数
n = 5

# 隣接行列の作成
A = np.array([[0, 1, 0, 0, 1], [1, 0, 1, 0, 0], [0, 1, 0, 1, 0], [0, 0, 1, 0, 1], [1, 0, 0, 1, 0]], dtype=np.float32)
A

array([[0., 1., 0., 0., 1.],
       [1., 0., 1., 0., 0.],
       [0., 1., 0., 1., 0.],
       [0., 0., 1., 0., 1.],
       [1., 0., 0., 1., 0.]], dtype=float32)

# グラフ信号 (各頂点に1次元の値)
x = np.array([1, 2, 3, 4, 5], dtype=np.float32)
x

array([1., 2., 3., 4., 5.], dtype=float32)

ここでは頂点数が5の単純なグラフを作成した。頂点同士が環状に接続されているイメージである。それぞれの頂点にスカラー値を持たせている。

3.2 空間領域での近傍集約

空間領域の手法では、隣接行列を用いて近傍の情報を集約する。特にGCNライクな平均化を行う簡単な例を示す。

# 次数行列 D の作成
D = np.diag(np.sum(A, axis=1))

# 正規化のために D の逆数平方根を計算
D_inv_sqrt = np.diag(1.0 / np.sqrt(np.sum(A, axis=1)))

# 自己ループを加える
A_hat = A + np.eye(n)
D_hat = np.diag(np.sum(A_hat, axis=1))
D_hat_inv_sqrt = np.diag(1.0 / np.sqrt(np.sum(A_hat, axis=1)))

# 空間領域における1層分の演算 (GCNの単純化)
W = np.ones((1, 1), dtype=np.float32)
Z = D_hat_inv_sqrt @ A_hat @ D_hat_inv_sqrt @ x.reshape(-1, 1) @ W
Z = Z.flatten()

print("Spatial Conv Result:", Z)

Spatial Conv Result: [2.66666667 2.         3.         4.         3.33333333]

ここでは単純化のため、重み $\mathbf{W}$ を全1とした。実際のGCN実装では、より高次元の行列となり、学習される。結果として得られる $Z_np$ は、空間領域での近傍情報集約後の信号である。

3.3 周波数領域での要素積

次に、周波数領域で実装する簡単な例を示す。実際にはフーリエ変換の計算で固有分解を要する。小規模グラフなので計算は可能だが、大規模になると困難となる。

# グラフラプラシアン
L = D - A

# 固有分解
Lambda, U = np.linalg.eig(L)
# Lambda は固有値のベクトル
# U の列が対応する固有ベクトル

# フーリエ変換
x_hat = U.T @ x

# ここでは適当に別のフィルタ h_hat を用意し、要素積を行う
# 例えばローパスフィルタを想定して、小さい固有値の成分だけ残す簡単な例
h_hat = np.exp(-0.5 * Lambda)  # e^(-0.5 * λ) という形のフィルタ

# 周波数領域での要素積
y_hat = x_hat * h_hat

# 逆フーリエ変換 (空間領域に戻す)
y = U @ y_hat

print("Frequency Domain Conv Result:", y.real)

Frequency Domain Conv Result: [1.9138175 2.111436  3.3199701 4.0403714 3.6144037]

この例では、固有値 $\lambda_i$ に基づいて減衰関数 $e^{-0.5 \lambda_i}$ を適用した。これによりローパスフィルタ的な効果を与えている。大規模グラフで同様の操作を行うには、固有分解がネックとなる。実際の企業利用においては、この周波数領域の計算コストを多項式近似などで削減する手法が多い。

4. 数学的詳細と計算例

4.1 畳み込みの定義域と写像

グラフ上の畳み込みは、信号空間から信号空間への写像として考えられる。すなわち、 $$ \text{Conv}: \mathbb{R}^n \times \mathbb{R}^n \rightarrow \mathbb{R}^n $$ という写像を定義するイメージである。ただし、ここでの $\mathbb{R}^n$ は頂点数が $n$ の場合のスカラー信号を表す。空間領域であれ周波数領域であれ、結果は再び $\mathbb{R}^n$ の信号となる。

4.2 フィルタ作用素の例

例えば、周波数領域でのフィルタ $g(\mathbf{L})$ を多項式展開するとき、 $$ g(\mathbf{L}) = \sum_{k=0}^{K} \alpha_k \mathbf{L}^k $$ となることがある。この $g(\mathbf{L})$ を信号ベクトル $\mathbf{x}$ に適用すると、 $$ g(\mathbf{L})\mathbf{x} = \sum_{k=0}^{K} \alpha_k \mathbf{L}^k \mathbf{x} $$ という形で空間領域の演算として表現できる。このようにすれば、大きな行列の固有分解を直接しなくても、$\mathbf{L}$ の累乗を計算すればよい。ただし、$\mathbf{L}$ を何度も乗算するため、計算量には注意が必要である。適度な近似階数 $K$ を選ぶことが肝要となる。

4.3 簡単な数値例

たとえば、頂点数 $n = 4$ のグラフを想定し、各辺が均一の重みをもつとする。次数がすべて2のような4頂点の環状グラフ(Cycle)を考える。ラプラシアン $\mathbf{L}$ は以下のようになる(自己ループなし)： $$ \mathbf{L} = \begin{bmatrix} 2 & -1 & 0 & -1 \\ -1 & 2 & -1 & 0 \\ 0 & -1 & 2 & -1 \\ -1 & 0 & -1 & 2 \end{bmatrix} $$

この $\mathbf{L}$ の固有値は $0, 2, 4, 2$ (順序は異なる場合がある)。対応する固有ベクトルを行列 $\mathbf{U}$ として集めるとする。ここで、ある信号 $\mathbf{x}$ に対してローパスフィルタ $g(\lambda)=e^{-0.5\lambda}$ を適用するとき、 $$ \mathbf{y} = g(\mathbf{L}) \mathbf{x} $$ となる。実際に数値を入れて計算すれば、周波数領域での減衰が確認できる。空間領域で同様の結果を得るには、多項式近似などで $\mathbf{L}$ の累乗と係数を使って実装できる。

5. コード実行例と数値実験

次に、もう少し踏み込んだ数値実験例を示す。具体的な計算例を示すことで、空間領域と周波数領域の違いを確認する。また、簡単な多項式近似によるフィルタ適用の例も示す。

5.1 ラプラシアンの多項式近似

例として、以下のようなPythonコードを挙げる。 Chebyshev近似を本格的に組むのは少し複雑なので、単純に $\mathbf{L}^k$ での多項式近似をする例を示す。

def polynomial_filter(L, x, coeff_list):
    """
    L^k の線形和でフィルタを適用する簡易関数
    coeff_list は [α0, α1, α2, ...] の係数リスト
    """
    y = np.zeros_like(x)
    L_power = np.eye(L.shape[0])
    for k_int, alpha_float in enumerate(coeff_list):
        if k_int > 0:
            L_power = L_power @ L
        y += alpha_float * (L_power @ x)
    return y


# サンプル: 多項式 g(L) = I - 0.1 * L  (K=1 の例)
coeff_list = [1.0, -0.1]
y_poly = polynomial_filter(L, x, coeff_list)
print("Polynomial Filter Result:", y_poly)

Polynomial Filter Result: [1.5 2.  3.  4.  4.5]

ここで、$g(\mathbf{L}) = \mathbf{I} - 0.1\mathbf{L}$ という一次多項式フィルタを適用した。これは周波数領域で考えれば、$g(\lambda) = 1 - 0.1\lambda$ というフィルタになる。実際に固有分解して比較してみれば、周波数領域での要素積と同等の結果が得られるはずである。

4.2 結果の比較と考察

空間領域GCN風の近傍集約: 頂点の局所構造(隣接行列)に対して正規化を行いながら情報を集約する。メモリアクセスや行列乗算が中心であり、大規模化しやすい。
周波数領域の要素積: 固有分解を用いて厳密な畳み込み演算を行う。小規模グラフでは高い表現力を発揮する。ただし、大規模グラフでは計算負荷がネックとなる。
多項式近似: 周波数領域の理論を空間領域で実装可能にするブリッジ的手法。フィルタ関数を多項式で近似し、$\mathbf{L}$ の累乗と係数の線形和で演算する。固有分解を回避できるが、近似精度は多項式の次数やグラフの構造に依存する。

5. 企業利用の観点と応用例

企業や研究所で利用されるケースをいくつか挙げる。たとえば、ユーザとアイテムの関係をグラフとして表現し、推薦システムを構築する。ユーザ同士の類似度や、アイテム同士の類似度を辺として定義する場合が多い。そこにGNNを適用すれば、ユーザの近傍ユーザが高評価しているアイテムをレコメンドする仕組みを自然に構築できる。このとき、空間領域の畳み込みによってユーザ-アイテム近傍関係を反映できる。一方で、もし大量の頂点があっても、周波数領域の固有分解を直接やるのは難しい。よって、実際の業務システムでは多くが空間領域ベースの近傍集約手法を採用している。また、Chebyshev多項式などによるスペクトル近似を組み合わせることで、さらに表現力を高めるケースもある。

グラフ理論を応用したSNSの分析でも、空間領域ベースのGNNが用いられる。友人関係(辺)を直接たどるような解釈で学習が進むため、説明もしやすい。広告配信などの分野でも、グラフ構造を利用したターゲティングが行われる場合、空間領域の畳み込みが用いられることが多い。

6. 実用上の注意点と展望

6.1 大規模グラフへの拡張

空間領域手法は、大規模グラフに対してスケールしやすいという特徴がある。周波数領域手法では、固有分解に $O(n^3)$ 程度の計算量が必要となる。よって、頂点数数百万レベルの実務環境では、スペクトル法の直接的な実装は困難である。そのため、多項式近似やサンプリングを組み合わせる手段が研究されている。

6.2 異種グラフや動的グラフへの対応

実際の産業応用では、頂点や辺の種類が複数にわたる異種グラフや、時間とともに辺が変化する動的グラフがしばしば登場する。こうしたグラフに対して畳み込みを定義する場合、空間領域・周波数領域いずれの手法でも工夫が必要となる。辺ごとに重みだけでなく、タイプ情報を扱う必要があり、近傍集約のルールが複雑になる。周波数領域手法では、ラプラシアン自体が複雑化するうえ、固有分解への負荷がさらに増大する。

6.3 今後の研究動向

GNN全体としては、Transformer構造との融合や、グラフ注意ネットワーク(Graph Attention Network)の拡張などが注目を集めている。また、大規模分散環境でのGNN学習手法や、高速サンプリング技術なども活発に研究されている。

一方、周波数領域に基づく厳密な解釈をどのように残しつつ、計算負荷を低減するか、という点も議論が続いている。

Chebyshev多項式以外の多項式近似や、ランダム特徴マッピングによる近似が検討されるケースもある。これらの研究は、企業の大規模推薦システムやSNS分析などに直結する可能性がある。

結論

この記事では、グラフ理論における空間領域の畳み込みと周波数領域での要素積の対応関係を概説した。空間領域の手法は、大規模グラフでの実装が容易であり、企業での実利用例が多い。周波数領域の手法は、厳密な畳み込み定義が可能だが、大規模化が困難である。両者を多項式近似によって橋渡しすることが、現実的な解法として注目を集めている。

参考文献:

Kipf, T. N., & Welling, M. (2017). “Semi-Supervised Classification with Graph Convolutional Networks.” ICLR.
Bruna, J., Zaremba, W., Szlam, A., & LeCun, Y. (2014). “Spectral Networks and Locally Connected Networks on Graphs.” ICLR.