天气降尺度一直是气象和机器学习交叉领域的热门话题,但大多数方法都停留在“用粗网格大尺度预报去插值出细网格局地天气”这个框框里。遇到山地、城市、海岸线这种地形复杂区域,插值出的结果往往平滑得失真,原因很简单:粗网格只能记住区域的均值,却丢掉了一个网格内部更小尺度的空间结构。最近有一类研究关注“亚网格描述符”,也就是想办法把网格内部的高分辨率信息提炼出来,辅助降尺度。而“Earth observation embeddings are effective sub-grid descriptors for probabilistic weather downscaling”这篇论文标题指向的,正是用地球观测嵌入来承担这个角色。
这篇文章想把这套思路讲透:什么是地球观测嵌入,为什么它能当亚网格描述符,概率性降尺度怎么和它结合,以及如果你想在自己的气象或遥感项目里复现类似框架,应该从哪些细节入手。
文章会先解释传统降尺度的困境,再拆解嵌入、描述符、概率预测三个概念,然后给出一套可落地的 PyTorch 示例流程,最后补充常见坑和工程建议。即使你没有气候科学背景,只要能写 Python、用过深度学习框架,也能理解这条技术路线的价值。
1. 降尺度任务里的核心矛盾:网格分辨率和局地信息
1.1 什么是降尺度
降尺度(downscaling)在气象学里做的事情很简单:全球气候模式或数值天气预报通常在一个很粗的网格上模拟大气,比如一个格子边长几十甚至上百公里。可你真正关心的是一个风电场、一个农场、一座城市的精细天气,这就要把粗网格信息“细化”到几公里甚至几百米的尺度。
但问题不在于“插值”。简单的双线性插值只是把大尺度的平滑场变成更密的平滑场,并不会产生真实的局地细节。例如,一个 30km 粗网格内可能包含山地、湖泊、城市,插值后这些地形影响完全体现不出来。要做到有效的降尺度,就必须引入额外信息,这个信息要能描述粗网格内部的小尺度状态,也就是“亚网格信息”。
1.2 网格均值掩盖了什么
任何数值模式,本质上是对真实大气的有限体积近似。一个网格内的物理量,比如温度、湿度、风速,模式只存一个平均值。但这个平均值不代表真实空间分布:山区迎风坡和背风坡差异很大,城市热岛和郊区差异很大,农田和森林的蒸散发差异也很大。
真实环境中,这些差异恰恰是决策者最关心的。如果只给一个网格平均的风速,风电场无法估算单台风机的出力;如果只给一个网格平均的降水,城市内涝预报无法知道哪一个街区会出现短时强降水。所以降尺度模型的基本需求是:在粗网格输入之外,找到一些能反映亚网格非均一性的特征。
1.3 传统替代方案的局限
传统统计降尺度会用经纬度、海拔、到海岸线距离等静态变量作为附加信息。这些变量确实捕捉了部分地形效应,但静态变量不包含时间变化。云的分布、土壤湿度、雪盖、植被状态这些随时间变化的亚网格因子,无法从静态变量中获得。也有工作尝试用高分辨率数值模拟来提供亚网格信息,但计算成本极高,很难在业务中推广。
因此,一个新的思路自然出现:从卫星等地球观测数据中学习一个紧凑的嵌入向量,用它作为动态的亚网格描述符。卫星影像每过几天甚至每小时就能覆盖全球,分辨率可以达到米级到公里级,既能反映静态地形也能反映动态地表状态。问题是,原始影像像素太多,不能直接塞进降尺度模型,必须压缩成一个低维、信息丰富的表示——这就是嵌入。
2. 地球观测嵌入、亚网格描述符、概率预测三者的关系
2.1 什么是嵌入
嵌入(embedding)这一概念在自然语言处理里用得最多,比如把“国王”和“男人”的关系映射到向量空间。在地球观测领域,嵌入通常指用一个神经网络把高维遥感数据压缩成低维稠密向量。例如,一张 256×256 的多光谱影像,可以通过编码器得到 256 维向量。这个向量要尽量保留对目标任务有用的信息,丢掉与任务无关的噪声。
从数学上看,嵌入就是学习一个映射:
[ z = f_\theta(x) ]
其中,(x) 是地球观测数据,(z) 是嵌入向量,(f_\theta) 是神经网络编码器。如果训练方式得当,(z) 中会包含云系统形态、地表纹理、空间异质性等对天气条件敏感的因子。这正是亚网格描述符需要的属性。
2.2 亚网格描述符需要什么性质
一个合格的亚网格描述符,在功能上应该满足三个条件:
- 可分辨亚网格结构:它要能区分两个大网格均值相同但内部空间结构完全不同的区域。
- 与目标气象要素相关:描述符里包含的信息必须对降水、温度、风等预测目标有解释力。
- 方便与粗网格数据融合:它不能只是高分辨率图像,而应该是一种可以和粗网格气象变量一起输入模型的向量或特征图。
地球观测嵌入天然满足前两点:遥感数据记录了地表和大气下垫面的真实状态,和天气过程高度相关;嵌入又是低维向量,方便和气象变量拼接。这里的关键是,嵌入不是人工设计的特征,而是从数据中自动学习得到的,因此它可能捕捉到那些很难用海拔或经纬度显式表达的模式。
2.3 概率性降尺度的价值
传统降尺度模型通常输出一个确定性的数值,比如“明天下午 3 点气温 26.5 摄氏度”。但天气预报本质上是概率问题,因为初始场误差、模式误差和混沌效应都导致结果有不确定性。对于高风险决策,用户更需要知道 26.5 度的不确定性有多大,比如 95% 置信区间是 24 到 29 度,还是 25 到 28 度。
概率性降尺度不再预测单一点位数值,而是预测一个概率分布。这个分布可以是参数化的,例如假设预测目标服从正态分布,模型输出均值和方差;也可以是非参数化的,比如使用分位数回归输出多个分位数,或者使用条件生成模型采样。
将地球观测嵌入引入概率性降尺度,不仅有助于提高确定性预测的精度,还能改善不确定性估计。因为嵌入补充了亚网格信息,模型可以区分出“这个区域的小尺度过程不确定性很大”和“这个区域过程相对稳定”,从而更合理地预测方差和分布形态。
3. 论文思路拆解:算法与架构层面发生了什么
从论文标题看,研究提出的方法不是简单的“把卫星图像像素跟气象数据拼接”,而是用了一个两阶段框架:
3.1 第一阶段:从地球观测数据学习嵌入
这一阶段的目标是获得一个高复用性的编码器。常见做法有两种:
- 自监督训练:使用大量卫星影像做对比学习,比如 SimCLR、MoCo 或 MAE,让模型学习影像中具备尺度不变性和旋转不变性的语义特征。这样学到的嵌入可以被下游任务直接用。
- 端到端训练:在降尺度任务中同步训练嵌入编码器和降尺度模型,嵌入因任务而定制。
论文标题里强调“嵌入有效”,说明重点可能放在验证预训练嵌入的迁移能力上。如果使用预训练编码器,下游训练成本和数据需求都会降低,对于样本稀缺的气象研究非常友好。
3.2 第二阶段:嵌入与粗网格气象数据融合
粗网格气象数据(比如 ERA5 再分析数据)以网格场形式输入降尺度模型,同时地球观测嵌入以特征向量或特征图的形式参与融合。具体架构可以有多种选择:
- 拼接:嵌入向量与气象网格特征在通道维度拼接,再送入解码器。
- 注意力融合:气象场作为 query,嵌入作为 key/value,让模型自适应地决定哪里需要侧重亚网格信息。
- 分层注入:在解码器的多个尺度注入嵌入,保证从全局到局部的信息都有嵌入的参与。
由于嵌入显式携带了亚网格空间结构信息,理论上模型不再需要完全依靠粗网格去“猜测”局地细节,因此降尺度结果的空间结构会更真实。
3.3 概率输出的建模
概率性预测在架构上通常体现为一个额外的输出头。最简单的做法是让神经网络输出两个通道:一个表示预测均值 (\mu),另一个表示预测方差 (\sigma^2),并假设观测服从高斯分布。训练时使用负对数似然损失:
[ \mathcal{L} = \frac{1}{N}\sum_i \left[ \frac{1}{2}\log(\sigma_i^2) + \frac{(y_i - \mu_i)^2}{2\sigma_i^2} \right] ]
如果想表达更复杂的分布,也可以输出分位数。论文标题中的“probabilistic”并不限定是哪种概率形式,但无论用哪种,评估时都需要使用概率性指标,如连续排序概率分数(CRPS)、区间覆盖率等。
3.4 为什么嵌入能让概率预测更好
很多人会忽视嵌入对不确定性校准的影响。概率预测最难的不是输出一个分布,而是让这个分布与真实观测的频率一致。嵌入提供的高分辨率地表信息,可以告诉模型哪些区域更“难以预测”。比如,有积云活动的区域,降水变化剧烈,不确定性应更大;晴空均一的区域,不确定性相对小。没有嵌入时,模型只能泛泛地预测一个大平均方差;有嵌入后,模型可以根据局部纹理和云特征动态预测空间变化的方差,这样得到的概率预测才更有校准价值。
4. 如果你想自己复现:数据和技术栈准备
4.1 需要准备哪些数据
要复现类似研究,通常需要三类数据:
- 粗网格气象场:例如 ERA5 再分析数据,分辨率约 25-30km。这是降尺度模型的输入,包含温度、湿度、风、气压等多个变量。
- 高分辨率观测或再分析目标:例如地面气象站观测、雷达降水反演或高分辨率区域模式输出。这些作为训练标签,分辨率一般在几公里或更细。
- 地球观测影像:例如 Sentinel-2(10-20m 分辨率)、MODIS(250-1000m)、Landsat(30m),需要覆盖研究区域和时间段。考虑到云遮挡和重访周期,也可以使用低轨气象卫星的高频数据。
数据预处理的核心挑战是空间对齐和时间匹配。粗网格的每个格点需要对应一个“邻域”的地球观测窗口,而不是只能取单个像素。这个邻域窗口就是我们想提炼成嵌入的区域。
4.2 技术栈与软件环境
推荐的 Python 技术栈如下:
- xarray:处理网格气象数据;rioxarray:读取遥感栅格。
- PyTorch / TensorFlow:实现嵌入编码器和降尺度模型。
- einops:方便张量维度变换。
- pytorch-lightning 或 Trainer:简化训练循环。
- omegaconf 或 yaml:管理实验配置。
在开始之前先建立一个最小环境:
conda create -n downscale python=3.10 conda activate downscale pip install torch xarray rioxarray netCDF4 rasterio einops pytorch-lightning这一套环境足够处理 ERA5 下载的 NetCDF 文件和卫星 GeoTIFF 文件。
4.3 需要注意的数据版权和授权
ERA5 和大部分卫星数据是开放使用的,但需要遵守相应许可。特别是如果以后要发表论文或商用,要确认数据的许可版本。例如,Sentinel 数据使用 Copernicus 计划开放协议。MODIS 数据来自 NASA,公开可用。无论如何,都要在文章中引用数据源。
5. 简化示例代码:从嵌入到概率性降尺度
下面给出一个简化但完整的 PyTorch 示例,演示如何把地球观测嵌入和粗网格气象输入结合,输出概率预测。代码不追求与某篇论文官方实现一致,只用于展示整体思路。
5.1 定义地球观测编码器
假设地球观测数据是小的图像块,尺寸为 C×H×W,我们用一个小型 CNN 编码器提取嵌入向量:
# 文件路径:eo_encoder.py import torch import torch.nn as nn class EOEncoder(nn.Module): """ 输入: 地球观测图像块 (B, C_EO, H, W) 输出: 嵌入向量 (B, D) """ def __init__(self, in_channels=4, embedding_dim=64): super().__init__() self.features = nn.Sequential( nn.Conv2d(in_channels, 32, kernel_size=3, stride=2, padding=1), nn.ReLU(), nn.Conv2d(32, 64, kernel_size=3, stride=2, padding=1), nn.ReLU(), nn.Conv2d(64, 128, kernel_size=3, stride=2, padding=1), nn.ReLU(), nn.AdaptiveAvgPool2d(1), ) self.fc = nn.Linear(128, embedding_dim) def forward(self, x): x = self.features(x) x = x.view(x.size(0), -1) return self.fc(x)这个编码器会输出一个 64 维向量。如果只输入单张影像,它会丢失位置信息;在实际系统中,你可能会把中心点的嵌入和相邻网格的嵌入序列一起输入,帮助模型感知空间上下文。
5.2 降尺度模型:融合气象网格与嵌入
粗网格气象场是一个 (B, C_met, H_grid, W_grid) 的张量,嵌入是 (B, D) 的向量。最简单的融合方式是把嵌入复制到每个网格位置,然后沿着通道维度拼接。
# 文件路径:downscale_model.py import torch import torch.nn as nn from eo_encoder import EOEncoder class ProbabilisticDownscaleModel(nn.Module): """ 输入: met_field: (B, C_met, H, W) 粗网格气象场 eo_obs: (B, C_EO, h, w) 地球观测图像块 输出: mu: (B, 1, H, W) 预测均值 log_var: (B, 1, H, W) 预测对数方差 """ def __init__(self, met_channels=8, eo_in_channels=4, embedding_dim=64): super().__init__() self.eo_encoder = EOEncoder(eo_in_channels, embedding_dim) # 粗气象场经过一个小型卷积编码 self.met_encoder = nn.Sequential( nn.Conv2d(met_channels, 64, kernel_size=3, padding=1), nn.ReLU(), nn.Conv2d(64, 128, kernel_size=3, padding=1), nn.ReLU(), ) # 融合嵌入后的解码器 self.decoder = nn.Sequential( nn.Conv2d(128 + embedding_dim, 128, kernel_size=3, padding=1), nn.ReLU(), nn.Conv2d(128, 64, kernel_size=3, padding=1), nn.ReLU(), ) self.mean_head = nn.Conv2d(64, 1, kernel_size=3, padding=1) self.logvar_head = nn.Conv2d(64, 1, kernel_size=3, padding=1) def forward(self, met_field, eo_obs): B, _, H, W = met_field.shape emb = self.eo_encoder(eo_obs) # (B, D) # 把嵌入特征铺到每个网格位置 emb_expanded = emb.view(B, -1, 1, 1).expand(B, -1, H, W) # 气象特征 met_feat = self.met_encoder(met_field) # 拼接 fused = torch.cat([met_feat, emb_expanded], dim=1) decoded = self.decoder(fused) mu = self.mean_head(decoded) log_var = self.logvar_head(decoded) return mu, log_var这个模型同时输出了均值和对数方差,方差通过 log_var 参数化,可以保证方差的非负性。真实项目中,你可能想把嵌入注入到多个尺度,但这里展示的是最小可行版本。
5.3 损失函数与训练流程
概率性降尺度最常用的损失是高斯负对数似然损失。为了数值稳定性,通常让模型输出对数方差。
# 文件路径:loss.py import torch import torch.nn.functional as F def gaussian_nll_loss(mu, log_var, target): """ mu, log_var: (B, 1, H, W) target: (B, 1, H, W) """ # 约束方差下限 log_var = torch.clamp(log_var, min=-10, max=10) var = torch.exp(log_var) # 负对数似然 loss = 0.5 * (log_var + ((target - mu) ** 2) / var) return loss.mean()训练过程与通常的回归任务几乎一样,只是标签不再是一个确定数值,而是观测值本身。我们可以把负对数似然作为损失函数,也可以用分位数损失。如果数据噪声较大,也可以先训练一个确定性模型,再单独拟合方差,但端到端联合训练通常效果更好。
5.4 完整训练循环示例
为了方便复现,使用 PyTorch Lightning 组织训练,避免写太多样板代码:
# 文件路径:train.py import pytorch_lightning as pl import torch import torch.nn as nn from torch.utils.data import DataLoader, Dataset from downscale_model import ProbabilisticDownscaleModel from loss import gaussian_nll_loss class DownscaleLightningModule(pl.LightningModule): def __init__(self, learning_rate=1e-3): super().__init__() self.model = ProbabilisticDownscaleModel() self.lr = learning_rate def training_step(self, batch, batch_idx): met_field, eo_obs, target = batch mu, log_var = self.model(met_field, eo_obs) loss = gaussian_nll_loss(mu, log_var, target) self.log('train_loss', loss) return loss def validation_step(self, batch, batch_idx): met_field, eo_obs, target = batch mu, log_var = self.model(met_field, eo_obs) loss = gaussian_nll_loss(mu, log_var, target) self.log('val_loss', loss, prog_bar=True) return loss def configure_optimizers(self): return torch.optim.Adam(self.parameters(), lr=self.lr)上面的代码省略了自定义 Dataset 的实现,在实际项目中你需要读取 NetCDF 和 GeoTIFF,并按经纬度裁剪对齐。流程并不复杂,但“对齐”这一步花费的调试时间往往比训练还长。
6. 如何验证概率性降尺度效果
很多初学者训练完概率模型后,还是习惯性地只算 RMSE。这在概率预测场景不够,因为 RMSE 只衡量均值预测的好坏,没有衡量不确定性是否可靠。
6.1 指标一:CRPS(连续排序概率分数)
CRPS 是最常用的概率预测评估指标。它度量预测分布的累计分布函数 (F) 与观测值 (y) 之间的差距:
[ CRPS = \int_{-\infty}^{\infty} \Big[ F(x) - \mathbb{1}(x \ge y) \Big]^2 dx ]
CRPS 值越低越好。它的优势是同时考虑了预测分布的校准性和锐度。对于高斯分布,CRPS 有解析表达式,计算很方便:
# 文件路径:evaluate.py import torch import torch.distributions as dist def gaussian_crps(mu, sigma, target): """ 计算高斯分布的 CRPS。 mu: 预测均值 sigma: 预测标准差 target: 观测值 """ dist_norm = dist.Normal(torch.zeros_like(mu), torch.ones_like(mu)) normalized = (target - mu) / sigma phi = torch.exp(dist_norm.log_prob(normalized)) Phi = dist_norm.cdf(normalized) crps = sigma * (normalized * (2 * Phi - 1) + 2 * phi - 1 / (torch.pi ** 0.5)) return crps.mean()这个函数返回的 CRPS 是数值越小越好。
6.2 指标二:区间覆盖率
你可以取预测分布的 90% 置信区间,统计真实观测落在区间内的比例。理想情况下应该接近 90%。如果覆盖率远低于 90%,说明模型过度自信;远高于 90%,说明预测分布过宽,信息量低。
6.3 指标三:确定性对比
仍然可以计算 RMSE 和 MAE,用于和确定性降尺度模型对比。一个好的概率模型应该在 RMSE 上不牺牲太多,同时得到更好的 CRPS。一般来说,使用嵌入的模型应该比不使用嵌入的模型 RMSE 更低、CRPS 更低、空间细节更好。
7. 这类方法适用什么场景
7.1 山地与复杂地形区域
对风速和降水的降尺度来说,地形是最重要的调制因子。地球观测影像中包含了山体阴影、植被分布、云阴影等信息,这些内容与近地面风场和降水关系密切。嵌入模型可以从这些影像中自动提取地形相关的上下文特征,从而让风电场微观选址的局地风速估计更加准确。
7.2 城市与土地利用边界
城市地表与周边农村在热力性质上差异巨大。卫星影像中的不透水面比例、冠层高度、夜间灯光等信息,都能够帮助模型预测城市热岛效应下的温度异常。降尺度模型如果只输入经纬度和海拔,很难区分不同城市形态的差异;嵌入向量可以解决这个问题。
7.3 无站点区域的降水估计
在很多发展中国家或山区,地面气象站稀疏,雷达覆盖不足。卫星降水估计产品往往分辨率较低,且依赖复杂算法。利用高分辨率光学或微波影像的嵌入特征,结合粗网格预报,有可能得到更可靠的局地降水概率分布,对洪水预警和农业生产有直接价值。
7.4 不适合什么场景
如果目标区域非常均一,比如大片开阔洋面或平坦沙漠,亚网格空间异质性很弱,地球观测嵌入的增益可能有限。此时增加嵌入只会增加计算开销,对预测提升不大。同样,如果地球观测数据长时间被云覆盖,缺失严重,嵌入提取会变得不可靠,效果可能不比使用地形变量好。
8. 常见问题与排查思路
下面整理在实现地球观测嵌入 + 概率性降尺度过程中可能遇到的高频问题。
| 问题现象 | 可能原因 | 排查方式 | 解决方案 |
|---|---|---|---|
| 训练 loss 不下降 | 气象与EO数据明显错位 | 可视化同一个样本的气象场和EO影像,检查对齐坐标 | 重新进行地理配准,统一重采样网格 |
| 嵌入向量在训练中崩塌(都一样) | 编码器或解码器用了过强的 dropout / 随机初始化不当 | 打印不同样本的嵌入向量范数和方差 | 改用残差连接,增大编码器容量,或先对编码器做预训练 |
| 预测方差几乎为常数 | 模型对不确定性建模不敏感 | 检查 log_var head 输出是否没参与有效梯度 | 用负对数似然损失替代 MSE,适当增加log_var的梯度 |
| 概率预测过宽或过窄 | 损失函数或评估指标使用不当 | 计算区间覆盖率 | 调整log_var下限;使用CRPS作为早停指标 |
| 训练和验证分布不对齐 | 数据划分存在空间泄漏 | 检查训练集与验证集是否同一片区域 | 按空间块划分数据(leave-location-out) |
| EO影像缺失过多 | 云遮挡或卫星重访周期长 | 统计有效像素比例 | 使用合成影像或多时相影像聚合,或采用云掩码 |
一个经常被忽略的点是:预训练嵌入和下游任务之间的分布偏移。如果你从 ImageNet 预训练模型或者遥感领域专用预训练模型得到嵌入,它是在某种传感器和光照条件下学习的。如果你的 EO 影像来自不同传感器或不同季节,嵌入分布可能漂移。解决办法是在下游数据上微调编码器,或者做通道统计归一化。
9. 工程实践建议与未来方向
9.1 不要从一开始就训练端到端大模型
复现这类研究,最容易犯的错误是一上来就设计一个很复杂的网络,然后因为训练时间过长而失去耐心。更稳妥的顺序是:
- 先实现不需要嵌入的基线模型,跑通训练和评估。
- 再设计一个简单的嵌入编码器,用预训练权重或自监督方法冻结提取特征。
- 验证嵌入特征确实带来提升后,再考虑端到端微调和多尺度融合。
这样每一阶段的结果都可控,容易定位问题。
9.2 建立数据版本与日志
气象数据和遥感数据非常庞杂,实验组合很多。建议给每个数据集建立版本,记录来源、时间范围、处理脚本和校验和。训练时用类似 wandb 或 mlflow 的工具记录每次实验的配置、指标和预测样本。否则过了两周,你可能完全忘了某个实验的嵌入是用哪颗卫星的哪个波段组合生成的。
9.3 重视空间泛化能力
在实际业务中,模型要在没有站点或没有高分辨率观测的地区使用。如果训练数据只集中在某几个区域,模型对完全陌生的区域效果可能很差。最好采用空间交叉验证,把研究区分成互相不重叠的块,轮流训练和验证。这种方法比随机采样更接近真实应用场景。
9.4 安全与可信度边界
天气降尺度模型无论做得多好,都不能替代官方气象机构的预报产品。如果要把这类方法用于工程决策,比如风机偏航控制、电力负荷预测、农业保险定价,要明确模型的输入依赖、失效模式和数据延迟。概率输出不是借口,模型必须经过充分的回算测试,尤其要检查极端事件的预测能力。模型可能在小尺度空间结构上提供非常有价值的信息,但在极端天气条件下,模式本身的可靠性仍是第一位的。
9.5 未来方向:基础模型与嵌入复用
地球观测领域正朝着大模型方向发展,一个在大规模卫星影像上预训练的编码器,可以同时服务于植被分类、水体检测、天气降尺度等多个任务。未来更有效的做法很可能不是为单个任务重新训练嵌入编码器,而是复用一套通用基础模型,逐任务微调。该论文标题强调“嵌入有效”,也说明通用嵌入作为亚网格描述符具有跨场景迁移的潜力。
另一个值得探索的方向是嵌入的可解释性。如果我们能知道嵌入向量的哪些维度代表了云结构、哪些维度代表了地表湿度,就能更好地与气象学知识结合,也能让模型预测更可信。这类工作可能需要结合归因分析或概念瓶颈模型。
如果你正在做气象数据与深度学习结合的工程项目,可以沿着本文思路,先在自己的数据上比较“有嵌入”和“无嵌入”两个版本的差异。不要被论文标题里的高级概念吓到,拆开来看,它不过是给降尺度模型加了一个从地球观测数据中学习的信息源,然后用概率化输出把不确定性表达出来。有了这个框架,后续要做的就是在特征融合、预训练策略和评估指标上做细致打磨了。