机器学习干旱预测实战：ml_drought项目解析-平芜编程栈

简介：这套代码是一套面向气候科学的端到端机器学习管道，源自 GitHub 上的 ml-clim/ml_drought 项目，核心目标是借助机器学习提高干旱预测能力并加深对干旱机制的理解。管道将数据标准化、训练集划分、模型训练、结果比较与评估整合为统一流程，各步骤按模块组织在 src 目录下，并用环境配置文件锁定依赖；使用者只需具备基础 Python 知识，即可沿管道完成从数据处理到模型评估的闭环实验，适合气候研究人员和希望将机器学习落地到气象场景的开发者使用。资源包为 zip 压缩格式，大小约 49.31MB，页面未单独列出文件总数与类型明细；从项目说明来看，包内包含 environment.yml 环境配置文件、Python 源码和用于演示的 Jupyter Notebook 样例，可通过 conda 一键创建 esowc-drought 环境并复现实验结果。已有 142 人学习。读者可以得到一套可复用的干旱预测工作流，学习如何模块化组织数据预处理、建模与评估环节，并方便地迁移到其他气候变量的机器学习研究中；项目结构和模块划分清晰，也适合作为后续扩展与二次开发的参考实现。

1. ml_drought 到底是什么：一个能复现的机器学习干旱预测项目

第一次打开ml_drought这个仓库时，我的第一反应是"又一个用机器学习打包的 demo"——直到我顺着它的 Notebook 跑完一遍，才发现它把干旱预测这件事拆成了完整的链路：从 Copernicus 的 ERA5 再分析数据下载，到标准化降水蒸散指数（SPEI）计算，再到随机森林训练与时间序列验证。项目挂在 github.com 上的ml-clim组织下，2019 年的代码，Jupyter Notebook 为主，不 fancy，但管线干净。它要解决的核心问题是：不靠物理数值模式，仅用机器学习对再分析数据做统计建模，能不能提前几个月识别干旱前兆。复现成本很低，一台普通笔记本就能跑通，适合刚接触气候数据、又想上手机器学习完整流程的人，也适合想知道"机器学习在气象里到底怎么落地"的从业者。这一篇我按拆过的资源来写，尽量把每一步都说得能直接操作。

2. 数据管线：Copernicus 再分析数据怎么变成训练集

2.1 你需要先装齐的工具清单

ml_drought 的 Notebook 依赖不算重，核心是xarray、netCDF4、pandas、scikit-learn，外加一个下载 ERA5 数据用的cdsapi。我第一次复现时直接pip install一把梭，结果版本冲突卡了半小时，后来固定了环境才顺畅。建议先建一个干净的虚拟环境再装：

python -m venv ml_drought_env source ml_drought_env/bin/activate # Windows 下用 ml_drought_env\Scripts\activate pip install pandas numpy xarray netCDF4 scikit-learn matplotlib pip install cdsapi # Copernicus CDS 下载专用

cdsapi是访问 Copernicus Climate Data Store 的官方 Python 客户端，没有它你没法自动拉 ERA5 数据；xarray负责读取 NetCDF 格式的再分析场，netCDF4是它的底层引擎。装完记得验证一下导入是否正常，缺底层库时xarray导入会报ImportError，这时候按报错信息补包就行。

2.2 CDS API 配置：下载 ERA5 最常见的卡点

下载 ERA5 不需要人工点网页，但需要先在 CDS 网站注册账号、拿到 API Key，然后写入本地配置文件。这个 Key 在~/.cdsapirc里，格式固定：

url: https://cds.climate.copernicus.eu/api key: <你的UID>:<你的API Key>

我一般会在下载脚本里加一个重试机制，因为 ERA5 的下载队列经常排队，单次请求失败很常见。下面的脚本是 ml_drought 这类项目里通用的取数写法，按月度降水、2m 气温、蒸散量三个变量拉 1979 年至今的数据：

import cdsapi c = cdsapi.Client() c.retrieve( "reanalysis-era5-single-levels-monthly-means", { "product_type": "monthly_averaged_reanalysis", "variable": [ "total_precipitation", "2m_temperature", "potential_evaporation" ], "year": [str(y) for y in range(1979, 2020)], "month": [f"{m:02d}" for m in range(1, 13)], "time": "00:00", "format": "netcdf", }, "era5_monthly.nc", )

这段代码的意图很直白：请求月度平均的 ERA5 单层数据，把三个变量一次性打包成一个 NetCDF 文件。variable列表里的三项决定了你后续能算什么指数——total_precipitation是降水项，potential_evaporation用来算蒸散，2m_temperature是温度项。year和month用列表推导式生成，比手写长字符串省事。如果网络环境不稳，建议外面再套一层for循环，按年逐次请求，失败则等 30 秒重试，避免整段任务因为一次超时全部作废。

2.3 把 NetCDF 读成表格：区域裁剪与时间维处理

ERA5 下载下来是全球 0.25° 网格，直接整场建模既慢又没必要。ml_drought 的做法是先裁剪目标区域，再聚合成以月为单位的 DataFrame。xarray的sel方法可以按经纬度切片，ds.mean(dim=...)可以算区域平均。这一步的细节决定后面特征矩阵长什么样：

import xarray as xr import pandas as pd ds = xr.open_dataset("era5_monthly.nc") # 裁剪地中海区域，lat 在前，lon 在后是 xarray 的默认维度顺序 region = ds.sel(lat=slice(30, 45), lon=slice(-10, 40)) # 按区域平均，把 3D 场压缩成时间序列 df = region.mean(dim=["lat", "lon"]).to_dataframe().reset_index() # 构造时间索引 df["time"] = pd.to_datetime(df["time"]) df = df.set_index("time").sort_index() df.to_csv("mediterranean_era5.csv")

sel(lat=slice(30, 45), lon=slice(-10, 40))是区域裁剪的常用写法，注意slice的顺序是"上界在前、下界在后"，反过来会得到空数据。mean(dim=["lat", "lon"])把整个区域的场平均成单条时间序列，这个操作很粗暴，但作为快速原型够用；如果后续要做空间分布预测，得保留经纬度维度，改成按格点单独建模。最后.to_dataframe()把xarray.Dataset拍平成pandas.DataFrame，time变量会自动变成索引。存成 CSV 之后，机器学习部分就跟普通表格数据没什么区别了。

注意：ERA5 的time变量在月度平均数据里是每月 1 号的 00:00，下游算滞后特征时直接用shift即可，不需要重采样对齐。

3. 特征工程与模型选型：干旱标签怎么来，模型怎么挑

3.1 SPEI 计算：把连续气象变量变成干旱等级

干旱不是一个可以直接测量的量，气象上常用标准化降水蒸散指数（SPEI）来定义"干不干、干到什么程度"。ml_drought 的做法是把降水与潜在蒸散的差值拟合成概率分布，再映射成标准正态分位数。SPEI 的经典范围是 -2 到 +2，正值偏湿，负值偏干。计算的核心是对差值序列做 gamma 分布拟合，scipy.stats.gamma足够处理：

from scipy.stats import gamma import numpy as np # 用降水与蒸散差构造水分盈余序列 df["water_balance"] = df["total_precipitation"] - df["potential_evaporation"] # 逐月拟合 gamma 分布，转换为标准化分位数 data = df["water_balance"].values # gamma 分布要求正值，做偏移保证输入为正 offset = -np.min(data) + 0.1 fit_data = data + offset shape, loc, scale = gamma.fit(fit_data) spi_values = gamma.cdf(fit_data, shape, loc=loc, scale=scale) df["spei"] = norm.ppf(spi_values) # 标准正态逆变换

这是 SPEI 的一种简化实现，正统做法还要按月份分组拟合，因为 1 月和 7 月的水分平衡分布并不同。gamma.fit返回三个参数：shape控制分布形态，loc是位置偏移，scale是尺度，这三者直接决定了 CDF 转换的结果。norm.ppf把累积概率映射到标准正态分位数，得到最终的 SPEI。实战中我更推荐直接用spei库或climate_indices库，但手动算一遍能帮你理解标签是怎么来的。得到 SPEI 后按阈值切分成类别，常见的分级如下：

SPEI 范围	分类	是否算干旱事件
>= 0.5	偏湿	否
-0.5 ~ 0.5	正常	否
-1.0 ~ -0.5	轻度干旱	弱干旱
-1.5 ~ -1.0	中度干旱	是
<= -1.5	重度干旱	是

3.2 滞后特征：让模型"看见过去"而不是猜未来

机器学模型不会自动理解"时间"，你要把前几个月的降水、气温、SPEI 显式地拼到当前样本里，模型才能学到"连续三个月降水偏低之后，第四个月出现干旱的概率上升"这类规则。pandas.shift是构造滞后特征最直接的工具：

# 构造 1、3、6 个月的滞后特征 for lag in [1, 2, 3, 6]: df[f"spei_lag{lag}"] = df["spei"].shift(lag) df[f"tp_lag{lag}"] = df["total_precipitation"].shift(lag) df[f"t2m_lag{lag}"] = df["2m_temperature"].shift(lag) # 删掉前 6 个月没有完整滞后特征的样本 df = df.dropna()

shift(lag)把整列下移，lag=1表示上个月的值，lag=6表示半年前的值。选择 1、2、3、6 这几个滞后窗口是气候预测里的常见做法：1 到 3 个月捕捉短期惯性，6 个月捕捉季节尺度信号。dropna()必不可少，否则前 6 行样本的滞后特征全是NaN，模型训练时会报错或者静默丢弃。你也可以加入滚动均值，比如rolling(3).mean()，但要注意滚动窗口和滞后长度的信息重复问题，加太多冗余特征会让模型过拟合。

3.3 为什么从随机森林和逻辑回归开始

ml_drought 仓库里最早用的是随机森林和逻辑回归，这是 2019 年前后气候机器学习项目的标准搭配。逻辑回归的优势是系数可以直接解释成"某个滞后特征每增加一个单位，干旱概率变化多少"，非常适合理解机制；随机森林则能捕捉非线性关系，对特征缩放不敏感，几乎不需要调参就能出基线。它们的训练成本都极低，在几千条样本上秒级完成，适合做基准对比：

from sklearn.ensemble import RandomForestClassifier from sklearn.linear_model import LogisticRegression feature_cols = [c for c in df.columns if "lag" in c] X = df[feature_cols] y = (df["spei"] <= -1.0).astype(int) # 1 表示干旱事件 models = { "logistic": LogisticRegression(max_iter=1000, class_weight="balanced"), "random_forest": RandomForestClassifier( n_estimators=300, max_depth=6, min_samples_leaf=10, class_weight="balanced", random_state=42, ), }

y的构造方式影响建模目标：把 SPEI 小于等于 -1.0 视为一次干旱事件，这是二分类；如果你想要分级预测，改成三分类或多分类即可。class_weight="balanced"是针对干旱样本偏少的常用手段，让模型在训练时给少数类更高的权重，逻辑回归里它等价于调整决策阈值，随机森林里它影响叶子节点的分裂权重。max_depth=6和min_samples_leaf=10是随手设置的保守参数，目的是防止树长太深、在干湿年份交替的数据上记忆噪声。

4. 训练与验证：Jupyter Notebook 里的全流程怎么组织

4.1 时间序列交叉验证：随机 shuffle 是第一个大坑

气候数据按时间排序，相邻月份高度相关。如果直接用train_test_split(shuffle=True)，训练集里混着未来的样本，模型会"偷看"答案，验证分数虚高。ml_drought 这类项目必须用时间序列交叉验证，常见做法是按年份切块，前 N 年训练、后 1 年验证，然后滚动：

import numpy as np from sklearn.model_selection import TimeSeriesSplit tscv = TimeSeriesSplit(n_splits=5) for train_idx, val_idx in tscv.split(X): X_train, X_val = X.iloc[train_idx], X.iloc[val_idx] y_train, y_val = y.iloc[train_idx], y.iloc[val_idx] # 每一折都用过去预测未来

TimeSeriesSplit的机制是把样本顺序切分成递增的训练块——第一折训练集占前 20%，第二折占前 40%，以此类推，测试块永远在训练块之后。这保证验证过程不会用到未来信息。n_splits=5意味着至少需要 6 年以上数据，数据太短时减到 3。这个环节的误用是"换汤不换药"：有人用了TimeSeriesSplit，但特征里却包含了未来月份的滞后值，交叉验证形同虚设，这一点我放到第五章细说。

4.2 训练、评估与混淆矩阵：别只看准确率

干旱事件在总样本里占比通常不到 20%，准确率这个指标会被"绝大多数月份不干旱"的多数类拉高，没有参考价值。ml_drought 的 Notebook 里用的是混淆矩阵加 F1 分数，配合精确率和召回率一起看：

from sklearn.metrics import classification_report, confusion_matrix rf = models["random_forest"] rf.fit(X_train, y_train) y_pred = rf.predict(X_val) print(classification_report(y_val, y_pred, target_names=["正常", "干旱"])) print(confusion_matrix(y_val, y_pred))

classification_report输出每一类的精确率、召回率、F1，干旱类的召回率尤其重要——漏报一次干旱事件，比误报一次正常月份代价高得多。随机森林还自带特征重要性，rf.feature_importances_可以直接画条形图，看哪些滞后特征起了作用。这一步属于"结果可解释"的起点，也是项目名里 "understand drought" 的含义所在。