天气模型排名:从预报对比到排行榜的Python实现

天气模型排名:从预报对比到排行榜的Python实现 天气模型排名指的是把不同数值天气预报模型放到同一批历史天气场景里用实际观测结果检验它们的预报误差最后按误差大小排出先后。这个做法在气象服务、能源调度、农业决策、物流规划和量化交易场景里非常有用业务方不需要关心模型内部有多少物理过程只需要知道过去一段时间里哪个模型的预报更接近真实天气从而决定后续引用哪条预报源。标题里的 Show HN: Ranking weather models by how their forecasts turned out 对应的正是这样一个系统把“事后验证”做成排行榜让模型质量可比较、可追踪、可复盘。下面会围绕如何搭建一套最小可用的天气模型排名系统展开。主线是从数据准备、对齐逻辑、评分指标到 Python 实现、异常排查和工程化扩展。所有代码和数据格式都是演示结构用于说明思路真实项目接入任何模型输出、观测站点或数据文件时需要按实际字段和授权调整。1. 天气预报模型排名要解决什么问题1.1 什么是天气模型排名为什么不能只信模型说明数值天气预报模型会用物理方程描述大气运动再通过超级计算机求解得到未来若干小时甚至十几天的气象要素预报。不同模型在物理方案、分辨率、资料同化方式上差别很大所以同一天、同一个地点、同一个变量的预报结果经常不一致。使用者最关心的问题是哪个模型更可信天气模型排名就是给这个问题提供一个可量化的答案。它不评价模型内部的物理设计只评价“历史输出和真实观测的差距”。通俗地说一次预报就是一次考试观测就是标准答案。模型在大量历史样本上考出来的平均得分就是它的排名依据。技术定义可以写成给定一组模型、一组起报时间、一组站点或格点、一组气象变量和一组验证时间逐个比较预报值与观测值计算统计算法指标再按指标排序。不能只信模型说明原因主要有三点。第一模型宣传材料通常强调分辨率、同化系统等能力却不一定提供与业务场景匹配的独立验证。第二模型在不同地区、季节和天气类型下的表现可能差异很大一个总体分数不能覆盖所有场景。第三模型版本会更新排名会变化只有持续验证才能跟踪真实变化。一个榜单系统的价值不是给出一次结论而是让结论可以被持续更新和质疑。1.2 排名系统的核心链路预报、观测、对齐、评分一套最小可用的天气模型排名系统核心链路只有五步。第一步收集多个模型在过去一段时间的预报。注意“预报”不是一句话而是结构化数据至少要包含起报时间、有效时间、站点或经纬度、气象变量和预报值。第二步收集同时段的观测数据。观测可以是气象站、网格分析场或卫星反演产品但要和预报使用同一套空间和时间口径。第三步对齐。这一步最容易被忽略却最影响结果。两个模型必须比较同一个有效时间、同一个地点、同一个变量的预报才公平。第四步计算误差指标。比如绝对误差、均方根误差、偏差、降水命中率等。第五步按模型聚合指标并生成排行榜再展示给业务方或下游系统。后面所有代码和配置都围绕这条链路展开。很多人一开始就把精力放在“把模型跑起来”或“把画图做漂亮”上其实排名系统最核心的工作量在前三步数据怎么来、字段怎么对应、口径怎么统一。数据和字段错了后面的指标再专业也没有意义。1.3 适用场景与读者定位这套方案适合以下几类读者。第一类是气象数据产品开发人员需要为业务方提供模型质量看板。第二类是数据工程师正在对接气象预报和观测数据需要设计统一验证流程。第三类是研究人员想快速验证新模型或新参数方案是否优于既有模型。第四类是后端开发人员需要把历史验证结果做成排行榜页面或 API。本文的例子使用 Python 和 pandas不依赖专业气象软件。这样做的目的是先把排名逻辑讲清楚再让读者根据自己的数据格式替换加载层。如果你已经在使用 xarray、cfgrib 或 GRIB2 文件只需要替换数据加载部分评分和排名逻辑可以原样复用。2. 数据准备历史预报和观测数据如何配对2.1 先统一时间口径预报发布时间与有效时间气象数据里至少有三个时间概念起报时间issue_time、预报时效lead_hour和有效时间valid_time。起报时间是模型开始计算的时刻通常是 00 时或 12 时预报时效是预测未来多少个小时后的大气状态有效时间是预报值对应的真实时刻它等于起报时间加上预报时效。排名系统比较的是“对同一个真实时刻的预报谁更准”所以最终必须以 valid_time 作为对齐键。如果两个模型起报时间相同但一个预报 24 小时另一个预报 48 小时它们对应的 valid_time 完全不同不能放在一起比较。同样如果两个模型 valid_time 相同但起报时间不同它们利用的初始资料可能不同可以比较但要注意时效差异。下面是一个典型字段示例model,issue_time,valid_time,lead_hour GFS,2025-01-01T00:00:00Z,2025-01-02T00:00:00Z,24 ECMWF,2025-01-01T00:00:00Z,2025-01-02T00:00:00Z,24在读取数据时所有时间字段都应该转成带时区的 UTC 时间。不要使用本地时间字符串做 join因为不同数据源可能使用不同时区容易出现“看起来相同实际差 8 小时”的问题。2.2 先统一空间口径格点与站点如何对齐天气预报模型通常输出网格数据比如 0.25 度、0.5 度分辨率的经纬度格点。气象观测站是离散的点可能落在两个格点之间。为了比较需要把模型格点插值到观测站点位置或者把观测值匹配到网格上。最简单可靠的做法是使用现成插值工具将模型值插值到站点然后把站点的经纬度转换成唯一的 station_id。在演示项目中为了聚焦排名逻辑可以先把数据整理成站点的 long 表。每条记录都包含 model、valid_time、station_id、variable、value。对齐时forecast 和 observation 都通过 valid_time、station_id、variable 三个字段连接。这样能避免插值实现干扰主流程。如果数据源来自网格加载时要注意经纬度表示方式和单位。部分数据用整数经纬度部分用浮点经纬度有的用 lon 表示 -180 到 180有的用 0 到 360。做空间对齐前需要先统一经纬度范围。2.3 用最小 CSV 结构承载评分数据为了能直接运行后面的 Python 代码这里定义两个最小 CSV 文件。第一个是历史预报文件 forecasts.csvmodel,issue_time,valid_time,station_id,variable,value GFS,2025-01-01T00:00:00Z,2025-01-02T00:00:00Z,S001,temperature_2m,4.2 ECMWF,2025-01-01T00:00:00Z,2025-01-02T00:00:00Z,S001,temperature_2m,5.0 ICON,2025-01-01T00:00:00Z,2025-01-02T00:00:00Z,S001,temperature_2m,3.9第二个是观测文件 observations.csvvalid_time,station_id,variable,value 2025-01-02T00:00:00Z,S001,temperature_2m,3.5这两个文件都很小但已经包含了排名系统需要的所有核心字段。变量名建议统一使用机器可读的小写字符例如 temperature_2m、wind_speed_10m、precipitation_24h。如果值有单位最好在表里增加 unit 字段或者先在预处理阶段统一单位避免不同模型输出摄氏度和华氏度造成误差计算错误。注意观测文件的粒度是“真实时刻的客观值”而预报文件还需要保留 model 和 issue_time。合并时只保留 valid_time 相同的记录否则样本对不齐排名会失真。2.4 数据质量检查清单在开始写评分代码前建议先对数据做一轮质量检查。下面这个清单可以直接复用检查项检查方式通过标准时间字段格式打印每个来源的 unique 时间样例全部使用 ISO8601 且带 Z 时区valid_time 覆盖范围检查 min/max预报和观测覆盖同一时间段station_id 是否一致计算集合差预报和观测站点集合一致variable 是否一致计算集合差变量名完全匹配是否存在重复记录按 key 分组统计每个 key 只有一条记录值是否缺测查看 value 为 NaN 或特殊值缺测数量低于阈值且被标记这个检查可以用 pandas 在预处理阶段自动完成。不要省略这些步骤因为时间或空间对齐一旦出错排名结果会出现系统性偏差而且肉眼很难发现。3. 用 Python 实现最小可运行的模型评分与排名流程3.1 项目结构和依赖在本地目录创建一个最小项目weather-leaderboard/ ├── data/ │ ├── forecasts.csv │ └── observations.csv ├── leaderboard.py ├── requirements.txt └── run.sh依赖只需要 pandas 和 numpy画图可以后面再加。安装命令pip install pandas numpy matplotlibrequirements.txt 可以写成pandas2.0 numpy1.24 matplotlib3.7这一步的目的是让示例可以快速运行。真实项目中如果使用 GRIB2 或 NetCDF需要额外引入 xarray、cfgrib 或 h5netcdf但核心评分逻辑不受影响。3.2 加载和对齐数据在 leaderboard.py 中先读取两个 CSV 文件import pandas as pd def load_data(forecast_path, observation_path): forecasts pd.read_csv(forecast_path, parse_dates[issue_time, valid_time]) observations pd.read_csv(observation_path, parse_dates[valid_time]) return forecasts, observations然后写对齐函数def align_forecast_observation(forecasts, observations): merged pd.merge( forecasts, observations, on[valid_time, station_id, variable], suffixes(_fc, _obs), howinner ) merged[error] merged[value_fc] - merged[value_obs] return merged这里使用内连接强制保留预报和观测都存在的样本。如果某个模型出现了数据缺口它的样本量会比其他模型小这时排名结果需要谨慎解读。使用 suffixes 区分预报值和观测值error 表示预报减观测的差值正偏差表示预报偏高。对齐之后需要检查一下样本数量def sample_summary(df): return df.groupby(model).size().reset_index(namesample_size)3.3 计算误差指标接下来定义三个连续变量误差指标平均绝对误差 MAE、均方根误差 RMSE、平均偏差 Bias。import numpy as np def mae(obs, fct): obs np.asarray(obs, dtypefloat) fct np.asarray(fct, dtypefloat) return float(np.mean(np.abs(fct - obs))) def rmse(obs, fct): obs np.asarray(obs, dtypefloat) fct np.asarray(fct, dtypefloat) return float(np.sqrt(np.mean((fct - obs) ** 2))) def bias(obs, fct): obs np.asarray(obs, dtypefloat) fct np.asarray(fct, dtypefloat) return float(np.mean(fct - obs))这三个函数都比较简单。MAE 直观反映平均误差大小RMSE 对大误差更敏感因为平方项放大了离群值的影响Bias 表示系统性误差正值偏暖或偏强负值偏冷或偏弱。排名时通常以 MAE 或 RMSE 为主指标Bias 作为辅助参考。3.4 按模型生成排行榜核心逻辑是按 model 分组对每组样本计算指标再排序def rank_models(merged, metricmae): records [] for model, group in merged.groupby(model): obs group[value_obs] fct group[value_fc] records.append({ model: model, sample_size: len(group), mae: mae(obs, fct), rmse: rmse(obs, fct), bias: bias(obs, fct), }) result pd.DataFrame(records) # 这里按误差指标从小到大排名bias 不作为主排名指标 result[rank] result[metric].rank(methodmin, ascendingTrue).astype(int) result result.sort_values([rank, model]).reset_index(dropTrue) return result这段代码先遍历每个模型再生成 DataFrame。metric 参数可以切换排名用的指标默认 mae。rank 列使用 pandas 的 rank 方法相同误差时并列名次。这里要注意样本量暂时只作为展示不参与排序。如果两个模型样本量差异很大需要在业务层决定是否过滤。主流程如下if __name__ __main__: forecasts, observations load_data(data/forecasts.csv, data/observations.csv) merged align_forecast_observation(forecasts, observations) if merged.empty: raise SystemExit(没有匹配到有效样本请检查时间、站点、变量字段) table rank_models(merged, metricmae) print(table.to_string(indexFalse))3.5 运行结果示例如果 forecasts.csv 中包含多个站点和多天数据运行后可能得到类似这样的输出model sample_size mae rmse bias rank ECMWF 2800 1.82 2.41 0.35 1 GFS 2800 2.05 2.78 0.10 2 ICON 2800 2.31 3.02 -0.42 3这个结果的业务