本文介绍了图神经网络（GNN）的数学基础，特别是GCN的实现，以及结合注意力机制、双向LSTM和编码器-解码器的多指标水质预测方法。作者还探讨了如何使用谱图理论和傅里叶变换进行图卷积操作，以及使用DGL库实现GCN进行图节点分类。
目录
摘要
本周结束了GNN的学习，通过对GCN进行数学推导，进一步加深了对GCN的理解，明白了图神经网络中运用卷积思想的数学意义，同时编写了一个简单的GCN模型对图节点进行预测。本周还阅读了一篇论文，学习到了一些时间序列数据预处理的方法，以及一些超参数优化算法。
ABSTRACT
This week, we have finished our study of GNN. Through mathematical derivations of GCN, we deepened our understanding of GCN, comprehending the mathematical significance of applying convolutional thinking in graph neural networks. Simultaneously, we developed a simple GCN model to predict graph nodes. Additionally, we read a paper this week and learned some preprocessing methods for time series data, as well as various parameter optimization algorithms beyond the BP algorithm.
1 文献阅读
摘要: 本文提出了一种基于注意力（Attention）辅助的双向LSTM（BiLSTM）和编码器-解码器（ED）的多指标水质预测方法，通过结合Savitzky-Golay滤波器（SG）、变分模态分解（VMD）、遗传模拟退火粒子群优化（GSPSO）等技术，实现了对水质时间序列的准确预测。
论文背景: 随着城市的发展和城市化的加速，城市水的健康循环是城市健康发展的重要基础，也是维护良好城市水环境的基本条件和保持健康城市生态环境的必要前提。水质的质量与人们的生命安全和财产安全密切相关，因此建立准确的水质预测系统以提高水质至关重要。水质预测是一个时间序列预测问题，通过高频水质传感器收集各种水质监测指标的时间序列数据，可以实时获取各种水质监测指标。
过去方案: 目前，时间序列预测方法主要分为传统统计方法和深度学习方法。传统的统计方法提取数据的线性关系，并开发了许多优秀的模型。例如，支持向量机（SVM）由于速度快，计算量少，被广泛用于时间序列预测。然而，很难应对大规模的训练样本。自回归积分移动平均线（ARIMA）结合了线性和非线性模型的优点，并提取了时间序列中的非线性关系。然而，水质变化受多种因素影响，其时间序列复杂且非线性。因此，仅靠传统的统计方法无法捕捉水质的细微变化。新兴的深度学习技术被广泛用于解决时间序列预测问题。递归神经网络（RNN）捕获长期相关性。作为一种典型的变体，长短期记忆（LSTM）解决了长序列训练过程中梯度消失和梯度爆炸的问题。然而，它无法得出给定结果的原因。人工神经网络模型也常用于时间序列预测，但它需要大量的参数和较长的学习时间。
论文方案: 该论文提出了一种混合模型SVABEG，该模型结合了SG滤波器、VMD、Attention机制、BiLSTM、ED结构和GSPSO。SVABEG首先采用SG滤波器和VMD分别去除原始时间序列中的噪声和非线性特征。然后，SVABEG结合BiLSTM、ED结构和Attention机制，分别捕获双向长期相关性，实现降维和提取关键信息。最后，SVABEG采用GSPSO来优化其超参数。
SG滤波器: SG滤波器是一种滤波技术，它采用时域中多项式最小二乘法的拟合。它可以平滑地对时间序列信号进行降噪，同时保留其原始特征。具体来说，它通过卷积过程拟合相邻数据点的连续子集。SG 滤波器有两个极其重要的参数，即滤波器窗口大小和多项式拟合阶数。窗口大小受卷积操作中的卷积系数数量的影响。实现代码与应用效果如下：
import numpy as np
from scipy.signal import savgol_filter
import matplotlib.pyplot as plt
# 创建一些虚拟的数据
x = np.linspace(0, 2 * np.pi, 100)
y = np.sin(x) + np.random.random(100) * 0.2 # 加入一些噪声
# 应用SG滤波器
yhat = savgol_filter(y, 15, 3) # 使用窗口大小为15和多项式阶数为3的SG滤波器
# 绘制原始数据和滤波后的数据
plt.plot(x, y, label='Noisy signal')
plt.plot(x, yhat, color='red', label='Filtered signal')
plt.legend()
plt.show()
VMD: VMD方法可以对原始时间序列进行分解。其采用非递归处理策略将原始复序列分解为多个相对稳定的子序列，从而减少噪声干扰，便于进一步预测。作为一种自适应时频分析方法，VMD广泛用于非线性和非平稳信号分解。它结合了维纳滤波算法、希尔伯特变换和频率混叠。它可以减少非线性和噪声时间序列的波动，并避免多模态混合问题。与其他信号分解方法相比，其分解模态分量的数量可以动态调整。这些分量稀疏且围绕中心频率波动很大，通过最小化每个模态分量的带宽总和来实现信号稳定性。实现代码与应用效果如下：
import numpy as np
import matplotlib.pyplot as plt
from vmdpy import VMD
# 测试信号及其参数
Fs = 1000 # 采样频率
N = 1000 # 采样点数
t = np.arange(1, N + 1) / N
fre_axis = np.linspace(0, Fs / 2, int(N / 2))
f_1 = 100
f_2 = 200
f_3 = 300
v_1 = (np.cos(2 * np.pi * f_1 * t))
v_2 = 1 / 4 * (np.cos(2 * np.pi * f_2 * t))
v_3 = 1 / 16 * (np.cos(2 * np.pi * f_3 * t))
v = [v_1, v_2, v_3] # 测试信号所包含的各成分
f = v_1 + v_2 + v_3 + 0.1 * np.random.randn(v_1.size) # 测试信号
# alpha惩罚系数，带宽限制经验取值为抽样点长度1.5-2.0倍
# 惩罚系数越小，各IMF分量的带宽越大，过大的带宽会使得某些分量包含其他分量言号
# a值越大，各IMF分量的带宽越小，过小的带宽是使得被分解的信号中某些信号丢失该系数常见取值范围为1000~3000
alpha = 2000
tau = 0 # tau 噪声容限，即允许重构后的信号与原始信号有差别。
K = 3 # K 分解模态（IMF）个数
DC = 0 # DC 若为0则让第一个IMF为直流分量/趋势向量
init = 1 # init 指每个IMF的中心频率进行初始化。当初始化为1时，进行均匀初始化。
tol = 1e-7 # 控制误差大小常量，决定精度与迭代次数
u, u_hat, omega = VMD(f, alpha, tau, K, DC, init, tol) # 输出U是各个IMF分量，u_hat是各IMF的频谱，omega为各IMF的中心频率
# 1 画原始信号和它的各成分
plt.figure(figsize=(10, 7))
plt.subplot(K + 1, 1, 1)
plt.plot(t, f)
for i, y in enumerate(v):
plt.subplot(K + 1, 1, i + 2)
plt.plot(t, y)
plt.suptitle('Original input signal and its components')
plt.show()
# 2 分解出来的各IMF分量
plt.figure(figsize=(10, 7))
plt.plot(t, u.T)
plt.title('all decomposed modes')
plt.show() # u.T是对u的转置
# 3 各IMF分量的fft幅频图
plt.figure(figsize=(10, 7), dpi=80)
for i in range(K):
plt.subplot(K, 1, i + 1)
fft_res = np.fft.fft(u[i, :])
plt.plot(fre_axis, abs(fft_res[:int(N / 2)]) / (N / 2))
plt.title('(FFT) amplitude frequency of IMF {}'.format(i + 1))
plt.show()
# 4 分解出来的各IMF分量的频谱
# print(u_hat.shape,t.shape,omega.shape)
plt.figure(figsize=(10, 7), dpi=80)
for i in range(K):
plt.subplot(K, 1, i + 1)
plt.plot(fre_axis, abs(u_hat[:, i][int(N / 2):]) / (N / 2))
plt.title('(VMD) amplitude frequency of the modes{}'.format(i + 1))
plt.tight_layout()
plt.show()
# 5 各IMF的中心频率
plt.figure(figsize=(12, 7), dpi=80)
for i in range(K):
plt.subplot(K, 1, i + 1)
plt.plot(omega[:, i]) # X轴为迭代次数，y轴为中心频率
plt.title('mode center-frequencies{}'.format(i + 1))
plt.tight_layout()
plt.show()
plt.figure(figsize=(10, 7))
plt.plot(t, np.sum(u, axis=0))
plt.title('reconstructed signal')
ED: ED是实现时间序列预测的常用框架，由编码器和解码器组成。编码器将输入序列视为具有语义的向量，解码器采用向量作为输入对目标序列进行解码。对于输入向量<X,Y>，我们的目标是通过ED框架预测目标Y，这是给定X输入序列的预测值。,
BiLSTM: BiLSTM不仅可以对自然语言处理任务中的上下文信息进行建模，还可以有效解决LSTM的反向编码问题（LSTM无法捕捉到输入信息中的反向关系，导致实际预测中缺少一些关键信息）。它可以捕获双向语义依赖关系，并在更细粒度的分类任务上实现更好的性能。
注意力机制: 当输入时间序列很长时，很难学习合理的向量表示。注意力机制可以打破ED框架在编解码中依赖于固定长度向量的限制。它使神经网络能够将有限的注意力集中在序列的重要信息上。通过这种方式，它提取了时间序列中最具代表性的特征，而较少关注序列中无关紧要的信息。通过在ED框架中加入注意力机制，可以对时间序列数据进行加权和转换，以提高系统性能。
GSPSO: 模型中包含大量的超参数，它们包括层数、批量大小、学习率、辍学率、权重衰减、时期数和序列长度等。为了使预测精度更高，大量研究证明，使用PSO调整超参数可以提高预测精度。本文选择了一种名为GSPSO的PSO改进版本（结合遗传算法、模拟退火算法与粒子群优化算法）来优化超参数。
2 Math of GCN
本节将探究图卷积神经网络（GCN）的数学原理，搞清楚GNN如何从数学的角度出发表示卷积操作的过程。
2.1 Basic of GCN
用于处理图结构数据的图神经网络（GNN）从本质上来说，就是给定图中每一个节点的特征以及图的邻接矩阵，通过输入这两种数据从而映射出下一层节点的特征。从节点的角度上来说，节点下一层的特征受到了节点上一层自身特征和周边邻居节点特征的影响。因此所有GNN的模型都是在构造映射的方式，例如GCN就代表了其中的一种。
2.2 Spectral Graph Theory
研究GCN免不了涉及到谱图理论，首先需要先了解两个与邻接矩阵相关的矩阵的特性。这两个矩阵，一个是拉普拉斯矩阵，另一个是拉普拉斯矩阵对应的对称规范化矩阵。之所以要研究这两个矩阵，是因为这两个矩阵都是实对称阵，都具有n个特征值和n个正交的特征向量，同时这两个矩阵都是半正定的，即它们的特征值都大于等于0。下图解释了为什么这两个矩阵是半正定的。
拉普拉斯矩阵对应的对称规范化矩阵还有一个最重要的特性是，它的特征值的范围是从0到2的（包括0和2）。证明如下：
2.3 Fourier Transformation
GCN也涉及到傅里叶变换的相关知识。什么是傅里叶变换？简单来说就是研究同一个事物在不同域（时域，频域）中的视角是什么样的，以及在不同域之间进行变换。通过傅里叶可以知道，任何周期函数，都可以看作是不同振幅，不同相位正弦波的叠加。为什么要用到傅里叶变换？因为傅里叶变换是无损的，而且有些操作在频域视角下更容易。比如说一段男女同时说话的音频中，从时域的视角下很难进行区分，但在频域的视角下，男生的音调频率通常较低，女生的音调频率通常较高，所以可以比较容易进行区分。将男生的声音去除后得到的女生声音再经过傅里叶逆变换就能得到女生的音频数据。下图是关于快速傅里叶变换 (Fast Fourier Transform)的例子。
很明显，图像可以直接使用卷积操作是因为图像具有比较规则的拓扑结构，但是图结构数据是有任意复杂的拓扑结构的，无法在空间域中给它一个特定形状的卷积核，因为每个节点可能有不同数量的邻居。根据傅里叶的启发，可以将图结构数据进行域的变换之后再进行卷积操作，卷积完成后再通过逆变换将域变回来。再回到前面所述的拉普拉斯矩阵，将其与一个列特征向量相乘，可以发现，拉普拉斯矩阵使得特征向量进行了一个类似聚合邻居的操作，每个特征都变成了与其邻居的一些交互，这种邻居的聚合其实就是卷积的本质。
理论上特征值分解的操作能够实现卷积操作，但其时间复杂度太高（）,并不能直接运用于实际。GCN对特征值分解的傅里叶变换做了一些限制，推导出了一种不需要特征值分解，复杂度与边的数量呈线性关系的一种卷积方法。
2.4 GCN
GCN的推导过程如下：
3 实现代码
基于DGL实现图节点的分类：
import dgl
import torch
import torch.nn as nn
import torch.nn.functional as F
import dgl.function as fn
# 创建一个简单图
G = dgl.DGLGraph()
G.add_nodes(5)
G.add_edges([0, 1, 1, 2, 2, 3, 3, 4], [1, 0, 2, 1, 3, 2, 4, 3])
# 定义一个GCN层
class GCNLayer(nn.Module):
# in_feats, out_feats 分别指定每个节点的输入特征和输出特征的维度
def __init__(self, in_feats, out_feats):
super(GCNLayer, self).__init__()
self.linear = nn.Linear(in_feats, out_feats)
def forward(self, g, feature):
# 创建一个局部作用域 with g.local_scope()，这是 DGL 提供的上下文管理器
# 用于确保在这个作用域内所做的修改不会影响全局图结构
with g.local_scope():
# 这一行将输入的节点特征 feature 存储到图 g 的节点数据中
# 'h' 是节点特征的名称，可以根据需要选择任何名称
g.ndata['h'] = feature
# 使用 update_all 函数执行消息传递过程
# 在消息传递过程中，首先从源节点 (src) 复制节点特征 'h'，并将其输出为 'm'，然后将这些消息进行汇总，计算邻居节点特征的和，并将其输出为 'h_neigh'
# 这个过程实现了 GCN 层的消息传递操作
g.update_all(fn.copy_u(u='h', out='m'),
fn.sum(msg='m', out='h_neigh'))
# 在消息传递完成后，从节点数据中提取两个重要的变量：h_neigh 表示邻居节点的特征总和，h_self 表示节点自身的特征
# 这些变量将用于计算 GCN 层的输出
h_neigh = g.ndata['h_neigh']
h_self = g.ndata['h']
# 将节点自身的特征和邻居节点的特征相加，得到了节点的总特征 h_total
# 这是 GCN 层的核心操作，它将节点的自身特征与邻居节点的信息相融合，以获得更丰富的表示
h_total = h_self + h_neigh
# 最后，将融合后的节点特征 h_total 通过线性变换 self.linear，得到最终的 GCN 层的输出
# 这个输出可以用于后续的分类或其他任务
return self.linear(h_total)
# 创建一个简单的GCN模型
class GCN(nn.Module):
def __init__(self, in_feats, hidden_size, num_classes):
super(GCN, self).__init__()
self.gcn1 = GCNLayer(in_feats, hidden_size)
self.gcn2 = GCNLayer(hidden_size, num_classes)
def forward(self, g, features):
x = F.relu(self.gcn1(g, features))
x = self.gcn2(g, x)
return x
# 设置模型超参数
in_feats = 5
hidden_size = 3
num_classes = 2
# 初始化模型和数据
model = GCN(in_feats, hidden_size, num_classes)
features = torch.eye(5)
labels = torch.LongTensor([0, 1, 0, 0, 1])
# 训练模型
optimizer = torch.optim.Adam(model.parameters(), lr=0.01)
criterion = nn.CrossEntropyLoss()
for epoch in range(100):
logits = model(G, features)
mask = torch.BoolTensor([True, True, False, False, False])
loss = criterion(logits[mask], labels[mask])
optimizer.zero_grad()
loss.backward()
optimizer.step()
# 进行预测
model.eval()
logits = model(G, features)
predicted = torch.argmax(logits, 1)
print("Predicted Class:", predicted)
运行结果：
总结
GNN是一种强大的模型，对处理图结构数据具有强大的作用。下一周我将开始Transformer相关知识的学习。
扫一扫
449
到【灌水乐园】发言
