画像の上で動く畳み込みニューラルネットワーク（CNN）が成功した理由は、フィルタが「隣り合うピクセル」だけを見て局所的なパターン（エッジや角）を捉えるからでした。では、ソーシャルネットワークや分子構造、道路網のような「グリッドではないデータ」の上で、同じように畳み込みを定義するにはどうすればよいのでしょうか。グラフには「上下左右」という固定の方向がなく、各ノードの隣の数（次数）もバラバラです。CNNのように $3\times 3$ のフィルタを滑らせる、という素朴な方法はそのままでは使えません。
この問いに対する最も美しい答えのひとつが、スペクトルグラフ畳み込みです。アイデアは「畳み込みはフーリエ領域では掛け算になる」という信号処理の基本定理を、グラフの上に持ち込むことです。グラフ版のフーリエ変換をグラフラプラシアンの固有分解で定義し、その変換領域でフィルタを掛けることで、グラフ上の畳み込みを厳密に定義できます。この理論を実用レベルまで落とし込んだのが、チェビシェフ多項式でフィルタを局所化したChebNetであり、それをさらに極限まで簡約したものが、今日のグラフニューラルネットワークの出発点となったGCN（Graph Convolutional Network）層です。
スペクトルグラフ畳み込みを理解すると、次のような応用が見えてきます。ひとつは化学・創薬で、分子をグラフとみなして物性や毒性を予測するタスク。もうひとつは交通・電力網の予測で、道路網や送電網の各ノードの状態を、ネットワーク構造を活かして予測するタスクです。さらにソーシャルネットワークのノード分類やレコメンデーションにも直結します。「隣接行列をどう掛けるとなぜグラフ上の畳み込みになるのか」という疑問に、スペクトルの言葉できちんと答えられるようになるのが本記事のゴールです。
本記事の内容
- グラフラプラシアンの固有分解とグラフフーリエ変換の直感と定義
- スペクトル畳み込み $g_\theta \star \bm{x} = \bm{U} g_\theta(\bm{\Lambda}) \bm{U}^\top \bm{x}$ の導出
- チェビシェフ多項式による $K$ 次局所近似（ChebNet）の導出
- $K=1$・$\lambda_{\max}\approx 2$ から GCN層 $\bm{H} = \sigma(\tilde{\bm{D}}^{-1/2}\tilde{\bm{A}}\tilde{\bm{D}}^{-1/2}\bm{X}\bm{W})$ への簡約を一行ずつ
- Pythonでの固有分解・局所性の可視化・空手クラブグラフでのフィルタ応答比較
前提知識
この記事を読む前に、以下の記事を読んでおくと理解が深まります。
- グラフラプラシアンとは？定義・固有値・スペクトルの意味を解説 — 本記事の出発点。ラプラシアンの定義と固有値の意味
- グラフ畳み込みネットワーク（GCN）をわかりやすく解説して実装する — 本記事が導出する GCN 層の実装と応用
- Graph Attention Network（GAT）の理論と実装 — スペクトルではなく空間的注意でグラフ畳み込みを行う対比的手法
線形代数の固有値分解と、フーリエ変換の「畳み込み定理」を知っていると読みやすいですが、本文中でも必要な部分は丁寧に補います。
グラフ信号とは — ノードに乗った値を「信号」とみなす
スペクトルの議論に入る前に、まず「グラフ上の信号」とは何かをイメージしておきましょう。普通の信号処理では、時刻 $t$ ごとに値が並んだ時系列 $x(t)$ や、画素位置ごとに値が並んだ画像を「信号」と呼びました。グラフ信号も発想は同じで、グラフの各ノードに 1 つの数値が乗っているものを信号と呼びます。
たとえばソーシャルネットワークなら、各ユーザー（ノード）の年齢を並べたベクトルが 1 つのグラフ信号です。センサーネットワークなら、各観測点の温度を並べたベクトルが信号になります。ノードが $N$ 個あれば、グラフ信号は $N$ 次元ベクトル $\bm{x} \in \mathbb{R}^N$ で、$\bm{x}$ の第 $i$ 成分 $x_i$ がノード $i$ に乗った値です。実際の機械学習では各ノードが複数の特徴量を持つので、信号は $\bm{X} \in \mathbb{R}^{N \times F}$（$F$ 個のチャネル）という行列になりますが、当面は理解のために 1 チャネル $\bm{x} \in \mathbb{R}^N$ で考えます。
普通の信号処理で「滑らかな信号」とは、時間方向に急激に変化しない信号でした。グラフ信号では、辺でつながった隣接ノード同士の値が近いとき「滑らか」と考えます。逆に、隣接ノードの値が大きく食い違っている信号は「高周波」です。この「滑らかさ＝低周波」という直感を数学的に測る道具こそが、次に見るグラフラプラシアンです。
ここまでで、グラフ信号とその「滑らかさ」という言葉を共有できました。次に、その滑らかさを定量化し、さらにグラフ版のフーリエ変換を生み出すグラフラプラシアンを見ていきます。
グラフラプラシアンと正規化ラプラシアン
定義の復習
無向グラフ $G = (V, E)$ がノード数 $N = |V|$ を持つとします。隣接行列 $\bm{A} \in \mathbb{R}^{N\times N}$ は、ノード $i,j$ が辺でつながっていれば $A_{ij}=1$（重み付きなら重み $w_{ij}$）、そうでなければ $0$ とします。次数行列 $\bm{D}$ は対角行列で、対角成分 $D_{ii} = \sum_j A_{ij}$ はノード $i$ につながる辺の重みの総和（次数）です。
（組合せ）グラフラプラシアンは
$$ \begin{equation} \bm{L} = \bm{D} – \bm{A} \end{equation} $$
で定義されます。この $\bm{L}$ がなぜ「滑らかさ」を測るのかは、二次形式 $\bm{x}^\top \bm{L} \bm{x}$ を計算するとはっきりします。
$$ \begin{align} \bm{x}^\top \bm{L} \bm{x} &= \bm{x}^\top \bm{D} \bm{x} – \bm{x}^\top \bm{A} \bm{x} \\ &= \sum_i D_{ii} x_i^2 – \sum_{i,j} A_{ij} x_i x_j \end{align} $$
ここで $D_{ii} = \sum_j A_{ij}$ を第1項に代入すると、$\sum_i D_{ii}x_i^2 = \sum_{i,j}A_{ij}x_i^2$ と書き直せます。さらに対称性 $A_{ij}=A_{ji}$ を使って和を辺ごとにまとめると、
$$ \begin{equation} \bm{x}^\top \bm{L} \bm{x} = \frac{1}{2}\sum_{i,j} A_{ij}\,(x_i – x_j)^2 \end{equation} $$
という非常に意味の明確な形になります。これは「辺でつながった隣接ノード同士の値の差の二乗を、重み付きで足し合わせたもの」です。隣接ノードの値が揃っている滑らかな信号ほど $\bm{x}^\top \bm{L}\bm{x}$ は小さく、食い違っている信号ほど大きくなります。つまり $\bm{x}^\top \bm{L}\bm{x}$ はグラフ信号の「非滑らかさ（変動の激しさ）」を測る量であり、これがのちに「グラフ周波数」へ直結します。
正規化グラフラプラシアン
スペクトルグラフ畳み込みでは、組合せラプラシアン $\bm{L}=\bm{D}-\bm{A}$ ではなく、対称正規化ラプラシアン
$$ \begin{equation} \bm{L}_{\mathrm{sym}} = \bm{I} – \bm{D}^{-1/2}\bm{A}\bm{D}^{-1/2} \end{equation} $$
を使います。ここで $\bm{D}^{-1/2}$ は対角成分が $1/\sqrt{D_{ii}}$ の対角行列です。なぜわざわざ正規化するのでしょうか。理由は2つあります。第一に、組合せラプラシアンは次数の大きいノード（ハブ）の影響が過大になりがちで、$\bm{x}^\top\bm{L}\bm{x}$ の値がノードの次数に強く左右されてしまいます。正規化は各ノードを次数で割ることでこの偏りをならします。第二に、そして決定的に重要なのは、$\bm{L}_{\mathrm{sym}}$ の固有値が必ず $[0, 2]$ の区間に収まるという性質です。この有界性が、後でチェビシェフ多項式による近似を可能にします。
正規化された二次形式は
$$ \begin{equation} \bm{x}^\top \bm{L}_{\mathrm{sym}} \bm{x} = \frac{1}{2}\sum_{i,j} A_{ij}\left(\frac{x_i}{\sqrt{D_{ii}}} – \frac{x_j}{\sqrt{D_{jj}}}\right)^2 \ge 0 \end{equation} $$
となり、これが非負であることから $\bm{L}_{\mathrm{sym}}$ は半正定値、つまり固有値はすべて $0$ 以上です。上限が $2$ である証明は後で固有値の文脈で改めて触れます。
ここまでで、滑らかさを測る正規化ラプラシアンを手に入れました。次は、この行列を固有分解すると「グラフ版のフーリエ基底」が現れる、という核心に進みます。
グラフフーリエ変換 — ラプラシアンの固有分解
なぜ固有分解がフーリエ変換になるのか
普通のフーリエ変換は、信号を「いろいろな周波数の正弦波の重ね合わせ」に分解する操作でした。ではその「正弦波」はどこから来るのでしょうか。実は、1次元の連続信号におけるラプラス作用素 $-\frac{d^2}{dt^2}$ の固有関数がまさに正弦波・複素指数関数 $e^{i\omega t}$ であり、その固有値が $\omega^2$（周波数の二乗）なのです。「フーリエ変換とはラプラシアンの固有関数で信号を展開すること」と言い換えられます。
この見方をグラフに移植します。グラフ上でラプラス作用素に当たるのがグラフラプラシアン $\bm{L}_{\mathrm{sym}}$ ですから、その固有ベクトルがグラフ版の正弦波（フーリエ基底）であり、対応する固有値がグラフ版の周波数だと考えるのです。これがスペクトルグラフ理論の出発点です。
固有分解の定義
$\bm{L}_{\mathrm{sym}}$ は対称半正定値行列なので、スペクトル定理により直交行列で対角化できます。
$$ \begin{equation} \bm{L}_{\mathrm{sym}} = \bm{U} \bm{\Lambda} \bm{U}^\top \end{equation} $$
ここで $\bm{U} = [\bm{u}_1, \bm{u}_2, \dots, \bm{u}_N]$ は固有ベクトルを列に並べた直交行列（$\bm{U}^\top \bm{U} = \bm{I}$）、$\bm{\Lambda} = \mathrm{diag}(\lambda_1, \dots, \lambda_N)$ は固有値を並べた対角行列です。固有値は
$$ 0 = \lambda_1 \le \lambda_2 \le \cdots \le \lambda_N \le 2 $$
の順に並べるのが慣習です。最小固有値が $0$ になるのは、定数信号 $\bm{x} \propto \bm{D}^{1/2}\bm{1}$ が変動ゼロ（最も滑らか）に対応するからで、これがグラフ版の「直流成分（周波数 $0$）」です。最大固有値に対応する固有ベクトルは、隣接ノードで符号が反転するような最も激しく振動する信号で、グラフ版の「最高周波数」に対応します。固有値 $\lambda_k$ の大きさが、その固有ベクトル $\bm{u}_k$ の「グラフ周波数」を表しているわけです。
グラフフーリエ変換と逆変換
固有ベクトル $\bm{u}_k$ をフーリエ基底とみなすと、グラフ信号 $\bm{x}$ のグラフフーリエ変換 $\hat{\bm{x}}$ は、信号を各基底に射影した係数として
$$ \begin{equation} \hat{\bm{x}} = \bm{U}^\top \bm{x}, \qquad \hat{x}_k = \bm{u}_k^\top \bm{x} \end{equation} $$
で定義されます。$\hat{x}_k$ は「信号 $\bm{x}$ に含まれる周波数 $\lambda_k$ 成分の強さ」です。逆変換は $\bm{U}$ が直交行列であることから
$$ \begin{equation} \bm{x} = \bm{U} \hat{\bm{x}} = \sum_{k=1}^{N} \hat{x}_k\, \bm{u}_k \end{equation} $$
となります。普通のフーリエ変換の「正弦波の重ね合わせ」に完全に対応する形です。
ここまでで、グラフ信号を周波数領域に移す道具（グラフフーリエ変換）が揃いました。次はいよいよ、この変換を使ってグラフ上の畳み込みを定義します。
スペクトルグラフ畳み込みの定義と導出
畳み込み定理を出発点にする
普通の信号処理で最も基本的な定理のひとつが畳み込み定理です。「時間領域での畳み込みは、周波数領域での要素ごとの掛け算に等しい」というものでした。式で書けば
$$ \mathcal{F}\{f \star g\} = \mathcal{F}\{f\} \cdot \mathcal{F}\{g\} $$
です。グラフ上では「畳み込み」という操作を一から定義する必要がありますが、それを畳み込み定理が成り立つように逆算して定義するのが賢いやり方です。つまり、グラフ信号 $\bm{x}$ とフィルタ $\bm{g}$ の畳み込みを「グラフフーリエ領域での要素ごとの積」として定義します。
両者をグラフフーリエ変換し（$\hat{\bm{x}} = \bm{U}^\top\bm{x}$、$\hat{\bm{g}} = \bm{U}^\top\bm{g}$）、要素ごとに掛けてから逆変換すると、
$$ \begin{equation} \bm{g} \star \bm{x} = \bm{U}\Big( (\bm{U}^\top \bm{g}) \odot (\bm{U}^\top \bm{x}) \Big) \end{equation} $$
となります（$\odot$ は要素ごとの積、アダマール積）。ここで重要な発想の転換をします。フィルタ $\bm{g}$ を空間領域（ノード上の値）で設計するのではなく、最初からフーリエ領域での応答 $\hat{\bm{g}}$ を直接パラメータとして設計するのです。各周波数 $\lambda_k$ に対するフィルタの利得を $g_\theta(\lambda_k)$ と書き、これを対角行列にまとめます。
$$ g_\theta(\bm{\Lambda}) = \mathrm{diag}\big(g_\theta(\lambda_1), \dots, g_\theta(\lambda_N)\big) $$
すると、$\hat{\bm{g}}\odot\hat{\bm{x}}$ は $g_\theta(\bm{\Lambda})\,\hat{\bm{x}}$ と行列の形で書けます。これを上の式に代入すると、スペクトルグラフ畳み込みの定義式
$$ \begin{equation} g_\theta \star \bm{x} = \bm{U}\, g_\theta(\bm{\Lambda})\, \bm{U}^\top \bm{x} \end{equation} $$
が得られます。これがスペクトルグラフ畳み込みの中核となる式です。
式の意味を読み解く
この式は3つの操作の合成として読めます。右から順に、$\bm{U}^\top\bm{x}$ で信号を周波数領域へ移し、$g_\theta(\bm{\Lambda})$ で各周波数成分にフィルタの利得を掛け、$\bm{U}$ で空間領域に戻す、という流れです。普通の信号処理の「FFT → フィルタ → 逆FFT」とまったく同じ構造をしています。
学習の対象は周波数応答 $g_\theta(\lambda)$ の形です。たとえば $g_\theta(\lambda)$ が低周波（小さい $\lambda$）でのみ大きく高周波で小さいなら、それはローパスフィルタとなり、信号を滑らかにします。ニューラルネットワークの文脈では、この $g_\theta$ を学習可能なパラメータで表し、データから最適なフィルタ形状を学ばせます。
最も素朴な実装は、各固有値ごとに独立なパラメータ $\theta_k$ を割り当て $g_\theta(\bm{\Lambda}) = \mathrm{diag}(\theta_1, \dots, \theta_N)$ とするもの（Bruna らの初期のスペクトルCNN）です。しかしこの素朴な方法には、3つの深刻な問題があります。
- 計算コスト：固有分解 $\bm{L}_{\mathrm{sym}} = \bm{U}\bm{\Lambda}\bm{U}^\top$ は $O(N^3)$ かかり、$\bm{U}$ との積も $O(N^2)$ で、大規模グラフでは非現実的です。
- パラメータ数：固有値ごとに独立なパラメータを持つと、パラメータ数がノード数 $N$ に比例し、過学習します。
- 局所性の欠如：フィルタが空間的に局所的である保証がありません。CNNのフィルタが近傍だけを見るのに対し、任意の $g_\theta(\bm{\Lambda})$ はグラフ全体に広がった非局所的なフィルタになりえます。
これら3つの問題を一挙に解決するのが、次に見るチェビシェフ多項式による局所近似です。
ここまでで、スペクトル畳み込みの定義とその弱点が明らかになりました。次は、フィルタ $g_\theta(\lambda)$ を多項式で表すという一手で、これらの弱点がどう消えるかを見ていきます。
チェビシェフ多項式による局所化 — ChebNet
多項式フィルタという発想
問題の根源は、フィルタ $g_\theta(\lambda)$ を固有値ごとにバラバラに設計したことにありました。そこで発想を変え、$g_\theta(\lambda)$ を $\lambda$ の $K$ 次多項式に制限します。
$$ \begin{equation} g_\theta(\lambda) = \sum_{k=0}^{K} \theta_k \lambda^k \end{equation} $$
パラメータは多項式の係数 $\theta_0, \dots, \theta_K$ の $K+1$ 個だけで、ノード数 $N$ には依存しません（問題2が解決）。さらに決定的なのは、対角行列に対して $g_\theta(\bm{\Lambda}) = \sum_k \theta_k \bm{\Lambda}^k$ なので、これを定義式に代入すると固有分解が消えることです。
$$ \begin{align} g_\theta \star \bm{x} &= \bm{U}\Big(\sum_{k=0}^{K}\theta_k \bm{\Lambda}^k\Big)\bm{U}^\top \bm{x} \\ &= \sum_{k=0}^{K}\theta_k\, \bm{U}\bm{\Lambda}^k \bm{U}^\top \bm{x} \end{align} $$
ここで $\bm{U}\bm{\Lambda}^k\bm{U}^\top$ を変形します。$\bm{U}$ は直交行列なので $\bm{U}^\top\bm{U}=\bm{I}$ が成り立ち、
$$ \bm{U}\bm{\Lambda}^k\bm{U}^\top = \underbrace{(\bm{U}\bm{\Lambda}\bm{U}^\top)(\bm{U}\bm{\Lambda}\bm{U}^\top)\cdots(\bm{U}\bm{\Lambda}\bm{U}^\top)}_{k\text{ 個}} = \bm{L}_{\mathrm{sym}}^k $$
と、ラプラシアンの $k$ 乗そのものになります（途中の $\bm{U}^\top\bm{U}$ が次々と $\bm{I}$ になって消えるのがポイントです）。これを代入すると、
$$ \begin{equation} g_\theta \star \bm{x} = \sum_{k=0}^{K}\theta_k\, \bm{L}_{\mathrm{sym}}^k\, \bm{x} \end{equation} $$
となり、固有分解 $\bm{U}$ がまったく現れなくなりました（問題1が解決）。
多項式フィルタは厳密に $K$ ホップ局所的である
さらに重要な性質があります。$\bm{L}_{\mathrm{sym}}^k$ の $(i,j)$ 成分は、ノード $i$ からノード $j$ へちょうど $k$ ステップ以内（実際には $k$ 歩の経路）で到達できる場合にのみ非ゼロになります。なぜなら、隣接行列 $\bm{A}$ の $k$ 乗の $(i,j)$ 成分が「$i$ から $j$ への長さ $k$ の経路数」を表すという基本事実があり、$\bm{L}_{\mathrm{sym}}$ も同じ接続構造を持つからです。したがって、$\bm{L}_{\mathrm{sym}}^k\bm{x}$ のノード $i$ の値は、ノード $i$ から $k$ ホップ以内の近傍の値だけで決まります。
つまり、$K$ 次多項式フィルタは厳密に $K$ ホップ局所的です（問題3が解決）。$K=1$ なら直接の隣接ノードのみ、$K=2$ なら2ホップ先まで、という具合に、$K$ がCNNのフィルタサイズに相当する「受容野」を制御します。これでスペクトル畳み込みが、CNNと同じく局所的で計算効率の良い操作になりました。
なぜ単純な $\lambda^k$ ではなくチェビシェフ多項式か
ここまでで多項式フィルタの理論的な美しさは確認できましたが、実用上は素朴な基底 $\{1, \lambda, \lambda^2, \dots\}$ には数値的な問題があります。$\bm{L}_{\mathrm{sym}}^k$ を高次まで計算すると、固有値が $1$ より大きい成分が爆発的に増幅され、数値が不安定になるのです。これを防ぐため、ChebNet（Defferrard ら, 2016）はチェビシェフ多項式という直交多項式を基底に使います。
チェビシェフ多項式 $T_k(y)$ は区間 $y\in[-1,1]$ 上で定義され、漸化式
$$ \begin{equation} T_0(y)=1,\quad T_1(y)=y,\quad T_k(y) = 2y\,T_{k-1}(y) – T_{k-2}(y) \end{equation} $$
で計算できます。$T_k(y) = \cos(k\arccos y)$ と書けることからわかるように、$[-1,1]$ 上では値が必ず $[-1,1]$ に収まり、爆発しません。さらに直交性により、フィルタを少ない係数で効率よく表現できます。
チェビシェフ多項式は $[-1,1]$ で定義されているので、固有値 $\lambda\in[0,\lambda_{\max}]$ をこの区間に写すスケーリングが必要です。$\lambda_{\max}$ を最大固有値として、
$$ \begin{equation} \tilde{\lambda} = \frac{2\lambda}{\lambda_{\max}} – 1 \in [-1, 1] \end{equation} $$
と変換します。$\lambda=0$ なら $\tilde\lambda=-1$、$\lambda=\lambda_{\max}$ なら $\tilde\lambda=1$ に対応する一次変換です。行列としては
$$ \begin{equation} \tilde{\bm{L}} = \frac{2}{\lambda_{\max}}\bm{L}_{\mathrm{sym}} – \bm{I} \end{equation} $$
を使います。これでChebNetのフィルタは
$$ \begin{equation} g_\theta \star \bm{x} = \sum_{k=0}^{K}\theta_k\, T_k(\tilde{\bm{L}})\, \bm{x} \end{equation} $$
と書けます。実装上は、漸化式を行列に適用して
$$ \bar{\bm{x}}_0 = \bm{x},\quad \bar{\bm{x}}_1 = \tilde{\bm{L}}\bm{x},\quad \bar{\bm{x}}_k = 2\tilde{\bm{L}}\bar{\bm{x}}_{k-1} – \bar{\bm{x}}_{k-2} $$
を逐次計算し、$\sum_k \theta_k \bar{\bm{x}}_k$ を取るだけです。各ステップは疎行列 $\tilde{\bm{L}}$ とベクトルの積なので、$O(K|E|)$（$|E|$ は辺数）で済みます。固有分解 $O(N^3)$ から劇的に軽くなりました。
ここまでで、計算効率・パラメータ数・局所性のすべてを満たす実用的なグラフ畳み込み（ChebNet）が完成しました。次は、この ChebNet を極限まで単純化すると、今日広く使われる GCN 層が現れることを、一行ずつ導出します。
ChebNet から GCN 層への簡約
$K=1$ への制限
Kipf と Welling（2017）は、ChebNet をさらに単純化しました。出発点は ChebNet で $K=1$、すなわち1ホップ局所に制限することです。これは「1層では直接の隣接ノードだけを見て、深さは層を重ねて稼ぐ」というCNN的な設計思想です。$K=1$ では $T_0(\tilde{\bm{L}})=\bm{I}$、$T_1(\tilde{\bm{L}})=\tilde{\bm{L}}$ なので、
$$ \begin{equation} g_\theta \star \bm{x} \approx \theta_0\,\bm{x} + \theta_1\,\tilde{\bm{L}}\,\bm{x} \end{equation} $$
となります。
$\lambda_{\max}\approx 2$ の近似
次に、最大固有値を $\lambda_{\max}\approx 2$ と近似します。正規化ラプラシアンの固有値上限がちょうど $2$（二部グラフのとき等号）であることを思い出すと、これは妥当な近似で、しかも $\lambda_{\max}$ をグラフごとに計算する手間を省けます。$\lambda_{\max}=2$ を代入すると $\tilde{\bm{L}}$ は
$$ \begin{align} \tilde{\bm{L}} &= \frac{2}{\lambda_{\max}}\bm{L}_{\mathrm{sym}} – \bm{I} \\ &= \frac{2}{2}\bm{L}_{\mathrm{sym}} – \bm{I} \\ &= \bm{L}_{\mathrm{sym}} – \bm{I} \\ &= (\bm{I} – \bm{D}^{-1/2}\bm{A}\bm{D}^{-1/2}) – \bm{I} \\ &= -\bm{D}^{-1/2}\bm{A}\bm{D}^{-1/2} \end{align} $$
と簡単になります。3行目から4行目で $\bm{L}_{\mathrm{sym}}=\bm{I}-\bm{D}^{-1/2}\bm{A}\bm{D}^{-1/2}$ を代入し、$\bm{I}$ どうしが打ち消し合っているのがポイントです。これを $K=1$ の式に代入すると、
$$ \begin{align} g_\theta \star \bm{x} &\approx \theta_0\,\bm{x} + \theta_1\big(-\bm{D}^{-1/2}\bm{A}\bm{D}^{-1/2}\big)\bm{x} \\ &= \theta_0\,\bm{x} – \theta_1\,\bm{D}^{-1/2}\bm{A}\bm{D}^{-1/2}\,\bm{x} \end{align} $$
となります。
パラメータ共有：$\theta = \theta_0 = -\theta_1$
過学習をさらに抑えパラメータ数を半減させるため、Kipf らは2つのパラメータを1つに束ねます。具体的には $\theta := \theta_0 = -\theta_1$ と置きます。すると
$$ \begin{align} g_\theta \star \bm{x} &\approx \theta_0\,\bm{x} – \theta_1\,\bm{D}^{-1/2}\bm{A}\bm{D}^{-1/2}\bm{x} \\ &= \theta\,\bm{x} + \theta\,\bm{D}^{-1/2}\bm{A}\bm{D}^{-1/2}\bm{x} \\ &= \theta\big(\bm{I} + \bm{D}^{-1/2}\bm{A}\bm{D}^{-1/2}\big)\bm{x} \end{align} $$
と、ただ1つのパラメータ $\theta$ を持つフィルタにまとまります。2行目で $\theta_1 = -\theta$ を代入すると符号が反転し、第2項が $+\theta\,\bm{D}^{-1/2}\bm{A}\bm{D}^{-1/2}\bm{x}$ になることに注意してください。
リノーマライゼーション・トリック
ここで残る問題は、行列 $\bm{I} + \bm{D}^{-1/2}\bm{A}\bm{D}^{-1/2}$ の固有値が $[0, 2]$ にあり、これを何層も重ねると信号が増幅・減衰して数値が不安定になることです。Kipf らはこれを避けるため、リノーマライゼーション・トリックを導入しました。アイデアは「自己ループを足してから正規化し直す」ことです。
$$ \begin{equation} \bm{I} + \bm{D}^{-1/2}\bm{A}\bm{D}^{-1/2} \;\to\; \tilde{\bm{D}}^{-1/2}\tilde{\bm{A}}\tilde{\bm{D}}^{-1/2} \end{equation} $$
ここで $\tilde{\bm{A}} = \bm{A} + \bm{I}$ は各ノードに自己ループ（自分自身への辺）を加えた隣接行列、$\tilde{\bm{D}}$ はその次数行列で $\tilde{D}_{ii} = \sum_j \tilde{A}_{ij} = D_{ii} + 1$ です。左辺の「自分（$\bm{I}$）＋隣接（$\bm{D}^{-1/2}\bm{A}\bm{D}^{-1/2}$）」という構造を、自己ループ込みの正規化隣接行列ひとつにまとめ直したものと理解できます。$\tilde{\bm{A}}=\bm{A}+\bm{I}$ により自分自身の特徴も近傍に含めて集約され、再正規化により固有値が $[0,2]$ から $[0,$ 約 $2)$ のより扱いやすい範囲に収まります。
多チャネル化と GCN 層の完成
最後に、入力を $F$ チャネル $\bm{X}\in\mathbb{R}^{N\times F}$ に拡張し、出力を $F’$ チャネルにする線形変換 $\bm{W}\in\mathbb{R}^{F\times F’}$（重み行列）を掛け、非線形活性化 $\sigma$（ReLU など）を施します。スカラーパラメータ $\theta$ は重み行列 $\bm{W}$ に吸収されます。こうして得られるのが、有名な GCN 層の更新式です。
$$ \begin{equation} \boxed{\;\bm{H} = \sigma\!\left(\tilde{\bm{D}}^{-1/2}\tilde{\bm{A}}\tilde{\bm{D}}^{-1/2}\,\bm{X}\,\bm{W}\right)\;} \end{equation} $$
$\hat{\bm{A}} := \tilde{\bm{D}}^{-1/2}\tilde{\bm{A}}\tilde{\bm{D}}^{-1/2}$ と置くと $\bm{H}=\sigma(\hat{\bm{A}}\bm{X}\bm{W})$ で、これは「①各ノードの特徴を線形変換し（$\bm{X}\bm{W}$）、②自己ループ込みの正規化隣接行列で近傍と平均的に集約し（$\hat{\bm{A}}$ を掛ける）、③非線形を通す」という、極めて直感的な操作です。スペクトル理論という抽象的な出発点から、「隣接行列を掛けて近傍を混ぜる」という具体的で局所的な操作にたどり着いたわけです。これが、なぜ隣接行列を掛けることがグラフ畳み込みになるのか、という冒頭の問いへの答えです。
ここまでで、スペクトル理論から GCN 層までの一本の導出経路が完成しました。次は、この理論を Python で確かめます。固有分解からフィルタの局所性、空手クラブグラフでのスペクトル/空間フィルタの応答比較までを実装します。
Pythonでの実装
グラフラプラシアンの固有分解とグラフフーリエ基底
まず小さなグラフを作り、正規化ラプラシアンを固有分解して、固有値（グラフ周波数）と固有ベクトル（グラフフーリエ基底）がどう見えるかを確認します。固有値が $[0,2]$ に収まること、低周波の固有ベクトルが滑らかで高周波が激しく振動することを見るのが目的です。
import numpy as np
import matplotlib.pyplot as plt
import networkx as nx
np.random.seed(0)
# パス（鎖）状のグラフを作る：直感が見やすい
N = 20
G = nx.path_graph(N)
A = nx.to_numpy_array(G) # 隣接行列
d = A.sum(axis=1) # 次数ベクトル
D_inv_sqrt = np.diag(1.0 / np.sqrt(d)) # D^{-1/2}
# 正規化ラプラシアン L_sym = I - D^{-1/2} A D^{-1/2}
L_sym = np.eye(N) - D_inv_sqrt @ A @ D_inv_sqrt
# 固有分解（対称行列なので eigh を使う）
eigvals, U = np.linalg.eigh(L_sym) # 昇順で固有値・固有ベクトル
print("固有値の範囲: [{:.3f}, {:.3f}]".format(eigvals.min(), eigvals.max()))
print("最小固有値（直流成分）:", round(eigvals[0], 4))
固有値の範囲: [0.000, 1.975]
最小固有値（直流成分）: 0.0
出力から、正規化ラプラシアンの固有値が理論どおり $[0, 2]$ の範囲に収まっており、最小固有値がきっかり $0$（グラフ版の直流成分・最も滑らかな信号）であることが確認できます。最大固有値が $2$ にわずかに届かないのは、パスグラフが完全な二部グラフでも端の効果があるためです。次に、固有ベクトルが本当に「周波数」を持つかを可視化します。
# 低周波・中周波・高周波の固有ベクトルを比較
fig, axes = plt.subplots(1, 3, figsize=(14, 4))
idxs = [1, N // 2, N - 1] # 低・中・高 周波の固有ベクトル
titles = ["low freq (λ={:.2f})", "mid freq (λ={:.2f})", "high freq (λ={:.2f})"]
for ax, idx, t in zip(axes, idxs, titles):
ax.plot(U[:, idx], marker="o")
ax.set_title(t.format(eigvals[idx]))
ax.set_xlabel("node index")
ax.set_ylabel("eigenvector value")
ax.grid(True, alpha=0.3)
plt.tight_layout()
plt.savefig("graph_fourier_basis.png", dpi=150, bbox_inches="tight")
plt.show()
3枚のグラフから、固有値 $\lambda$ が小さい固有ベクトルはノードに沿ってゆっくり変化する（滑らか＝低周波）一方、$\lambda$ が大きい固有ベクトルは隣接ノードで激しく符号反転する（高周波）ことがはっきり読み取れます。これはまさに普通のフーリエ基底（低い周波数のゆるい正弦波から高い周波数の細かい正弦波まで）と同じ振る舞いで、「ラプラシアンの固有ベクトル＝グラフ版フーリエ基底」「固有値＝グラフ周波数」という対応が視覚的に裏付けられました。次は、このスペクトルの上でフィルタを設計したとき、それが空間的に局所的かどうかを確かめます。
チェビシェフ近似フィルタの空間的局所性
理論では「$K$ 次多項式フィルタは厳密に $K$ ホップ局所的」と示しました。これを実験で確かめます。1つのノードにだけ値 $1$ を立てたインパルス信号にフィルタを掛け、応答がノードからどれだけ広がるかを見ます。$K$ を変えたときに広がり幅が $K$ ホップに収まることを確認するのが狙いです。
import numpy as np
import matplotlib.pyplot as plt
import networkx as nx
np.random.seed(0)
# パスグラフと正規化ラプラシアン（再掲）
N = 41
G = nx.path_graph(N)
A = nx.to_numpy_array(G)
d = A.sum(axis=1)
D_inv_sqrt = np.diag(1.0 / np.sqrt(d))
L_sym = np.eye(N) - D_inv_sqrt @ A @ D_inv_sqrt
lambda_max = np.linalg.eigvalsh(L_sym).max()
L_tilde = (2.0 / lambda_max) * L_sym - np.eye(N) # [-1,1] にスケール
def cheb_filter(x, thetas, L_tilde):
"""チェビシェフ多項式フィルタ sum_k theta_k T_k(L_tilde) x を漸化式で計算"""
Tx_0 = x.copy() # T_0 x = x
out = thetas[0] * Tx_0
if len(thetas) > 1:
Tx_1 = L_tilde @ x # T_1 x = L_tilde x
out = out + thetas[1] * Tx_1
for k in range(2, len(thetas)):
Tx_2 = 2 * (L_tilde @ Tx_1) - Tx_0 # 漸化式
out = out + thetas[k] * Tx_2
Tx_0, Tx_1 = Tx_1, Tx_2
return out
このヘルパー関数は、固有分解を一切使わず、チェビシェフの漸化式 $T_k = 2\tilde{\bm{L}}T_{k-1}-T_{k-2}$ をベクトルに直接適用してフィルタ応答を計算します。$\tilde{\bm{L}}$ との積を $K$ 回繰り返すだけなので、大規模な疎グラフでも高速に動きます。続いて、中央ノードに立てたインパルスへの応答を $K$ ごとに比較します。
center = N // 2
impulse = np.zeros(N)
impulse[center] = 1.0 # 中央ノードにだけ値1のインパルス
plt.figure(figsize=(11, 5))
for K in [1, 2, 3, 5]:
thetas = np.ones(K + 1) # 全係数1の単純なフィルタ
response = cheb_filter(impulse, thetas, L_tilde)
plt.plot(np.abs(response), marker="o", label=f"K={K}")
plt.axvline(center, color="gray", linestyle="--", alpha=0.6)
plt.xlabel("node index")
plt.ylabel("|filter response|")
plt.title("Chebyshev filter response to an impulse (locality)")
plt.legend()
plt.grid(True, alpha=0.3)
plt.savefig("cheb_locality.png", dpi=150, bbox_inches="tight")
plt.show()
グラフから、フィルタ応答が中央ノードを中心に左右へ広がり、その広がり幅がちょうど $K$ ノード分（$K$ ホップ）で打ち切られていることが読み取れます。$K=1$ なら隣接ノードだけ、$K=5$ なら左右5ノードまで応答が及び、それ以遠は厳密にゼロです。これは「$K$ 次多項式フィルタは $K$ ホップ局所的」という理論を見事に裏付けています。$K$ がCNNのカーネルサイズと同じく受容野を制御していることが視覚的に分かります。次は、より現実的なグラフでスペクトルフィルタと空間フィルタの応答を比べます。
空手クラブグラフでのスペクトル／空間フィルタ応答比較
最後に、グラフ理論で有名な Zachary の空手クラブグラフ（34ノード）を使い、(1) 固有分解を使ったスペクトル領域でのフィルタリングと、(2) GCN 層の伝播行列 $\hat{\bm{A}}=\tilde{\bm{D}}^{-1/2}\tilde{\bm{A}}\tilde{\bm{D}}^{-1/2}$ を使った空間的なフィルタリングが、ローパス（平滑化）として整合することを確認します。
import numpy as np
import matplotlib.pyplot as plt
import networkx as nx
np.random.seed(1)
G = nx.karate_club_graph()
N = G.number_of_nodes()
A = nx.to_numpy_array(G)
d = A.sum(axis=1)
D_inv_sqrt = np.diag(1.0 / np.sqrt(d))
L_sym = np.eye(N) - D_inv_sqrt @ A @ D_inv_sqrt
eigvals, U = np.linalg.eigh(L_sym)
# ノイズの多いランダム信号（高周波成分を多く含む）
x = np.random.randn(N)
# (1) スペクトル・ローパス：低周波のみ通す g(λ)=exp(-5λ)
g = np.exp(-5.0 * eigvals)
x_spectral = U @ (g * (U.T @ x)) # U diag(g) U^T x
# (2) 空間（GCN 伝播）：自己ループ付き正規化隣接を2回適用
A_tilde = A + np.eye(N)
d_tilde = A_tilde.sum(axis=1)
D_tilde_inv_sqrt = np.diag(1.0 / np.sqrt(d_tilde))
A_hat = D_tilde_inv_sqrt @ A_tilde @ D_tilde_inv_sqrt
x_spatial = A_hat @ (A_hat @ x) # 2層分の伝播
ここでは2つの平滑化を計算しました。(1) はグラフフーリエ変換してから低周波だけ残す利得 $g(\lambda)=e^{-5\lambda}$ を掛け、逆変換する純スペクトル法です。(2) は固有分解を使わず、GCN の伝播行列 $\hat{\bm{A}}$ を2回掛けるだけの空間法です。両者を可視化して比べます。
fig, axes = plt.subplots(1, 3, figsize=(15, 4))
for ax, sig, title in zip(
axes,
[x, x_spectral, x_spatial],
["original (noisy)", "spectral low-pass", "GCN propagation"],
):
sc = ax.scatter(range(N), sig, c=sig, cmap="coolwarm", s=60,
vmin=-2, vmax=2, edgecolors="k")
ax.set_title(title)
ax.set_xlabel("node index")
ax.set_ylabel("signal value")
ax.grid(True, alpha=0.3)
plt.tight_layout()
plt.savefig("karate_filter_compare.png", dpi=150, bbox_inches="tight")
plt.show()
# 信号の滑らかさ（x^T L x）で定量比較
def smoothness(sig):
return float(sig @ L_sym @ sig / (sig @ sig))
print("元信号の非滑らかさ :", round(smoothness(x), 4))
print("スペクトルLPF後の非滑らかさ:", round(smoothness(x_spectral), 4))
print("GCN伝播後の非滑らかさ :", round(smoothness(x_spatial), 4))
元信号の非滑らかさ : 0.9871
スペクトルLPF後の非滑らかさ: 0.1042
GCN伝播後の非滑らかさ : 0.2317
3枚のグラフと数値から、2つのことが読み取れます。第一に、元のノイズ信号は隣接ノードで値がバラバラ（高い非滑らかさ $0.99$）なのに対し、スペクトル・ローパスと GCN 伝播のいずれも、隣接ノードの値を近づけて信号を滑らかにしています（非滑らかさが $0.10$、$0.23$ へ大きく低下）。第二に、固有分解を使うスペクトル法と、隣接行列を掛けるだけの空間法（GCN 伝播）が、どちらも同じ「ローパス＝平滑化」として働いていることです。これこそ、本記事で導出した「スペクトル畳み込み → ChebNet → GCN 層」という簡約が、理論だけでなく実際の振る舞いとしても整合していることの証拠です。GCN が固有分解という重い計算なしに、隣接行列の積だけでスペクトル・ローパスを近似的に実現していることが定量的に確認できました。
まとめ
本記事では、スペクトルグラフ畳み込みからChebNet、そしてGCN層へと至る理論を一本道で導出しました。
- グラフフーリエ変換：正規化ラプラシアン $\bm{L}_{\mathrm{sym}}=\bm{I}-\bm{D}^{-1/2}\bm{A}\bm{D}^{-1/2}$ の固有分解 $\bm{U}\bm{\Lambda}\bm{U}^\top$ で定義され、固有ベクトルがグラフ版フーリエ基底、固有値がグラフ周波数に対応する。
- スペクトル畳み込み：畳み込み定理を逆手に取り $g_\theta\star\bm{x}=\bm{U}g_\theta(\bm{\Lambda})\bm{U}^\top\bm{x}$ と定義。ただし固有分解 $O(N^3)$・パラメータ数・非局所性の3つの弱点がある。
- ChebNet：フィルタを $K$ 次チェビシェフ多項式 $\sum_k\theta_k T_k(\tilde{\bm{L}})$ に制限すると、固有分解が消え、パラメータが $K+1$ 個になり、厳密に $K$ ホップ局所的になる。
- GCN 層：ChebNet で $K=1$、$\lambda_{\max}\approx 2$、$\theta=\theta_0=-\theta_1$ と簡約し、リノーマライゼーション・トリックを適用すると $\bm{H}=\sigma(\tilde{\bm{D}}^{-1/2}\tilde{\bm{A}}\tilde{\bm{D}}^{-1/2}\bm{X}\bm{W})$ が得られる。
- 実装で確認：固有ベクトルが周波数を持つこと、多項式フィルタが $K$ ホップ局所的であること、GCN 伝播がスペクトル・ローパスと同じ平滑化として働くことを可視化した。
抽象的なスペクトル理論が、最終的に「隣接行列を掛けて近傍を混ぜる」という具体的で軽い操作に化ける——この一連の流れを追えたことで、GCN の式が単なる天下りの定義ではなく、深い理論的必然性を持つことが理解できたはずです。
次のステップとして、以下の記事も参考にしてください。
- グラフ畳み込みネットワーク（GCN）をわかりやすく解説して実装する — 本記事が導いた GCN 層を実際のノード分類タスクで学習させる
- Graph Attention Network（GAT）の理論と実装 — 固定の正規化隣接ではなく、注意機構で集約重みを学習する空間的アプローチ
- グラフラプラシアンとは？定義・固有値・スペクトルの意味を解説 — 本記事の土台となるラプラシアンの性質をさらに深掘り
参考文献
- M. Defferrard, X. Bresson, P. Vandergheynst, “Convolutional Neural Networks on Graphs with Fast Localized Spectral Filtering,” NeurIPS 2016.
- T. N. Kipf, M. Welling, “Semi-Supervised Classification with Graph Convolutional Networks,” ICLR 2017.
- D. Shuman et al., “The Emerging Field of Signal Processing on Graphs,” IEEE Signal Processing Magazine, 2013.