ARTICLE DETAIL

资讯详情

深耕郑州网站建设与运营推广的一线实战洞察。

TOA深度学习反演PM2.5:从数据准备到模型训练的完整指南

TOA深度学习反演PM2.5:从数据准备到模型训练的完整指南 简介这是一份基于Python的遥感毕业设计项目聚焦TOA深度学习反演PM2.5面向计算机、人工智能、遥感、环境等专业的在校学生、教师及科研人员也可作为毕业设计、课程设计或项目演示的参考原型。压缩包内共6个文件主要为4个Python脚本、1个Jupyter Notebook和1个TXT说明文档Python脚本涵盖数据准备、特征提取、PM2.5数值提取及深度学习模型训练等环节Notebook则以交互方式呈现TOA反演PM2.5的完整流程TXT文档提供项目说明与环境配置提示。整个资源约12KB代码精简便于快速阅读与二次开发。该项目已有243人学习代码经测试运行成功作者毕业答辩评审平均分达到94.5分质量有保障。读者可以借此掌握遥感参数反演中的数据处理思路、深度学习模型构建方法以及从数据导入到结果输出的全链路实现基础较好的开发者还能在现有代码上修改扩展为其他气溶胶反演或多因子预测任务非常适合作为论文实验或工程落地的起点。1. 毕设选题“TOA深度学习反演PM2.5”到底在做什么第一次看到这个标题的读者多半会卡在“TOA”三个字母上。TOA是Top of Atmosphere大气层顶。传统遥感反演PM2.5的路线是先做大气校正拿到地表反射率再反演气溶胶光学厚度AOD最后用AOD和气象要素建统计模型回归PM2.5。这条路线中间任何一环出错结果都会偏。而“TOA深度学习反演PM2.5”这个方向思路是跳过大isnt校正和AOD反演这两个中间步骤直接把卫星接收到的TOA反射率、观测几何角度、气象辅助数据一起喂给神经网络端到端回归地面站点浓度。这个思路在近几年的论文里被反复验证是可行的因为深度学习模型有能力从TOA信号里隐式学习气溶胶和地表贡献的耦合关系。对毕业设计来说它的价值在于不需要手写大气校正流程不需要依赖MOD04等气溶胶产品数据链条短、可视化效果好、可解释性工作也好做。本文面向正在做遥感方向毕设的学生以及想快速上手的从业者按照“原理—数据—建模—避坑—出图”的顺序把这条路完整走一遍。2. TOA反射率凭什么直接喂给神经网络物理基础和选型理由2.1 TOA反射率不是“没校正的图”它是完整的物理观测很多教材把TOA反射率描述成“需要被校正的原始数据”这个说法容易误导。TOA反射率是卫星传感器在大气层顶接收到的表观反射率它天然包含三部分贡献地表反射、大气分子散射瑞利散射、气溶胶散射与吸收。传统大气校正的目标是把后两者剥离出去但TOA本身是地表与大气的耦合结果信息量比地表反射率更大。在PM2.5反演这个任务里气溶胶信号恰恰是我们想要的。如果先做大气校正等于先把气溶胶信息当成噪声去掉再想办法从AOD产品里找回来——这本身就是一条弯路。深度学习时代的一个常用做法是直接使用TOA数据让模型在训练过程中自己学习“哪些波段组合对气溶胶敏感哪些波段用来估计地表贡献”。从信息论角度看端到端方案少了一个有损压缩步骤。具体到数据形式上常见做法是用Landsat 8/9 OLI、Sentinel-2 MSI或MODIS的TOA产品。Landsat和Sentinel-2的空间分辨率在10~30米适合做城市尺度的精细反演MODIS分辨率250~500米适合做大区域覆盖。毕设场景里我一般建议优先选Landsat因为空间分辨率高、出图好看、和地面站点匹配时定位误差更小。2.2 深度学习如何隐式学习气溶胶—地表解耦无需显式反演AOD的端到端回归传统反演PM2.5的物理链条可以写成TOA → 大气校正 → 地表反射率暗像元 → AOD → 经验模型 → PM2.5。这条链条的每一个环节都有独立的误差来源大气校正假设地表为朗伯体、AOD反演依赖气溶胶模型假设城市型/大陆型/海洋型、经验模型又受气象条件影响。误差逐级累积。深度学习端到端方案的链条是TOA 角度 气象 → 神经网络 → PM2.5。模型内部自行构造特征它可能在第几层等价地实现了一个“暗像元地表反射率估计”在更深的层里实现“AOD代理变量提取”最后回归。这种隐式分解不需要人工指定气溶胶类型训练数据里有多少种大气状况模型就有机会学到多少种“规则”。在实际操作里输入特征通常是这样的组合TOA反射率所有波段如果用的是Landsat 8就是海岸波段、蓝、绿、红、近红外、两个短波红外、卷云波段、太阳天顶角、观测天顶角、相对方位角、以及可选的辅助变量边界层高度、温度、风速。模型输出是站点位置的PM2.5浓度数值。这里要特别提醒不要只给模型TOA反射率而不给角度参数因为TOA反射率随观测角度变化很大没有角度信息模型很难分清“反射率变化是气溶胶变化还是几何变化导致的”。2.3 为什么不用传统的遥感随机森林或PROSAIL反演方案有读者会问既然要端到端用随机森林回归、支持向量机也可以为什么非要深度学习答案是随机森林对特征交互的建模能力有限。PM2.5浓度和TOA反射率的关系是高度非线性的并且依赖空间邻域信息——某像元的TOA值可能没变但它周围像元的纹理变化可能意味着云影或气溶胶团块。卷积神经网络能提取这种空间上下文随机森林只能逐像元处理。同理像PROSAIL反演LAI这类基于物理模型的方法是正向建模的思路参数化地表和大气然后最小化模拟与观测的差异。这种方案对先验知识要求极高需要知道气溶胶类型、地表BRDF参数PM2.5又没有直接的物理光学信号——它通过气溶胶消光间接影响TOA且受湿度影响很大。用PROSAIL反演LAI很成熟但没有任何一款物理模型能直接从TOA反射率正向模拟出PM2.5浓度。深度学习在这里不是“黑匣子替代物理模型”而是物理模型压根覆盖不了的场景下唯一可行的反演范式。3. 数据准备从遥感影像到归一化训练样本3.1 数据源选择Landsat 8/9、Sentinel-2与站点数据的时间窗口匹配策略毕设反演PM2.5最理想的数据源是Landsat系列。Landsat重访周期16天单景覆盖范围185公里空间分辨率30米全色15米对PM2.5这种在几十公里尺度上有空间自相关的大气污染物来说足够用。如果做赛季级别的分析可以叠加Landsat 8和Landsat 9两颗卫星重访周期缩短到8天。Sentinel-2重访周期5天10米分辨率但只有13个波段、部分波段的信噪比在气溶胶敏感区表现不如Landsat稳定。地面PM2.5站点数据国内能拿到的是环境监测总站的逐小时浓度。和卫星匹配时要注意时间窗口卫星过境是当地时间上午10:00左右站点数据取卫星过境时刻前后1小时的平均值。不要直接用整点值因为地面浓度小时变化剧烈10:30过境的卫星和10:00整点浓度可能差异很大。取前后1小时平均能平滑这一误差。如果站点数据是日均值那就只能和卫星做粗略对齐模型的R²会明显下降。还有一个毕设常见的问题是站点数量不够。一个省的国控站点大约100~200个匹配到晴空无云的Landsat影像上一次过境能匹配到的有效样本可能只有几十到一百条。这是正常现象不需要因为样本少就放弃。解决办法在后面第4章会提到——用数据增强和小模型。3.2 用GEE批量导出TOA数据和角度数据避免本地大气校正的坑推荐用Google Earth EngineGEE的Python API来批量获取TOA反射率而不是本地下载Landsat Level-1产品再自行处理。GEE里Landsat 8/9 Collection 2 Level-1数据集已经完成了辐射定标和几何校正直接取反射率波段即可不需要自己写辐射定标公式。更关键的是GEE导出的数据包含SR_BTOA反射率以及太阳/观测角度波段不需要额外下载元数据文件来解析角度。import ee import geemap import pandas as pd ee.Initialize() # 以某站点坐标为中心做一个100m x 100m的缓冲区 lon, lat 116.40, 39.90 # 示例北京某站点 point ee.Geometry.Point([lon, lat]).buffer(100) dataset ee.ImageCollection(LANDSAT/LC08/C02/T1_L2) \ .filterBounds(point) \ .filterDate(2020-01-01, 2021-12-31) \ .filter(ee.Filter.lt(CLOUD_COVER, 20)) def extract_toa(img): # 取TOA反射率波段Collection 2 Level-1中为SR_B开头的波段 # 注意L2产品里是SR_BL1产品里才是SR_B的原始TOA # 这里用L1产品的SR_B并配合对应的辐射定标参数更常规 # 但GEE的L2产品比较直接配合角度波段即可 return img.select( [SR_B1, SR_B2, SR_B3, SR_B4, SR_B5, SR_B6, SR_B7, SR_B9, SR_B10, SR_B11], [coastal, blue, green, red, nir, swir1, swir2, cirrus, therm1, therm2] ).addBands(img.select([SAA, SZA, VAA, VZA])) \ .set(system:time_start, img.get(system:time_start)) collection dataset.map(extract_toa) # 导出每个站点位置的TOA向量 def sample_point(img): return img.sample(point, scale30).first() samples collection.map(sample_point).getInfo()以上代码的逻辑分三步先按站点坐标过滤影像集合然后重命名波段并追加太阳和观测角度波段最后用sample在缓冲区中心提取该位置上的TOA值。需要注意的点是scale参数设为30和Landsat分辨率一致sample返回的是一个Feature不是逐像元数组后续需要转成DataFrame。如果你想要一个空间patch比如以站点为中心画5x5的窗口需要用sampleRectangle而不是sample。还有一个容易出错的地方Landsat Collection 2 Level-1的TOA产品在GEE里的波段名是SR_B开头的但Level-2地表反射率产品波段名也是SR_B开头两者看起来一样、数值含义完全不同。一定要确认自己用的数据集ID是LANDSAT/LC08/C02/T1_L1才是TOA产品T1_L2是地表反射率。上面的示例代码里用的是L2产品如果你严格要TOA需要改成L1并使用multiply缩放系数。# 如果使用Level-1产品做TOA需要手动施加缩放系数 import numpy as np def apply_scale(img): # L1产品的TOA反射率需要乘以0.00001才能得到真实反射率 toa_bands img.select( [SR_B1, SR_B2, SR_B3, SR_B4, SR_B5, SR_B6, SR_B7, SR_B9] ).multiply(0.00001).float() angles img.select([SAA, SZA, VAA, VZA]).float() return toa_bands.addBands(angles).set( system:time_start, img.get(system:time_start) )参数说明0.00001是Landsat Collection 2 Level-1产品的TOA反射率缩放因子。热红外波段SR_B10、SR_B11在这个产品里是绝对辐射亮度单位W/(m²·sr·μm)需要另一个公式转换一般不放进反演PM2.5的模型输入因为热红外波段对气溶胶不敏感反而引入地表温度干扰。cirrus波段卷云波段建议保留它能有效识别薄卷云污染。3.3 站点匹配与数据集划分按时间切分防止时空泄漏拿到每个站点的TOA向量之后要和地面PM2.5站点做匹配然后用一句话概括最重要的数据集划分原则按时间划分。如果随机划分训练集和测试集同一天内各站点的样本有强烈的空间自相关——模型记住站点位置就能“猜对”测试集R²会虚高到0.9以上看起来完美实际完全不可用。import pandas as pd import numpy as np from sklearn.model_selection import GroupShuffleSplit # df_merged: 每行是一个站点某次卫星过境的TOA向量PM2.5浓度 # 关键列: date(日期), station_id(站点ID), pm25(浓度), feature列 # 按日期分组保证同一天的样本全进同一个集合 gss GroupShuffleSplit(n_splits1, test_size0.2, random_state42) train_idx, test_idx next(gss.split(df_merged, groupsdf_merged[date])) df_train df_merged.iloc[train_idx] df_test df_merged.iloc[test_idx] # 验证训练集和测试集的日期不重叠 assert set(df_train[date]).isdisjoint(set(df_test[date])), 日期重叠 # 对特征做标准化注意只用训练集统计量 feature_cols [coastal, blue, green, red, nir, swir1, swir2, cirrus, SZA, VZA, SAA, VAA] mean df_train[feature_cols].mean() std df_train[feature_cols].std() df_train[feature_cols] (df_train[feature_cols] - mean) / std df_test[feature_cols] (df_test[feature_cols] - mean) / std这里有两个容易被忽略的细节。一个是标准化必须只用训练集的均值和标准差不能混入测试集信息否则就是数据泄漏。二是把太阳天顶角SZA、观测天顶角VZA也一起标准化角度变化范围大不标准化会让神经网络训练不稳定。另外PM2.5浓度本身建议做对数变换因为浓度分布是右偏的直接用原始值会导致模型对高浓度样本过拟合而忽略中低浓度样本。做完匹配和划分后常规做法是把特征存成npy格式方便后续模型加载特征形状是(n_samples, n_features)。如果你想做空间patch输入那就在GEE里用sampleRectangle导出一个(n_samples, patch_size, patch_size, n_bands)的四维数组python侧再另外处理。4. 模型搭建与训练CNN回归与关键超参数4.1 选型为什么用一维CNN或浅层二维CNN而不是Transformer和MLP样本量决定了模型规模。一个省三年Landsat数据、过境晴空样本大约2000~3000条Transformer在这种数据量下必过拟合。MLP多层感知机能建模特征交互但无法利用空间邻域信息。一个务实的选型是基础方案用MLP做快速baseline进阶方案使用一维CNN把各波段当成序列建模或者用浅层二维CNN处理站点邻域patch。二维CNN的方式更适合“以站点为中心取patch”的数据组织方式。这样子输入不再是把每个站点压成一个向量而是保留空间结构——模型能看到站点周边大气状况的空间分布。常见做法是取32x32或16x16的patch输入通道是10个波段输出是PM2.5标量。为了控制参数规模只用4~6个卷积层每层卷积核数量从16开始逐层翻倍通道数不超过128。为什么不能直接抄图像分类的ResNet结构因为图像分类有ImageNet千万级预训练数据而你的遥感样本只有几千条模型一旦超过5百万参数就必然过拟合。同时PM2.5浓度与像素值的关系是平滑的回归关系不需要特别深的网络来提取高度抽象特征。一个反直觉的经验是在这个任务上3层卷积的模型比18层ResNet效果更好因为深层网络把小样本里的噪声也当成了特征。import torch import torch.nn as nn class PM25CNN(nn.Module): def __init__(self, n_bands10, patch_size32): super().__init__() # 输入: (batch, n_bands, patch_size, patch_size) self.features nn.Sequential( nn.Conv2d(n_bands, 16, kernel_size3, padding1), nn.ReLU(), nn.BatchNorm2d(16), nn.MaxPool2d(2), # 16x16 nn.Conv2d(16, 32, kernel_size3, padding1), nn.ReLU(), nn.BatchNorm2d(32), nn.MaxPool2d(2), # 8x8 nn.Conv2d(32, 64, kernel_size3, padding1), nn.ReLU(), nn.BatchNorm2d(64), nn.AdaptiveAvgPool2d(1) # 全局平均池化 - (batch, 64, 1, 1) ) self.head nn.Sequential( nn.Flatten(), nn.Linear(64, 32), nn.ReLU(), nn.Dropout(0.3), nn.Linear(32, 1) ) def forward(self, x): x self.features(x) return self.head(x).squeeze(-1) model PM25CNN(n_bands10, patch_size32)代码里的两个关键设计值得展开说。BatchNorm2d放在卷积和ReLU之后作用是缓解内部协变量偏移让小样本训练更稳定没有它模型在前500轮训练里loss会剧烈震荡。AdaptiveAvgPool2d把任意尺寸的特征图压成1x1这样模型对patch大小不敏感如果你想在测试时改用更大的patch来预测比如64x64网络结构不需要改只改输入尺寸就行。Dropout(0.3)是专门用来对抗过拟合的。在小样本回归任务里Dropout比L2正则化更有效因为它同时做了模型集成。不过注意不要加太高0.5以上的Dropout在样本量很小时反而会让训练不收敛。4.2 训练循环与评估指标用RMSE而不是只看R²训练过程的组织比模型结构更容易被忽视。毕设里最常犯的错误是只用一个固定学习率跑几百轮然后报告一个好看的训练集R²。下面这段代码是一个带学习率衰减和早停的完整训练循环通用性很强。import torch.optim as optim from torch.utils.data import DataLoader, TensorDataset def train_model(model, X_train, y_train, X_val, y_val, epochs200): # X: (n, n_bands, patch, patch), y: (n,) train_ds TensorDataset(torch.FloatTensor(X_train), torch.FloatTensor(y_train)) val_ds TensorDataset(torch.FloatTensor(X_val), torch.FloatTensor(y_val)) train_loader DataLoader(train_ds, batch_size32, shuffleTrue) val_loader DataLoader(val_ds, batch_size32, shuffleFalse) optimizer optim.AdamW(model.parameters(), lr1e-3, weight_decay1e-4) scheduler optim.lr_scheduler.ReduceLROnPlateau( optimizer, modemin, factor0.5, patience10 ) criterion nn.MSELoss() best_val_loss float(inf) best_state None patience_counter 0 for epoch in range(epochs): model.train() train_loss 0.0 for xb, yb in train_loader: optimizer.zero_grad() pred model(xb).squeeze(-1) loss criterion(pred, yb) loss.backward() optimizer.step() train_loss loss.item() * xb.size(0) model.eval() val_loss 0.0 with torch.no_grad(): for xb, yb in val_loader: pred model(xb).squeeze(-1) loss criterion(pred, yb) val_loss loss.item() * xb.size(0) train_loss / len(train_ds) val_loss / len(val_ds) scheduler.step(val_loss) # 早停连续20轮验证loss不降就停止 if val_loss best_val_loss: best_val_loss val_loss best_state {k: v.clone() for k, v in model.state_dict().items()} patience_counter 0 else: patience_counter 1 if patience_counter 20: print(fEarly stop at epoch {epoch}) break model.load_state_dict(best_state) return model这里有两个参数决定成败。第一个是AdamW里的weight_decay设1e-4它是权重衰减和L2正则化本质一样防止权重变大导致过拟合。第二个是ReduceLROnPlateau的patience参数设为10的意思是在验证loss连续10轮不降后把学习率减半比固定学习率衰减策略灵活得多——它只在真正“卡住”时才调整。评估指标上除了R²一定要报告RMSE均方根误差。R²高不代表模型好如果测试集浓度范围小R²天然就低反过来样本量小、范围跨度大R²可能虚高。RMSE单位是μg/m³能直观反映预测误差。另外建议把预测值和真实值的散点图加上1:1线比任何数值指标都有说服力。4.3 数据增强小样本遥感反演的后悔药当有效样本只有几百条时数据增强不是可选项是必选项。遥感TOA数据的增强方式和图像分类不同不能用随机裁剪、旋转90度这些会破坏物理意义的操作——旋转后太阳方位角信息就错了。适合的做法是添加噪声、通道扰动和微小几何偏移。def augment_toa(batch_x, batch_y, noise_std0.01): 对TOA反射率patch添加高斯噪声和通道扰动 x batch_x.clone() # 1. 高斯噪声模拟传感器噪声 noise torch.randn_like(x) * noise_std x x noise # 2. 通道缩放扰动模拟不同大气条件下波段间相对强度变化 scale torch.empty(x.size(0), x.size(1), 1, 1).uniform_(0.95, 1.05) x x * scale # 3. 随机翻转只在空间维度翻转不改变物理含义 if torch.rand(1) 0.5: x torch.flip(x, dims[2]) if torch.rand(1) 0.5: x torch.flip(x, dims[3]) return x, batch_y这里的核心思想是加噪声而不是改变几何结构。高斯噪声的std设为0.01是参考Landsat TOA反射率的典型传感器噪声水平。通道缩放的幅度在±5%以内超过这个范围会破坏波段间比例关系模型学到的是被污染的分布。随机翻转对部分方向性特征有影响但因为patch很小32x32空间纹理信息权重不高实际影响有限。在训练时调用增强要注意频率不是每个epoch都对全部数据做增强那样模型可能见过太多变形版本反而学不到本征特征。常见做法是每轮以50%概率对每个batch做一次增强其余时间用原始数据。5. 常见问题排查5个让反演结果翻车的细节5.1 测试集R²高达0.9但站点对比曲线完全对不上现象模型指标很漂亮训练集R²0.95、测试集R²0.88但画出测试集每天的平均预测浓度曲线和真实浓度曲线明显偏离尤其在重污染时段。原因这是典型的时空泄漏。随机划分数据集时同一天的站点数据部分进了训练集、部分进了测试集。因为同一天内PM2.5浓度在空间上有强相关性模型实际上在“记忆”当天的浓度水平而不是学习TOA和PM2.5的关系。跨日期的泛化能力一检验就露馅。解决严格按日期分组划分数据已经在前文代码中展示。更严谨的做法是拿一整年的数据做测试——训练集用前两年测试集用第三年全年这样检验的是模型对新日期、新气象条件的泛化能力论文里审稿人也认。5.2 云和冰雪的TOA特征被模型误判为高PM2.5现象预测结果里云覆盖区域的PM2.5浓度明显偏高冬天有积雪时整片区域都显示出虚假的污染高值。原因云和雪的TOA反射率在可见光波段远高于气溶胶的贡献模型学到了“反射率高→浓度高”的粗糙规则。如果GEE导出数据时只用CLOUD_COVER小于20%这个整体阈值筛选影像但站点patch内部仍有零星云块没被识别就会混入噪声样本。解决在GEE端用QA_PIXEL波段做逐像元云掩膜把掩膜后的像素置为NaN后续建模时剔除。同时在特征里把质量波段也保留下来。下面是GEE端云掩膜的完整做法。def mask_clouds(img): qa img.select(QA_PIXEL) # Bit 3: cloud, Bit 4: cloud shadow cloud_mask qa.bitwiseAnd(1 3).eq(0).And( qa.bitwiseAnd(1 4).eq(0) ) return img.updateMask(cloud_mask) # 使用方式 collection dataset.map(lambda img: mask_clouds(img).select( [SR_B2, SR_B3, SR_B4, SR_B5, SR_B6, SR_B7, SR_B9, SAA, SZA, VAA, VZA] ))要注意的是云掩膜会让patch内部产生空洞有效像元数减少。所以在导出patch时统计有效像元比例低于60%的样本直接丢弃否则模型输入里大量填充值训练会被干扰。5.3 模型在白天晴空场景一切正常有薄霾时预测崩盘现象能见度好的日子里预测值和实测值吻合一旦出现中度霾AQI150模型预测值严重偏低甚至出现负值。原因薄霾场景下TOA反射率变化有双重效应——气溶胶增加会增大可见光波段的TOA反射率但也增加了大气吸收可能降低特定波段的信号。训练数据里中度和重度污染样本占比很小通常不到10%模型根本没机会学到这个区间。解决一是对PM2.5浓度做对数变换后再回归缓解极端值权重问题。二是使用加权损失函数给高浓度样本更大的权重迫使模型关注这些少数样本。代码里把MSELoss换成加权版本权重设为浓度值的平方根浓度越高的样本权重越大。def weighted_mse_loss(pred, target): weight torch.sqrt(target 1) # 浓度高权重大 loss weight * (pred - target) ** 2 return loss.mean()这样做能显著改善重污染时段的拟合但不能根治——根治需要补充更多重污染样本。毕设阶段如果实在凑不齐可以在论文讨论里诚实地说明该限制比模型假装能预测所有情况要更好。5.4 Landsat和Sentinel-2数据混用时模型训练loss反复震荡现象把Landsat 8和Sentinel-2的数据混合在一起训练模型loss在训练过程中忽高忽低验证集指标比只用单源数据还差。原因两代传感器的波段设置不一致。Landsat 8有海岸波段433-453nm和卷云波段1360-1390nmSentinel-2没有海岸波段但有B1443nm和B101375nm且两者的波段半高宽不同。模型输入通道语义不一样网络表层会试图同时适应两种分布在特征空间里造成冲突。解决只保留两者共有的波段蓝、绿、红、近红外、两个短波红外统一用Landsat的分辨率重采样并用公式校正波段差异。更省事的办法是不混用数据源直接用单源数据做完整实验把混合数据作为对比实验来说明模型的泛化能力。5.5 初始patch大小选错导致模型对站点位置偏移极度敏感现象模型训练好了但把预测patch从站点中心偏移5个像素预测值就波动30%以上周边有河流或山体时更明显。原因patch太小比如8x8模型依赖的中心像元权重过大邻域信息不足patch太大比如128x128则引入了过多远处的地表信息而这些地表信息与站点浓度无关。解决32x32是经过实验验证的折中选择。如果站点周边地表类型单一比如农田16x16也能用如果站点在城市中周边几公里内地表类型变化大就需要64x64。判断patch是否合适的直观标准把模型预测结果做成浓度空间分布图看站点周围是否出现以站点为中心的不自然“斑点”——如果有说明patch大小不匹配空间相关尺度。6. 落地验证与出图把反演结果做成年份序列和空间分布图模型训练完成只是第一步毕业设计答辩时真正拉开差距的是验证和出图。建议至少做三张图预测值vs实测值散点图带1:1参考线、按月份聚合的时间序列对比图、以及一张城市尺度的PM2.5空间分布反演图。散点图用matplotlib绘制1:1线加上去可以直观看到模型是否系统性偏低或偏高。时间序列对比图更有说服力把2020年12个月的站点平均实测浓度和模型预测浓度画在同一个坐标系里用两条曲线对比。如果两条曲线趋势一致但峰值有滞后或超前说明模型对气象条件变化敏感度有问题。空间分布图是最费功夫也最出效果的一步。使用GEE导出的采样点TOA数据作为输入把研究区范围划分成32x32的网格逐网格预测PM2.5浓度再插值成连续分布图。这里有个经验不要对每个网格独立预测后直接拼接因为相邻网格的预测值在边界处会突变——先对研究区影像做重采样和切块预测完所有块后用高斯滤波做一次轻微平滑可以降低这种“拼图效应”。import numpy as np import matplotlib.pyplot as plt from scipy.ndimage import gaussian_filter # pred_grid: 预测出的二维浓度场假设是 (rows, cols) pred_grid_smoothed gaussian_filter(pred_grid, sigma1.0) fig, axes plt.subplots(1, 2, figsize(12, 5)) im0 axes[0].imshow(pred_grid, cmapRdYlBu_r, vmin10, vmax120) axes[0].set_title(Raw model prediction) im1 axes[1].imshow(pred_grid_smoothed, cmapRdYlBu_r, vmin10, vmax120) axes[1].set_title(After Gaussian smoothing) plt.colorbar(im1, axaxes[1], labelPM2.5 (μg/m³)) plt.tight_layout() plt.show()sigma1.0是常用经验值具体大小要看像元对应的地面距离。每个像元对应30米时sigma1.0的平滑相当于做了约90米范围的近邻平均如果像元对应500米MODIS反演sigma取1.5才合适。最后一章想分享一个我自己的习惯模型训练结束后不急着跑全流程出图而是先随机挑5个站点单独画它们各自的实测vs预测时间序列。这个动作几乎每次都能暴露问题——整个测试集指标漂亮不代表单站点时间尺度上可靠。如果你按这套流程做到这一步把站点、数据和代码整理好论文的“讨论与展望”部分就专门写模型对极端污染事件捕捉能力的局限导师认可度会明显高于只报R²的应付性结论。希望这个方向的经验能帮到你的选题和实现环节。本文还有配套的精品资源点击获取
返回列表