【图算法】图卷积的演变-从谱图卷积到GCN

基础

傅里叶变换

傅里叶级数是对周期为T的确定性信号做展开，而傅里叶变换将周期推广到无穷，能对具有任意长度的信号做展开。

傅里叶级数和傅里叶变换是什么关系？
如下为傅里叶变换公式：

\[\hat{f}(t)={\int}f(x){\exp}^{-iwt}dx = {\int}f(x) \left(cos(wx) + isin(wx) \right)dx \]

用欧拉公式将\({\exp}^{-iwt}\)展开后，可发现相当于用不同频率的正弦和余弦信号作为基向量和\(f(x)\)做内积，从而将其从时域变到频域空间。

图傅里叶变换

要将傅里叶变换推广到图上，其关键是找到图信号的基函数。

拉普拉斯算子(Laplacian operator)\(∆\)的物理意义是空间二阶导数，其准确定义是：标量梯度场中的散度（梯度的散度，二阶偏导之和），可用于描述物理量的流入流出，例如热传播。

传统傅里叶的基函数\({\exp}^{-iwt}\)可视为拉普拉斯算子的特征向量，频率为特征值：

\[∆{e}^{-iwt}=\frac{{\partial}^2}{{\partial}t^2}{e}^{-iwt}={-w^2}{e}^{-iwt} \]

其可视为广义的特征方程，和矩阵的特征向量类似（在经过变换后只是改变大小，而方向不变）。

图拉普拉斯矩阵\(\mathbf{L}\)是拉普拉斯算子的在图（离散空间）上的的推广，其可用于衡量图上信号的平滑程度。

那么，\(\mathbf{L} \in \mathbb{R}^{N \times N}\)的特征向量就可以类比基函数 \({\exp}^{-iwt}\)，作为图傅里叶变换的基函数(向量)。

\[\mathbf{L} = \mathbf{D} - \mathbf{A} = \mathbf{U} \mathbf{\Lambda} \mathbf{U}^{T}, 拉普拉斯矩阵特征分解 \]

\[\mathbf{U} = (\mathbf{u}_1, \mathbf{u}_2, {\cdots}, \mathbf{u}_N) \in \mathbb{R}^{N \times N}, 拉普拉斯矩阵特征向量 \]

\[ \mathbf{U}^{-1} =\mathbf{U}^T, \mathbf{U} \mathbf{U}^T = \mathbf{I} \]

\[ \mathbf{\Lambda}, 拉普拉斯矩阵特征值对角矩阵 \]

\[ \phi_l = \mathbf{u}_l^{T} \mathbf{x},基向量上的分量 \]

一个包含\(N\)个节点图，图上每个节点\(i\)有一个标量信号，这时图信号\(\mathbf{x} \in \mathbb{R}^{N}\)都可表示为拉普拉斯矩阵特征向量（基向量）的线性组合。

\[ \mathbf{\Phi} = \left[ \begin{matrix} \phi_1 \\ ... \\ \phi_N \end{matrix} \right]=\mathbf{U}^T \mathbf{x}， 图傅里叶变换 \]

\[ \mathbf{x} ={\sum}_{l} \phi_l \mathbf{u}_l = \mathbf{U} \mathbf{\phi}，图傅里叶逆变换 \]

图傅里叶变换，在这里就是将图信号\(\mathbf{x}\)投影（内积计算分量）到\(\mathbf{L}\)的特征向量构成的基向量上。就是将\(\mathbf{x}\)从原始空间变到新的空间-频域。

\[ \mathbf{\Phi}=\mathbf{U}^T \mathbf{x}= \mathbf{U}^T (\mathbf{U} \mathbf{\Phi}), \mathbf{x} = \mathbf{U} \mathbf{\Phi} = \mathbf{U} (\mathbf{U}^T \mathbf{x}) \]

图傅里叶变换(Graph Fourier Transformation）就是基于图拉普拉斯矩阵，将图信号从空域(顶点上)转换到谱域(频域）的一种方法。

第一代：Spectral Network

卷积

卷积运算: 两个实值函数的卷积运算可以理解成，以其中一个函数为权重，对另一个函数做加权平均的操作（这样可以令函数平滑降噪）。

卷积运算的目的不限于此。其操作可以得到一个新的函数，相当于是这两个函数内积的结果。

例子，给定\(f(t)\)表示时刻\(t\)的测量值。由于测量值可能存在噪声，且时间上越近的测量结果越相关，可使用加权方法对最近测量值赋予高权重，来获得加权平滑结果。

可用一个加权函数\(g(a)\)来实现，其中\(a\)代表测量值距当前时刻的时间间隔：

\[s(t) = (f * g)(t) = \sum_{a=-\infty}^{+\infty} f(a) g(t-a)， 离散情况 \]

若\(t\)为目标时刻，当\(a\)取\(t\)时，间隔为\(t-a=0\)，有\(f(t) g(0)\)，如下表所示（实际大于\(t\)取不到）：

\(a\)	\(t-3\)	\(t-2\)	\(t-1\)	\(t\)
\(f(a)\)	\(f(t-3)\)	\(f(t-2)\)	\(f(t-1)\)	\(f(t)\)
\(g(t-a)\)	\(g(3)\)	\(g(2)\)	\(g(1)\)	\(g(0)\)

\[s(t) =(f * g)(t) = \int f(a) g(t-a) d{a}，连续情况 \]

卷积定理：卷积的傅里叶变换等于傅里叶变换的乘积（时域卷积，等于在频域做乘积）

\[ F\{f*g\} = F[f] {\odot} F[g]] \]

通过傅里叶逆变换可以得到:

\[f*g = F^{-1}[F[f]{\odot}F[g]] \]

图上的卷积

在图上做，图信号和滤波器g的卷积：

\[输入 \mathbf{x} \in \mathbb{R}^{N}，每一个节点有一个标量 \]

\[ \mathbf{g} \in \mathbb{R}^{N}，滤波器向量 \]

那么，图上的卷积可以定义为：

\[\mathbf{x} \star \mathbf{g} = \mathbf{U} \left( (\mathbf{U}^{T}\mathbf{x}) {\odot} (\mathbf{U}^{T}\mathbf{g})\right) = \mathbf{U} (\mathbf{U}^{T}\mathbf{x} {\odot} \mathbf{\theta}) \]

\[把 \mathbf{U}^{T}\mathbf{g}统一视为一个，\mathbf{\theta} \in \mathbb{R}^{N} (傅里叶变换后的滤波器 \mathbf{g}) \]

传统滤波器需根据经验设定，在这里可将滤波器视为：可参数化的卷积核

\[(\mathbf{U}^{T}\mathbf{x}) {\odot} \mathbf{\theta}，\mathbf{\theta} {\odot} (\mathbf{U}^{T}\mathbf{x})，交换顺序不影响 \]

卷积运算中的乘法为element-wise product，即在频域的乘法。在这里，其直观意义就是：

用卷积核的参数对频域信号的每个分量进行加权操作，来实现滤波（不同的频率分量有不同的权重系数，例如可对高频分量施以更低权重）

那么将卷积核向量展开为对角矩阵形式（行变换），有：

\[\mathbf{g}_{\theta} = diag{(\mathbf{\theta}) } = \left[ \begin{matrix} {\theta}_1 & ... & 0 \\ ... & ... & ... \\ 0 & ... & {\theta}_N \end{matrix} \right] \]

最后，可得到：

\[\mathbf{x} \star \mathbf{g} = \mathbf{U} (\mathbf{U}^{T}\mathbf{x} {\odot} \mathbf{\theta} ) \\ = \mathbf{U} (\mathbf{\theta} {\odot} \mathbf{U}^{T}\mathbf{x} ) \\ = \mathbf{U} \mathbf{g}_{\theta} \mathbf{U}^{T} \mathbf{x} \]

假设每个节点有\(d\)维的特征，即通道数为\(d\)（\(d\)个图信号）:

\[\mathbf{X} = \left[ \begin{matrix} {x_{11}} & {x_{12}} & ... & {x_{1d}} \\ ... & ... & ... \\ {x_{n1}} & {x_{n2}} & ... & {x_{nd}} \end{matrix} \right] = \left[ \begin{matrix} {\mathbf{x}_1} & {\mathbf{x}_2} & ... & {\mathbf{x}_d} \end{matrix} \right] \]

注意，\(\mathbf{X} \in \mathbb{R}^{N \times d}\), 每一个通道可使用多个卷积核（类似CNN，拓展通道数）。

对于第\(l\)层谱图卷积，通道数为\({d_l}\)：

\[假设第 {l} 和{l+1}层的节点状态为： \mathbf{X}^{(l)} \in \mathbb{R}^{N \times d_l}， \mathbf{X}^{(l+1)} \in \mathbb{R}^{N \times d_{l+1}} \]

\[\mathbf{X}^{(l)}_{:i} =\mathbf{x}^{(l)}_{i} \in \mathbb{R}^{N} \]

使用\(d_l * d_{(l+1)}\)个卷积核，每次在全部通道分别用\(d_l\)个卷积核并将结果求和，重复\(d_{(l+1)}\)次，得到输出特征通道：

\[\mathbf{x}^{l+1}_j={\sigma}(\mathbf{U} {\sum}_{i=1}^{d_l} \mathbf{\Theta}^l_{i,j} \mathbf{U}^T \mathbf{x}^l_i), (j = 1, \dots, d_{(l+1)}) \]

\[\mathbf{\Theta}^l_{i,j}直接视为模型参数， \mathbf{U} \mathbf{\Theta}^l_{i,j} \mathbf{U}^T对应CNN中的卷积核(复杂度 O(n^2)) \]

Spectral Graph Convolution操作定义为：

• 计算图拉普拉斯(graph Laplacian)的特征值分解，得到特征向量
• 将图信号进行图傅里叶变换， 然后使用卷积核进行滤波，然后再进行图傅里叶逆变换

缺点：

• 图拉普拉斯特征分解\(O(n^3)\)复杂度， 前向传播\(O(n^2)\)。
• 卷积核参数量大： \(N * d_l * d_(l+1)\), 易过拟合($N $ 为节点数量)
• 在空域上没有明确定义，不能局部化到节点上

基于谱图卷积的频域方法，学到的滤波器都是基于拉普拉斯特征分解，也就是取决于图的结构。这也就意味着，在一个特定结构上训练得到的模型，并不能直接应用到另外一个结构不同的图上。

第二代：ChebNet

切比雪夫网络实现了：快速局部化和低复杂度

\[\mathbf{g} \star \mathbf{x} = \mathbf{x} \star \mathbf{g} = (\mathbf{U} \mathbf{g}_{\theta} \mathbf{U}^{T}) \mathbf{x}， 谱图卷积 \]

从图信号分析的角度考虑，希望这个过滤函数\(\mathbf{g}\)能有较好的局部化（只影响节点的局部邻居点）。

故可把\(\mathbf{g}\)定义成\(\mathbf{L}\)的函数\(\mathbf{g}_{\theta}(\mathbf{L})\)， 例如\(\mathbf{L}\)的多项式。(这里就有从频域转向空域的意思)

因为作用一次拉普拉斯矩阵\(\mathbf{L}\), 相当于在图上把信息扩散到1阶邻居。

图信号被这个滤波器过滤后 (拉普拉斯矩阵乘法仅与特征值相关)，得到：

\[\mathbf{y} = \mathbf{g}_{\theta} (\mathbf{L})\mathbf{x} = \mathbf{g}_{\theta} (\mathbf{U} \mathbf{\Lambda} \mathbf{U}^{T}) \mathbf{x} = \mathbf{U} \mathbf{g}_{\theta} (\mathbf{\Lambda}) \mathbf{U}^{T} \mathbf{x} \]

也就是说，可把谱域图卷积中的卷积核, 看作拉普拉斯矩阵特征值\(\mathbf{\Lambda}\)的函数。通常，可选择使用一个多项式卷积核:

\[\mathbf{g}_{\theta}(\mathbf{\Lambda}) = \sum_{k=0}^{K} \mathbf{\theta_{k}} \mathbf{\Lambda}^{k} \]

其中，参数\(\mathbf{\theta_{k}}\)是多项式的系数。通过这个定义，我们现在只需要\(K+1个\)参数（\(K远小于N\)）这大大降低了参数学习过程的复杂度。就相当于:

\[\mathbf{g}_{\theta}(\mathbf{L}) = \sum_{k=0}^{K} \mathbf{\theta_{k}} \mathbf{L}^{k} \]

因此信息最多在每个节点传播\(K\)步，即即卷积的局部化。

ChebNet进一步提出了加速方案，把 \(\mathbf{g}_{\theta}(\mathbf{\Lambda})\) 近似为\(K\)阶切比雪夫多项式的：

\[\mathbf{g}_{\theta}(\mathbf{\Lambda}) = \sum_{k=0}^{K} \theta_{k} T_{k}(\tilde{\mathbf{\Lambda}}) \]

其中，Tk是k阶切比雪夫多项式。

\[\tilde{\mathbf{\Lambda}} = 2 \mathbf{\Lambda}_n / \lambda_{max} - \mathbf{I}_n是一个对角阵 \\ 主要将特征值对角阵映射到[-1，1]区间 \\ \lambda_{max}是\mathbf{L} 最大的特征值，\theta_{k} \in \mathbb{R}^{K}为切比雪夫系数向量 \]

之所以采用切比雪夫多项式，是因为考虑到它具有很好的性质，可以循环递归求解:

\[T_{k}(\mathbf{x})=2 \mathbf{x} T_{k-1}(\mathbf{x})-T_{k-2}(\mathbf{x}) \]

\[从初始值 T_{0}(\mathbf{x})=1, T_{1}(\mathbf{x})=\mathbf{x}开始,采用递归公式，可求得k阶T_k的值 \]

为了避免特征值分解，将式（3.8）写回为L的函数:

\[\begin{aligned} \mathbf{y} =\boldsymbol{g} * \mathbf{x} & = \mathbf{U} \mathbf{g}_{\theta} (\mathbf{\Lambda}) \mathbf{U}^{T} \mathbf{x} \\ & = \mathbf{U} \left( \sum_{k=0}^{K} \theta_{k} T_{k}(\tilde{\mathbf{\Lambda}}) \right) \mathbf{U}^{T} \mathbf{x} \\ & = \sum_{k=0}^{K} \theta_{k} \left(\mathbf{U} T_{k}(\tilde{\mathbf{\Lambda}}) \mathbf{U}^{T}\right) x \\ &=\sum_{k=0}^{K} \theta_{k} T_{k}(\tilde{\mathbf{L}}) \mathbf{x} \end{aligned} \]

\[其中, \tilde{\mathbf{L}}=\frac{2}{\lambda_{\max }} \mathbf{L}-\mathbf{I}_{N}。这个式子是拉普拉斯矩阵的K次多项式。 \]

因此，它仍然保持\(K\)-局部化(节点仅被其周围的\(K\)阶邻居节点所影响）。可以看到，ChebNet本质上已经是在用多阶次的\(\mathbf{L}\)对图信号进行处理。

ChebNet要学习的参数就是切比雪夫多项式的权重系数，同时还需确定切比雪夫多项式的阶数\(K\)。

第三代：GCN

GCN进一步对ChebNet进行了局部化来限制卷积操作的范围，从而来减缓过拟合的问题。

具体地，它将切比雪夫多项式的项数设为\(K=1\)，它还近似了\(\lambda_{\max } \approx 2\)，最后简化的方程如下：

\[ \mathbf{g}_{\theta^{\prime}} \star \mathbf{x} \approx \theta_{0}^{\prime} \mathbf{x}+\theta_{1}^{\prime}\left(\mathbf{L}-\mathbf{I}_{N}\right) \mathbf{x}=\theta_{0}^{\prime} \mathbf{x}-\theta_{1}^{\prime} \mathbf{D}^{-\frac{1}{2}} \mathbf{A} \mathbf{D}^{-\frac{1}{2}} \mathbf{x} \]

使用两个无限制的参数\(\theta'_0\)和\(\theta'_1\)。

在通过设置\(\theta=\theta_{0}^{\prime}=-\theta_{1}^{\prime}\)来限制参数的数量之后，可得到以下表达式：

\[ \mathbf{y} =\boldsymbol{g} * \mathbf{x} = \mathbf{g}_{\theta} \star \mathbf{x} \approx \theta\left(\mathbf{I}_{N}+\mathbf{D}^{-\frac{1}{2}} \mathbf{A} \mathbf{D}^{-\frac{1}{2}}\right) \mathbf{x} \]

值得一提的是，叠加使用这个操作会导致数值不稳定性以及梯度爆炸或消失(因为不断地乘以同一个矩阵)。因此，该论文里面使用了重规范化操作(renormalization)：

\[ \mathbf{I}_{N}+\mathbf{D}^{-\frac{1}{2}} \mathbf{A} \mathbf{D}^{-\frac{1}{2}}\rightarrow{}\tilde{\mathbf{D}}^{-\frac{1}{2}} \tilde{\mathbf{A}} \tilde{\mathbf{D}}^{-\frac{1}{2}} \]

其中， 自环邻接矩阵\(\tilde{\mathbf{A}}=\mathbf{A}+\mathbf{I}_{N}\)，\(\tilde{\mathbf{D}}_{i i}=\sum_{j} \tilde{\mathbf{A}}_{i j}\)。

然后，论文将模型扩展为含有\(C\)个输入通道的信号，\(\mathbf{X} \in \mathbb{R}^{N \times C}\)以及\(F\)个滤波器来用于提取特征：

\[ \mathbf{Z}=\tilde{\mathbf{D}}^{-\frac{1}{2}} \tilde{\mathbf{A}} \tilde{\mathbf{D}}^{-\frac{1}{2}} \mathbf{X} \boldsymbol{\Theta} \]

其中，\(\Theta \in \mathbb{R}^{C \times F}\)是滤波器参数矩阵，\(\mathbf{Z} \in\mathbb{R}^{N \times F}\)是卷积信号矩阵。

此时，GCN已经可以很好的联系到空域。GCN可以看成一层网络，而每一个GCN层相当于对目标节点的一阶邻居即自己做加权求和。具体推导可看之前博客： Pytorch-geometric: Creating Message Passing Networks 构建消息传递网络教程 中第三节关于GCN空域解释。

至此，单层的GCN变成了一个一阶模型，它每次卷积只能处理图上的1阶邻居信息。若要处理K阶邻居，需通过堆叠\(K\)个上述GCN层，来扩大图卷积地感受野。

实际上后续地改进，通常也是从这个角度出发，如果将$ \tilde{\mathbf{D}}^{-\frac{1}{2}} \tilde{\mathbf{A}} \tilde{\mathbf{D}}^{-\frac{1}{2}} $ 视为加了自环的转移矩阵\(\mathbf{P}\) ，那么完全可以预先计算其幂\(\mathbf{P}^{k}\)或邻接矩阵的幂\(\mathbf{A}^{k}\)来直接获取高阶邻居信息。