简介:这份资源面向环境科学专业学生、水务工程技术人员及相关领域研究人员,提供一套融合图卷积网络与长短期记忆网络的区域级多井地下水位时空预测方案,用于解决多井空间分布差异大、水位关联性强条件下的同步预测难题。压缩包内共1个PDF文件,约3.47MB,完整呈现GCN-LSTM模型的建模思路与实验过程。内容围绕成都城区56口观测井五年历史数据展开,详细给出空间自相似矩阵与属性自相似矩阵的构建方式、邻接矩阵组织图结构的流程,以及编码器-解码器架构下GCN提取空间关联、LSTM挖掘时间依赖的实现细节,并对比仅考虑水位特征与引入温度、降雨等气象因素后的预测精度差异。读者可据此理解时空联合建模的完整技术路线,掌握多井同步预测的验证方法与结果分析思路,为城市水务管理、防灾减灾辅助决策等场景提供可复现的参考。目前已有436人学习。
1. 区域级地下水位预测:为什么单井 LSTM 一定会翻车
如果你手上有几十上百口监测井,想用 LSTM 做地下水位预测,大概率会遇到一个很尴尬的局面:单井模型在训练集上拟合得漂漂亮亮,一到枯水期或周边新增取水井,预测曲线就开始离谱。原因不复杂——地下水位从来不是一口井自己的事,相邻井之间存在明显的水力联系,抽水、补给、隔水层分布都会让水位在空间上互相牵制。只喂时间序列、不喂空间关系,模型学到的只是"这口井过去怎么波动",而不是"这片区域的水往哪走"。
GCN-LSTM 就是冲着这个痛点来的:用图卷积神经网络(GCN)在每一时刻聚合相邻井的水位信息,把空间依赖压进节点特征;再用 LSTM 沿时间轴建模,捕捉滞后响应和季节性。两者串起来,输入是"图结构 + 多井时间序列",输出是区域内每口井未来若干天的水位。这套方案适合做区域尺度地下水资源评估、矿区/灌区水位预警、以及需要跨井推理的场合。下面按我实际搭过的一版流程,从建图、造数据、写模型到排查坑,一步步拆开讲。
2. 把监测井变成图:GCN-LSTM 的输入到底长什么样
2.1 图结构怎么建,邻接矩阵不是随便连
GCN 的核心输入是邻接矩阵 A 和节点特征矩阵 X。地下水位场景里,节点就是监测井,边代表井间的水力联系。常见做法有三种:一是按地理距离阈值连边,比如两井直线距离小于 2 km 就连;二是用反距离权重,把 1/d 作为边权;三是结合水文地质单元,同一含水层组内才连边。我一般先用距离阈值快速跑通,再根据残差图微调。
import numpy as np from scipy.spatial.distance import cdist # wells_xy: (N, 2) 每口井的经纬度或投影坐标 # dist_threshold: 连边距离阈值,单位与坐标一致 def build_adjacency(wells_xy, dist_threshold=2000): dist = cdist(wells_xy, wells_xy, metric='euclidean') # 距离小于阈值且非自身,置为反距离权重 A = np.where((dist < dist_threshold) & (dist > 0), 1.0 / (dist + 1e-6), 0.0) # 对称化,避免有向图带来的信息不对称 A = np.maximum(A, A.T) # 行归一化,防止度数大的节点主导聚合 deg = A.sum(axis=1, keepdims=True) deg[deg == 0] = 1.0 A_norm = A / deg return A_norm A = build_adjacency(wells_xy, dist_threshold=2000) print(A.shape) # (N, N)这段代码做了三件事:算距离、按阈值连边并赋反距离权重、行归一化。参数dist_threshold是最需要调的,设太小图会碎成孤岛,GCN 退化成单井全连接;设太大所有井互相连,空间信息被平均掉。我的经验是先看井距分布的中位数,取 1.5 到 2 倍中位数试。归一化那步别省,否则高度数节点在聚合时会淹没邻居。
2.2 时间滑窗与特征工程:LSTM 吃的是什么
建完图,接下来把多井水位序列切成监督学习样本。假设有 N 口井、T 个时间步,用过去 L 步预测未来 H 步,输入张量形状是(样本数, L, N, 特征数)。特征除了水位本身,通常还会加降雨、气温、开采量。注意:GCN 是在每个时间步上对节点做空间聚合,所以数据要按时间步组织好。
def make_sliding_windows(series, lookback=30, horizon=7): # series: (T, N, F) 时间, 井, 特征 X, Y = [], [] T = series.shape[0] for t in range(T - lookback - horizon + 1): X.append(series[t:t+lookback]) # (L, N, F) Y.append(series[t+lookback:t+lookback+horizon, :, 0]) # 只预测水位 return np.array(X), np.array(Y) X, Y = make_sliding_windows(series, lookback=30, horizon=7) print(X.shape, Y.shape) # (样本, 30, N, F) (样本, 7, N)lookback=30对应约一个月的历史窗口,horizon=7是未来一周。地下水位滞后性强,窗口太短学不到补给响应,太长则引入噪声且显存吃紧。我一般从 30 天起步,枯水期可拉到 60。标签只取水位通道,是因为降雨等外生变量未来值往往拿不到,强行预测会引入泄漏。
2.3 模型结构:GCN 和 LSTM 谁先谁后
有两种串法:一种是每个时间步先过 GCN 做空间聚合,再把聚合后的序列喂给 LSTM;另一种是先用 LSTM 抽时间特征,再对每个时间步做 GCN。前者更符合"先空间后时间"的物理直觉,也是我推荐的默认结构。核心代码如下:
import torch import torch.nn as nn class GCNLayer(nn.Module): def __init__(self, in_dim, out_dim): super().__init__() self.linear = nn.Linear(in_dim, out_dim) def forward(self, x, A): # x: (B, N, F), A: (N, N) h = self.linear(x) # (B, N, out_dim) h = torch.einsum('nn,bnf->bnf', A, h) # 邻居聚合 return torch.relu(h) class GCNLSTM(nn.Module): def __init__(self, feat_dim, gcn_hidden=32, lstm_hidden=64, horizon=7): super().__init__() self.gcn = GCNLayer(feat_dim, gcn_hidden) self.lstm = nn.LSTM(gcn_hidden, lstm_hidden, batch_first=True) self.head = nn.Linear(lstm_hidden, horizon) def forward(self, x, A): # x: (B, L, N, F) B, L, N, F = x.shape x = x.reshape(B * L, N, F) h = self.gcn(x, A) # (B*L, N, gcn_hidden) h = h.reshape(B, L, N, -1) # 对每口井独立跑 LSTM,共享权重 h = h.permute(0, 2, 1, 3).reshape(B * N, L, -1) out, _ = self.lstm(h) out = out[:, -1, :] # 取最后时间步 out = self.head(out) # (B*N, horizon) return out.reshape(B, N, -1)einsum('nn,bnf->bnf', A, h)就是图卷积的聚合操作,A 是归一化邻接矩阵。LSTM 部分把井维度和批次维度合并,实现权重共享——每口井用同一套时序参数,但输入已经带上了邻居信息。gcn_hidden和lstm_hidden是主要容量参数,井数少可以调小防过拟合,井数上百则适当加大。输出层对每口井独立预测 horizon 步。
3. 训练与评估:从数据划分到指标怎么看
3.1 数据划分不能随机打乱
时间序列最忌讳随机划分。正确做法是按时间切:前 70% 训练,中间 15% 验证,最后 15% 测试。如果用随机划分,未来信息会泄漏到训练集,指标虚高,上线必崩。另外,标准化参数只能用训练集统计量,再应用到验证和测试。
def temporal_split(X, Y, ratios=(0.7, 0.15, 0.15)): n = len(X) n_train = int(n * ratios[0]) n_val = int(n * ratios[1]) return (X[:n_train], Y[:n_train], X[n_train:n_train+n_val], Y[n_train:n_train+n_val], X[n_train+n_val:], Y[n_train+n_val:]) Xtr, Ytr, Xva, Yva, Xte, Yte = temporal_split(X, Y) # 标准化:只用训练集 mean, std = Xtr.mean(axis=(0,1,2), keepdims=True), Xtr.std(axis=(0,1,2), keepdims=True) Xtr, Xva, Xte = (Xtr-mean)/std, (Xva-mean)/std, (Xte-mean)/std标准化维度要对齐(1,1,N,F),即每口井每个特征单独标准化,因为不同井的水位量级可能差几十米。用全局均值会把浅井和深井混在一起。
3.2 损失函数与评估指标
回归任务默认 MSE,但地下水位预测里我更关注峰值和枯水期,所以常用组合损失:MSE 加一项对高水位的加权。评估指标除了 RMSE、MAE,一定要看 NSE(纳什效率系数),它对水文序列的拟合优度更敏感。
def nse(y_true, y_pred): # y_true, y_pred: (样本, horizon, N) 或展平 y_true = y_true.reshape(-1) y_pred = y_pred.reshape(-1) return 1 - np.sum((y_true-y_pred)**2) / np.sum((y_true-y_true.mean())**2) # 训练循环关键片段 criterion = nn.MSELoss() optimizer = torch.optim.Adam(model.parameters(), lr=1e-3) for epoch in range(200): model.train() pred = model(torch.FloatTensor(Xtr), torch.FloatTensor(A)) loss = criterion(pred, torch.FloatTensor(Ytr)) optimizer.zero_grad(); loss.backward(); optimizer.step()NSE 大于 0.75 通常算可用,大于 0.9 算优秀。如果测试集 NSE 远低于训练集,先查是不是图结构把不相关的井连在了一起,或者标准化用了全局统计量。
3.3 一个可跑通的最小训练配置
把上面拼起来,给一组我常用的起步超参:lookback=30、horizon=7、gcn_hidden=32、lstm_hidden=64、lr=1e-3、batch_size=32、epochs=200,优化器 Adam,早停耐心 20 轮。井数在 50 到 200 之间时这套配置基本能跑出合理结果。显存不够就减 batch_size 或 gcn_hidden,别一上来就堆大模型——地下水位数据量通常撑不起。
4. 避坑与排查:这五个坑我基本都踩过
4.1 现象:验证集 loss 正常,测试集 NSE 为负
原因:时间划分时标准化用了全量数据,或者滑窗跨越了训练/测试边界,未来信息泄漏。解决:严格按时间切分,标准化统计量只从训练集算,滑窗生成后再切分,不要先切分再滑窗。
4.2 现象:模型对所有井输出几乎相同的值
原因:邻接矩阵归一化后过于平滑,GCN 把每口井的特征都平均成了区域均值。解决:检查dist_threshold是否过大,改用反距离权重并保留自环,或在 GCN 后加残差连接保留节点自身信息。
4.3 现象:训练 loss 下降但预测曲线滞后一步
原因:LSTM 学到了"复制上一时刻"的捷径,尤其在水位变化平缓时。解决:在损失里加大变化量的权重,或把输入做一阶差分,让模型学增量而非绝对值。
4.4 现象:新增监测井后模型无法推理
原因:GCN 的邻接矩阵维度固定,井数一变整个图就变了。解决:预留可扩展的图构建流程,新井按同样规则连边后重新归一化;或者训练时用掩码,让模型见过不同规模的子图。
4.5 现象:枯水期预测系统性偏高
原因:训练集里枯水期样本少,模型偏向丰水期分布。解决:对枯水期样本过采样,或在损失里按季节加权,也可以引入降雨滞后项作为外生输入。
5. 进阶技巧:用残差图反哺图结构,让 GCN 自己学会连边
跑通基础版之后,最值得投入的一步是让图结构从数据里学出来,而不是全靠距离阈值拍脑袋。我的做法是:先用基础模型跑一遍,计算每口井的预测残差,残差相关性高的井对,说明它们之间存在模型没捕捉到的共同驱动,就在下一轮把这条边加上或加大权重。这相当于用残差图做一次图结构微调,通常能把测试集 NSE 再抬 0.03 到 0.08。
# 残差相关性驱动的边权增强 residuals = Yte - pred_te # (样本, horizon, N) res_corr = np.corrcoef(residuals.reshape(-1, residuals.shape[-1]).T) # (N, N) # 残差正相关的井对,增强其边权 A_enhanced = A + 0.3 * np.where(res_corr > 0.5, res_corr, 0) np.fill_diagonal(A_enhanced, 0) deg = A_enhanced.sum(axis=1, keepdims=True); deg[deg==0]=1 A_enhanced = A_enhanced / deg0.3是增强系数,0.5是相关性阈值,这两个参数按区域调整。残差相关性高不一定代表水力连通,也可能是共同受降雨影响,所以增强后要用验证集确认,别直接上测试集调。
验证方法上,我习惯留出最近一个完整水文年做滚动预测:每次预测未来 7 天,然后窗口前移 7 天,累计看全年 NSE 和峰值误差。这比单次划分更接近实际使用。另外,把预测结果按井位画成空间热力图,能直观看出模型在哪些区域系统性偏差,往往比看数字更快定位问题。
最后说个习惯:我每次改图结构或损失函数,都会固定随机种子跑三遍取均值,单次结果好看不算数。地下水位数据噪声大,一次跑赢可能是运气。这套 GCN-LSTM 方案值不值得做,取决于你是否有跨井的空间依赖需要建模——如果井距很远、水力联系弱,老老实实单井 LSTM 加外生变量可能更省事。希望帮到你。
本文还有配套的精品资源,点击获取