Laplacian行列

今回はLaplacian行列について記事を書きます.Laplacian行列は機械学習でもデータのクラスタリングなどにしばしば用いられるので,情報系の学生にとってはそこそこ馴染みのあるものだと思いますが,Laplacian行列をちゃんと(導出とか物理学的解釈とか)説明している日本語の教科書やweb記事は私の知る限りほぼ存在しないように思います.おかげで私もB4の頃は勉強にとても苦労しました...... 今回は工学(特に,情報工学や数理工学)の学生を対象にしているので万人向けではないですが,理工系の学部生程度の知識があれば十分読めると思います.タブン.

前口上はさておき,早速説明に移りましょう.

Laplacian行列はグラフラプラシアン,キルヒホッフ行列とも呼ばれますが,ここではLaplacian行列で統一します.ここでいう"グラフ"というのは,散布図やヒストグラムなどのグラフではなく,グラフ理論で扱われるグラフ,つまり事物や人の関係性を記述する数学概念です.

グラフはノードとリンクによって構成されます.ノード集合を V=\{1,2,\dots,n\},リンク集合を E \subset V \times Vとすると,グラフはG=G(V,E)と書くことができます.さらに,ノード i-j間のリンク(i,j)の方向を区別する場合は有向グラフ,区別しない場合は無向グラフと呼ばれます.以下では,グラフGは無向グラフとして話を進めます.

グラフを数学的に取り扱うため,グラフ理論ではグラフ構造を行列を用いて表現します.グラフ構造の行列表現として,国内のグラフ理論の教科書では,接続行列や隣接行列が挙げられることが多いですが,一般に接続行列は非正方行列となり扱いづらいので,実用的には隣接行列の方がよく用いられます.グラフGの隣接行列 \boldsymbol{A}=[A_{ij}]の各成分 A_{ij}はリンクの有無を表し,ノード i-j間にリンクがあれば A_{ij}=1,リンクがなければ A_{ij}=0となります.

さて,これで準備が整ったのでいよいよLaplacian行列の導出に移りましょう.このために,次のようなグラフ上の拡散過程を考えます.

ノード i上の何らかの分布量 \psi_{i}がリンクを介してグラフ上を拡散する状況を考えます.なお,ノード jからノード iへの流入は差 \kappa(\psi_{j}-\psi_{i}) (\kappa>0) に比例すると仮定します(差 \kappa(\psi_{j}-\psi_{i})が負になる場合はノード iからノード jへの流出を意味します).このとき,ノード iの分布量 \psi_{i}の時間発展は以下のように記述できます.

 \frac{d\psi_{i}}{dt} = \kappa\sum_{j=1}^{n}A_{ij}(\psi_{j}-\psi_{i}) = \kappa\sum_{j=1}^{n}(A_{ij}-\delta_{ij}d_{i})\psi_{j}

ここで, \delta_{ij}クロネッカーのデルタ, d_{i}:=\sum_{j=1}^{n}A_{ij}はノード iの次数を表します.上式をベクトルで書くと

 \frac{d\boldsymbol{\psi}}{dt} = \kappa(\boldsymbol{A}-\boldsymbol{D})\boldsymbol{\psi}

となります.ただし, \boldsymbol{D}:=\mathrm{diag}(d_{1},d_{2},\dots,d_{n})は次数行列です.ここで \boldsymbol{L}:=\boldsymbol{D} - \boldsymbol{A} とすると,グラフ上の拡散方程式として

 \frac{d\boldsymbol{\psi}}{dt} = -\kappa \, \boldsymbol{L} \, \boldsymbol{\psi}

以下が得られます.上式を連続空間における拡散方程式 \frac{du}{dt}=\kappa \, \Delta \, uと比較すると,行列\boldsymbol{L}ラプラス演算子 \Deltaに対応することがわかります.これが行列 \boldsymbol{L}がLaplacian行列と呼ばれる由縁です.

ところで,グラフ上の拡散方程式の右辺にはマイナスが付いていますが,これは単に定義の問題で本質には影響しません.論文によっては \boldsymbol{L}:=\boldsymbol{A} - \boldsymbol{D} と定義しているものもあります.一般に \boldsymbol{L}:=\boldsymbol{D} - \boldsymbol{A} とされているのは,このように定義すると \boldsymbol{L}が非負正定値行列になって便利というだけの理由だと思います.

さて,ここまでの議論でグラフ上の拡散現象がLaplacian行列によって記述できることがわかりました.しかし,Laplacian行列が本質的にグラフにおけるラプラス演算子であることを考えれば,拡散の他にも様々な現象を表現できそうです.私の知る限り,グラフ上の拡散方程式の他に,反応拡散方程式,回路方程式,振動方程式などがLaplacian行列によって記述できたと記憶しています.以下では,その一例としてグラフ上の振動方程式を紹介したいと思います.

 n個のノードがバネで結合されており,隣接ノード同士で相互に影響を与え合いながら振動する連成振動系を考えます.時刻 tにおけるノード iの平衡点からの変位を x_{i}=x_{i}(t)とし,復元力はノード間の変位の差に比例すると仮定します.ノードの質量を m,バネ係数を kとし,一般化運動量を p_{i}=p_{i}(t)とすると,この系のハミルトニアン Hは以下のように書けます.

 H = \sum_{i \in V}\frac{(p_{i})^{2}}{2m}+\frac{k}{2}\sum_{(i,j) \in E}(x_{i}-x_{j})^{2}= \sum_{i \in V}\frac{(p_{i})^{2}}{2m}+\frac{k}{2}( {}^{t}\boldsymbol{x} \boldsymbol{L} \boldsymbol{x} )

上式右辺の第1項は運動エネルギー,第2項はポテンシャルエネルギーを表しています.

さて,ハミルトニアンから正準方程式を導出すると

 \frac{dp_{i}}{dt} = -\frac{\partial H}{\partial x_{i}} = - k \sum_{j=1}^{n}L_{ij}x_{j}

 \frac{dx_{i}}{dt} = \frac{\partial H}{\partial p_{i}} = \frac{p_{i}}{m}

となります. 上式をまとめて p_{i}を消去すると,ノード iに関する振動方程式

 m \frac{d^{2}x_{i}}{dt^{2}}=- k \sum_{j=1}^{n}L_{ij}x_{j}

が得られます.これをベクトル形式で表現するとグラフ上の振動方程式

 m \frac{d^{2}\boldsymbol{x}}{dt^{2}}=- k \boldsymbol{L} \, \boldsymbol{x}

 を得ます.ここでは,面倒くさかったため簡単のため質量やバネ係数を全て同一としましたが,ノード毎に質量が異なる場合やリンク毎にバネ係数が異なる場合にも容易に拡張可能です.

グラフ上の振動方程式の応用先として代表的なものにグラフの可視化があります.確かPythonのNetworkXモジュールのグラフレイアウトの一つとして利用されていた気がします.また,分子動力学の分野においてタンパク質間相互作用を近似的に記述するためにちょうど上記のような振動モデルを導入していたと思います.

 

改めて読み返すと説明の拙い部分が随所に見られますが,ちゃんと伝わるでしょうか.この記事を読んだ方々がLaplacian行列について少しでも理解を深めていただけたなら僥倖です.にゃん.