当前位置：首页 > news >正文

深度学习与遥感入门（五）｜GAT 构图消融 + 分块全图预测：更稳更快的高光谱图分类（PyTorch Geometric 实战）

news 2025/8/12 5:30:32

系列回顾：

（一）CNN 基础：高光谱图像分类可视化全流程，链接：https://mp.weixin.qq.com/s/P4IOG0WTDuoBEprfGSWMTQ
（二）HybridNet（CNN+Transformer）：提升全局感受野，链接：https://mp.weixin.qq.com/s/Zlev4Z0b3VE7a6jOOpzaAA
（三）GCN 入门实战：基于光谱 KNN 的图卷积分类与全图预测，链接：https://mp.weixin.qq.com/s/Vo5QNA7gkqbYhYg10krbnQ
（四）空间–光谱联合构图的 GCN：RBF 边权 + 自环 + 早停，得到更稳更自然的全图分类结果，链接：https://mp.weixin.qq.com/s/G7VnMzhby4Fvmjwh_R_z7A
本篇（五）：在（四）的基础上，加入 GAT（注意力图网络）与构图消融，并实现分块全图预测防 OOM，确保在显存有限的设备上也能稳定跑完全图。

一、这篇要解决什么问题？

构图到底该怎么选？
只用光谱？只用空间？还是融合？我们做一键消融：pure_spectral / pure_spatial / fused。
GAT vs. GCN？
GCN需要显式边权（连续相似度）；GAT让模型自己学习“看谁更重要”。本篇支持一键切换：MODEL='GCN' | 'GAT'。
全图预测容易爆显存？
一次性为整幅图构 KNN 图很吃内存。本篇提供分块局部构图 + 内核拼接，用 BLOCK_SIZE / OVERLAP 控制，显著降低内存占用。

二、这篇都做了哪些事情？

构图消融：光谱KNN、空间KNN、空间–光谱融合（RBF边权 + 对称化 + 无向图 + 自环）
模型双模：GCN（显式传 edge_weight） / GAT（连通性图即可）
早停 + 固定种子：可复现且节省时间
分块全图预测：局部构图与推理，仅写回“内核区域”，减少边界效应与显存压力
可视化：混淆矩阵、各类别准确率条形图、全图分类图

三、方法解释

3.1 构图（含消融）

光谱：在 PCA 特征空间做 KNN，距离 → 中位数归一 → RBF 相似
空间：在像素坐标平面做 KNN，距离 → 中位数归一 → RBF 相似
融合：S = α·S_spec + (1-α)·S_spat，然后对称化并转无向图，最后加自环
GCN：使用连续边权；GAT：使用连通性图（不需要权重）

3.2 分块全图预测

将图像划成重叠块（BLOCK_SIZE，OVERLAP），对每个块单独构图并推理
仅把去掉 overlap 的内核区域写回全图，保证拼接平滑
极端未覆盖像素兜底为 0 类（可改成最近邻填补）

3.3 结果展示

在这里插入图片描述

四、完整可运行代码（含详细注释）

直接复制运行。按需修改 DATA_DIR / X_FILE / Y_FILE。
默认数据：KSC。依赖：torch, torch_geometric, scikit-learn, matplotlib, seaborn, scipy.

# -*- coding: utf-8 -*-
"""
Series (V): GAT & Graph Construction Ablation for HSI Classification (PyTorch Geometric)
- 解决内存溢出问题（全图分块处理：局部构图 + 内核拼接）
- 邻接矩阵强制对称化、无向化；GCN用连续RBF边权，GAT用连通性
- 完善早停逻辑与复现性
- 构图消融：pure_spectral / pure_spatial / fused
- 可视化：混淆矩阵、各类别准确率、全图预测
"""import os
import numpy as np
import scipy.io
import torch
import torch.nn as nn
import torch.nn.functional as F
from torch_geometric.data import Data
from torch_geometric.nn import GCNConv, GATConv
from torch_geometric.utils import add_self_loops, to_undirectedfrom sklearn.preprocessing import StandardScaler
from sklearn.decomposition import PCA
from sklearn.metrics import accuracy_score, classification_report, confusion_matrix
from sklearn.neighbors import kneighbors_graph
from sklearn.model_selection import train_test_splitimport matplotlib
import matplotlib.pyplot as plt
import seaborn as sns# ----------------- 全局样式 -----------------
matplotlib.rcParams['font.family'] = 'SimHei'
matplotlib.rcParams['axes.unicode_minus'] = False
plt.rcParams['figure.dpi'] = 120
sns.set_theme(context="notebook", style="whitegrid", font="SimHei")# ----------------- 可配置超参（可一键消融&切换） -----------------
DATA_DIR = r"your_path"
X_FILE = "KSC.mat"
Y_FILE = "KSC_gt.mat"# PCA
PCA_DIM = 15                 # 若设为0，则自动选择解释方差≥PCA_VAR_THRESHOLD的维度
PCA_VAR_THRESHOLD = 0.95# 构图消融：'pure_spectral' | 'pure_spatial' | 'fused'
GRAPH_MODE = 'fused'# 模型：'GCN' | 'GAT'
MODEL = 'GAT'# KNN与RBF
K_SPEC = 6
K_SPAT = 6
ALPHA   = 0.7
SIGMA_SPEC = 1.0
SIGMA_SPAT = 1.0
SELF_LOOP_VALUE = 1.0# 训练
HIDDEN = 64
DROPOUT = 0.5
LR = 0.01
WD = 5e-4
MAX_EPOCHS = 200
PATIENCE = 15
TRAIN_RATIO = 0.3
SEED = 42# GAT参数
GAT_HEADS = 4
GAT_CONCAT = True# 分块全图预测
BLOCK_SIZE = 128
OVERLAP = 32# ----------------- 复现性 -----------------
def set_seeds(seed=42):import randomrandom.seed(seed)np.random.seed(seed)torch.manual_seed(seed)torch.cuda.manual_seed_all(seed)torch.backends.cudnn.deterministic = Truetorch.backends.cudnn.benchmark = Falseset_seeds(SEED)
device = torch.device('cuda' if torch.cuda.is_available() else 'cpu')
print(f"当前设备：{device}")# =========================
# 数据加载与预处理
# =========================
def load_data(x_path, y_path):"""加载高光谱数据和标签，自动识别键名"""x_data = scipy.io.loadmat(x_path)y_data = scipy.io.loadmat(y_path)x_key = [k for k in x_data.keys() if not k.startswith('__')][0]y_key = [k for k in y_data.keys() if not k.startswith('__')][0]return x_data[x_key], y_data[y_key]print("加载数据...")
X_image, y_image = load_data(os.path.join(DATA_DIR, X_FILE), os.path.join(DATA_DIR, Y_FILE))
h, w, bands = X_image.shape
print(f"图像尺寸: {h} x {w}, 波段数: {bands}")# 展平 & 标签处理
X_flat = X_image.reshape(-1, bands)
y_flat = y_image.reshape(-1).astype(np.int32)
mask_labeled = y_flat != 0
X_labeled = X_flat[mask_labeled]
y_labeled = (y_flat[mask_labeled] - 1).astype(np.int64)
num_classes = int(len(np.unique(y_labeled)))
print(f"有效样本: {len(y_labeled)}，类别数: {num_classes}")# 空间坐标（全图&标注）
rows_all = np.repeat(np.arange(h), w)
cols_all = np.tile(np.arange(w), h)
coords_all = np.stack([rows_all, cols_all], axis=1).astype(np.float32)
coords_labeled = coords_all[mask_labeled]# 标准化 + PCA（仅在标注样本上fit；全图transform）
scaler = StandardScaler()
X_scaled = scaler.fit_transform(X_labeled)if PCA_DIM == 0:pca_probe = PCA(n_components=None, random_state=SEED)pca_probe.fit(X_scaled)cum_var = np.cumsum(pca_probe.explained_variance_ratio_)PCA_DIM = int(np.argmax(cum_var >= PCA_VAR_THRESHOLD) + 1)print(f"自动选择PCA维度: {PCA_DIM}（解释方差: {cum_var[PCA_DIM-1]:.2f}）")pca = PCA(n_components=PCA_DIM, random_state=SEED)
X_pca = pca.fit_transform(X_scaled)# =========================
# 构图工具（对称化+无向图）
# =========================
def _normalize_dist(d):"""将KNN距离按全局中位数归一，抗异常值"""data = d.data.copy()if np.any(data > 0):med = np.median(data[data > 0])else:med = 1.0d_norm = d.copy()d_norm.data = data / (med + 1e-12)return d_normdef _rbf_sim_from_dist(d_norm, sigma):"""RBF相似度：exp(-d^2/(2*sigma^2))"""sim = d_norm.copy()sim.data = np.exp(-(sim.data ** 2) / (2.0 * (sigma ** 2) + 1e-12))return simdef build_graph(feats, coords, mode='fused',k_spec=6, k_spat=6,alpha=0.7, sig_spec=1.0, sig_spat=1.0,to_connectivity_for_gat=False):"""构建无向图（确保对称化）返回: edge_index (LongTensor[2,E]), edge_weight (FloatTensor[E]或None)"""if mode == 'pure_spectral':d = kneighbors_graph(feats, n_neighbors=k_spec, mode='distance', include_self=False)d = _normalize_dist(d)sim = _rbf_sim_from_dist(d, sig_spec)mat = (sim + sim.T) / 2elif mode == 'pure_spatial':d = kneighbors_graph(coords, n_neighbors=k_spat, mode='distance', include_self=False)d = _normalize_dist(d)sim = _rbf_sim_from_dist(d, sig_spat)mat = (sim + sim.T) / 2else:  # fusedd_spec = kneighbors_graph(feats, n_neighbors=k_spec, mode='distance', include_self=False)d_spec = _normalize_dist(d_spec)sim_spec = _rbf_sim_from_dist(d_spec, sig_spec)sim_spec = (sim_spec + sim_spec.T) / 2d_spat = kneighbors_graph(coords, n_neighbors=k_spat, mode='distance', include_self=False)d_spat = _normalize_dist(d_spat)sim_spat = _rbf_sim_from_dist(d_spat, sig_spat)sim_spat = (sim_spat + sim_spat.T) / 2mat = sim_spec * alpha + sim_spat * (1.0 - alpha)row, col = mat.nonzero()weight = mat[row, col]if hasattr(weight, "A1"):weight = weight.A1edge_index = np.vstack([row, col]).astype(np.int64)edge_weight = weight.astype(np.float32)# 转无向（双向边）edge_index = to_undirected(torch.tensor(edge_index, dtype=torch.long))if to_connectivity_for_gat:edge_weight = Noneelse:edge_weight = torch.tensor(edge_weight, dtype=torch.float32)return edge_index, edge_weight# =========================
# 构建训练图（仅标注像素）
# =========================
print(f"构建训练图：mode={GRAPH_MODE}, model={MODEL}")
edge_index_train, edge_weight_train = build_graph(X_pca, coords_labeled, mode=GRAPH_MODE,k_spec=K_SPEC, k_spat=K_SPAT, alpha=ALPHA,sig_spec=SIGMA_SPEC, sig_spat=SIGMA_SPAT,to_connectivity_for_gat=(MODEL == 'GAT')
)# 自环
if edge_weight_train is None:edge_index_train, _ = add_self_loops(edge_index_train, num_nodes=X_pca.shape[0])
else:edge_index_train, edge_weight_train = add_self_loops(edge_index_train, edge_attr=edge_weight_train,fill_value=SELF_LOOP_VALUE, num_nodes=X_pca.shape[0])# PyG数据
x_train = torch.tensor(X_pca, dtype=torch.float32)
y_train = torch.tensor(y_labeled, dtype=torch.long)
data_train = Data(x=x_train, edge_index=edge_index_train, y=y_train)
if edge_weight_train is not None:data_train.edge_weight = edge_weight_train# 分层划分
train_idx, test_idx = train_test_split(np.arange(len(y_labeled)),test_size=1 - TRAIN_RATIO,random_state=SEED,stratify=y_labeled
)
train_idx = torch.tensor(train_idx, dtype=torch.long)
test_idx = torch.tensor(test_idx, dtype=torch.long)
print(f"训练样本: {len(train_idx)}，测试样本: {len(test_idx)}")# =========================
# 模型定义
# =========================
class GCNNet(nn.Module):def __init__(self, in_dim, hidden, num_classes, dropout=0.5):super().__init__()self.conv1 = GCNConv(in_dim, hidden, add_self_loops=False, normalize=True)self.conv2 = GCNConv(hidden, num_classes, add_self_loops=False, normalize=True)self.dropout = dropoutdef forward(self, data):x, ei = data.x, data.edge_indexew = getattr(data, 'edge_weight', None)x = self.conv1(x, ei, edge_weight=ew)x = F.relu(x)x = F.dropout(x, p=self.dropout, training=self.training)x = self.conv2(x, ei, edge_weight=ew)return x
#重点学习这个部分
class GATNet(nn.Module):def __init__(self, in_dim, hidden, num_classes, heads=4, concat=True, dropout=0.5):super().__init__()self.gat1 = GATConv(in_dim, hidden, heads=heads, concat=concat,add_self_loops=True, dropout=dropout)out_in = hidden * heads if concat else hiddenself.gat2 = GATConv(out_in, num_classes, heads=1, concat=False,add_self_loops=True, dropout=dropout)self.dropout = dropoutdef forward(self, data):x, ei = data.x, data.edge_indexx = self.gat1(x, ei)x = F.elu(x)x = F.dropout(x, p=self.dropout, training=self.training)x = self.gat2(x, ei)return x# 初始化
if MODEL == 'GCN':model = GCNNet(PCA_DIM, HIDDEN, num_classes=num_classes, dropout=DROPOUT).to(device)
elif MODEL == 'GAT':model = GATNet(PCA_DIM, HIDDEN, num_classes=num_classes,heads=GAT_HEADS, concat=GAT_CONCAT, dropout=DROPOUT).to(device)
else:raise ValueError("MODEL must be 'GCN' or 'GAT'")data_train = data_train.to(device)
train_idx = train_idx.to(device)
test_idx = test_idx.to(device)optimizer = torch.optim.Adam(model.parameters(), lr=LR, weight_decay=WD)
criterion = nn.CrossEntropyLoss()# =========================
# 训练 & 评估 & 早停
# =========================
def train_one_epoch(model, data, train_idx, opt, criterion):model.train()opt.zero_grad()out = model(data)loss = criterion(out[train_idx], data.y[train_idx])loss.backward()opt.step()return loss.item()@torch.no_grad()
def eval_model(model, data, train_idx, test_idx):model.eval()out = model(data)pred = out.argmax(dim=1)tr_acc = accuracy_score(data.y[train_idx].cpu().numpy(), pred[train_idx].cpu().numpy())te_acc = accuracy_score(data.y[test_idx].cpu().numpy(),  pred[test_idx].cpu().numpy())return tr_acc, te_acc, pred.cpu().numpy()best_acc = 0.0
best_state = None
no_improve_epochs = 0for ep in range(1, MAX_EPOCHS + 1):loss = train_one_epoch(model, data_train, train_idx, optimizer, criterion)tr_acc, te_acc, _ = eval_model(model, data_train, train_idx, test_idx)if te_acc > best_acc:best_acc = te_accbest_state = {k: v.detach().cpu().clone() for k, v in model.state_dict().items()}no_improve_epochs = 0else:no_improve_epochs += 1if ep % 10 == 0 or ep == 1:print(f"[{MODEL}|{GRAPH_MODE}] Epoch {ep:03d} | 训练损失: {loss:.4f} | 训练Acc: {tr_acc:.4f} | 测试Acc: {te_acc:.4f}")if no_improve_epochs >= PATIENCE or ep >= MAX_EPOCHS:print(f"停止训练于Epoch {ep:03d}，最佳测试Acc={best_acc:.4f}")breakif best_state is not None:model.load_state_dict(best_state)
model.eval()# =========================
# （1）测试集可视化
# =========================
with torch.no_grad():out_train_graph = model(data_train)pred_train_graph = out_train_graph.argmax(dim=1).cpu().numpy()y_true_test = data_train.y[test_idx].cpu().numpy()
y_pred_test = pred_train_graph[test_idx.cpu().numpy()]cm = confusion_matrix(y_true_test, y_pred_test)
class_names = [f"类{i+1}" for i in range(num_classes)]plt.figure(figsize=(10, 7))
sns.heatmap(cm, annot=True, fmt='d', cmap="Blues",cbar=False, square=True,xticklabels=class_names, yticklabels=class_names)
plt.xlabel("预测标签", fontsize=12)
plt.ylabel("真实标签", fontsize=12)
plt.title(f"{MODEL}（{GRAPH_MODE}）测试集混淆矩阵", fontsize=14, weight='bold')
plt.tight_layout()
plt.show()print("分类报告：")
print(classification_report(y_true_test, y_pred_test, target_names=class_names))class_acc = cm.diagonal() / cm.sum(axis=1).clip(min=1)
plt.figure(figsize=(10, 4.5))
ax = sns.barplot(x=class_names, y=class_acc, edgecolor='black')
plt.ylim(0, 1.0)
for i, v in enumerate(class_acc):ax.text(i, v + 0.02, f"{v:.2f}", ha='center', fontsize=10)
plt.title(f"各类别准确率（测试集）— {MODEL}/{GRAPH_MODE}", fontsize=13, weight='bold')
plt.xlabel("类别"); plt.ylabel("准确率")
sns.despine()
plt.tight_layout()
plt.show()# =========================
# （2）分块全图预测（局部构图 + 内核拼接）
# =========================
def predict_full_image_tiled(model,X_all_pca,h, w,block_size=128,overlap=32,mode='fused',k_spec=6, k_spat=6,alpha=0.7, sig_spec=1.0, sig_spat=1.0,is_gat=False):"""分块进行局部构图并预测，减少内存占用。- 对每个块，使用块内像素的 PCA 特征与坐标构图并预测；- 仅将“内核区域”（去掉 overlap 的中心区域）写回全图，避免边界效应；"""N = h * wpred_all = -1 * np.ones(N, dtype=np.int64)# 全图坐标（一次构造，块内切片）rows_all = np.repeat(np.arange(h), w)cols_all = np.tile(np.arange(w), h)step = max(1, block_size - 2 * overlap)for r0 in range(0, h, step):for c0 in range(0, w, step):r1 = min(r0 + block_size, h)c1 = min(c0 + block_size, w)rr = np.arange(r0, r1)cc = np.arange(c0, c1)RR, CC = np.meshgrid(rr, cc, indexing='ij')idx_block = (RR * w + CC).reshape(-1)X_block = X_all_pca[idx_block]coords_block = np.stack([rows_all[idx_block], cols_all[idx_block]], axis=1).astype(np.float32)# 局部构图edge_index_b, edge_weight_b = build_graph(X_block, coords_block, mode=mode,k_spec=k_spec, k_spat=k_spat, alpha=alpha,sig_spec=sig_spec, sig_spat=sig_spat,to_connectivity_for_gat=is_gat)# 自环if edge_weight_b is None:edge_index_b, _ = add_self_loops(edge_index_b, num_nodes=X_block.shape[0])else:edge_index_b, edge_weight_b = add_self_loops(edge_index_b, edge_attr=edge_weight_b,fill_value=SELF_LOOP_VALUE, num_nodes=X_block.shape[0])data_b = Data(x=torch.tensor(X_block, dtype=torch.float32, device=device),edge_index=edge_index_b.to(device))if edge_weight_b is not None:data_b.edge_weight = edge_weight_b.to(device)with torch.no_grad():out_b = model(data_b)pred_b = out_b.argmax(dim=1).cpu().numpy()# 仅写入内核区域inner_r0 = r0 if r0 == 0 else r0 + overlapinner_c0 = c0 if c0 == 0 else c0 + overlapinner_r1 = r1 if r1 == h else r1 - overlapinner_c1 = c1 if c1 == w else c1 - overlapirr = np.arange(inner_r0, inner_r1)icc = np.arange(inner_c0, inner_c1)if len(irr) == 0 or len(icc) == 0:continueIRR, ICC = np.meshgrid(irr, icc, indexing='ij')inner_idx_global = (IRR * w + ICC).reshape(-1)inner_rr_local = IRR - r0inner_cc_local = ICC - c0inner_idx_local = (inner_rr_local * (c1 - c0) + inner_cc_local).reshape(-1)pred_all[inner_idx_global] = pred_b[inner_idx_local]# 兜底：未覆盖到的像素填0类（可改近邻）if np.any(pred_all < 0):pred_all[pred_all < 0] = 0return pred_allprint("对全图像素构图并预测（分块）... 可能稍耗时")
X_all_scaled = scaler.transform(X_flat)
X_all_pca = pca.transform(X_all_scaled)pred_all = predict_full_image_tiled(model,X_all_pca=X_all_pca,h=h, w=w,block_size=BLOCK_SIZE,overlap=OVERLAP,mode=GRAPH_MODE,k_spec=K_SPEC, k_spat=K_SPAT,alpha=ALPHA, sig_spec=SIGMA_SPEC, sig_spat=SIGMA_SPAT,is_gat=(MODEL == 'GAT')
)def show_pred_map(pred, h, w, title=None):if title is None:title = f"{MODEL} 全图预测（{GRAPH_MODE}，块={BLOCK_SIZE}, 重叠={OVERLAP}）"pred_map = pred.reshape(h, w)n_colors = int(np.max(pred)) + 1 if pred.size else 1cmap = matplotlib.colormaps.get_cmap('tab20').resampled(max(n_colors, 1))plt.figure(figsize=(10, 7.5))im = plt.imshow(pred_map, cmap=cmap, interpolation='nearest')cbar = plt.colorbar(im, ticks=range(n_colors), shrink=0.85)cbar.set_label('类别', rotation=90)plt.title(title, fontsize=14, weight='bold')plt.axis('off')plt.tight_layout()plt.show()show_pred_map(pred_all, h, w)
print("Done.")