基于时空深度学习模型的数值降水预报后处理

郑超昊，尹志伟，曾钢锋，许月萍，周鹏，刘莉

(1.浙江大学 建筑工程学院，浙江 杭州 310058；2.台州市水利局 浙江 台州 318000；3.台州市水利水电勘测设计院有限公司 浙江 台州 318000)

降水是水循环的重要组成部分，也是地表水文过程的基本驱动元素[1].降水预报是进行流域水文预报、水资源规划以及防洪防汛等研究的基础资料.降水时空差异性大，因此如何获取高精度、高分辨率的降水预报是水文领域一大难题[2].确定性降水预报无法完全模拟大气过程和初始场资料观测误差的影响，降雨预报的模式逐渐从确定性预测转向多成员集合预报，以集合的形式来描述未来降雨发生的概率情况[3].由于初值扰动、模式设计、集合成员数等方面的不完善因素，导致集合预报产品存在预报误差、集合成员离散度低的问题[4]而无法直接使用，必须运用后处理方法提升集合预报的可靠性和准确性，才能保证高质量的降水预报信息供给[5].

在对数值降水预报的统计后处理方法研究上，代刊等[6]系统介绍了定量降水集合预报的统计学后处理技术；Zhao等[7]评估了分位数映射方法对季节性降雨预测的校正；Andrew等[8-9]发现，参数化后处理方法贝叶斯联合概率建模的校正效果较全面.深度学习方法的学习和泛化能力出色[10]，具有非线性关系数据处理能力和数据稀疏地区较好表现的独特优势，在偏差纠正后处理研究中被广泛应用[11]：Qiu等[12]以多任务卷积神经网络考虑时间序列和站点之间的相关关系来进行降水预报，Krasnopolsky等[13]使用人工神经网络(artificial neural network, ANN)对降水集合预报进行后处理，皇甫江等[14]利用深度学习开展偏振雷达定量降水估测研究.在降水预报后处理的应用研究中，CNN和LSTM模型鲜有提及.

本研究提出结合CNN和LSTM的时空深度后处理模型(CNN-LSTM).该模型具有对时间序列数据进行分析的能力[15]，能够有效提取数据的空间信息特征特点，还能够同时考虑预报降水数据空间特征和时间相关性.本研究以浙江省椒江流域为对象，将不同预报时效的欧洲中期天气预报中心（European centre for medium-range weather forecasts,ECMWF）数值预报网格数据作为输入数据，以CMA-CMORPH实测雨量数据作为参照，对后处理前后降水资料进行全面评估.

1 研究区域与数据

1.1 研究区域概况

椒江流域地处浙江中东部，东临东海，介于东经120°21′～121°47′，北纬28°30′～29°31′之间，主要位于台州市境内，属于亚热带季风气候，多年平均气温为17.1 ℃，多年平均降雨量为1 634.9 mm.该区域的降水具有明显的季节性，全年降水主要集中在4—6月梅雨季和7—10月台风期；径流量年内分配主要集中在4—10月.永安溪和始丰溪是位于椒江上游的2条最主要河流[16].永安溪为椒江的源头，长144 km，流域面积为2 704 km2.始丰溪是椒江流域二级支流，为椒江水系最大支流，主流长132.7 km，流域面积为1 610 km2.如图1为椒江流域所处浙江省位置分布图，以及永安溪、始丰溪流域位置以及水系分布图.

图1 椒江流域的流域分布、水系、降水格点数据示意图Fig.1 Schematic diagram of basin distribution, water system and precipitation grid data in Jiaojiang river basin

1.2 研究数据

1.2.1 ECMWF降水预报数据 使用的预报数据为ECMWF数值降水预报产品，从ECMWF的官网（www.ecmwf.int）下载得到.该数据由51个预报成员组成，空间分辨率为0.5°×0.5°，选取协调世界时零点（北京时间08:00）起预报的未来0～72 h逐6 h的降水预报进行后处理.考虑到降水预报可能存在的落区误差，选取的经纬度范围为东经120°～122°，北纬27°30′～29°30′，该范围完全覆盖椒江流域及其周边区域.以2008—2018年汛期4—10月逐6 h降水资料作为输入数据，椒江流域在2012、2015和2016年遭受典型台风事件，并且这些年份的降水数据完整性高，无需插补，故选定2012、2015和2016年作为测试期，其余时间的数据作为率定期.

1.2.2 CMA-CMORPH融合数据 CMA-CMOPRH融合降水数据作为验证数据，从国家气象信息中心下载得到，为中国30 000个自动站与CMORPH降水产品融合的逐时降水量网格数据集，拥有2008—2018年分辨率为0.1°的逐小时降水数据.融合数据能够有效减少系统偏差和随机误差，CMA-CMOPRH融合降水数据总体误差水平不超过10%，优于同类型产品在研究区域的精度[17].椒江流域47个实测雨量基本站点逐小时降水数据不连续且不充分，不足以支持深度学习模型的建模，而CMA-CMOPRH融合降水数据精度较高，数据较长，因此利用该降水产品替代实测数据对ECMWF数据进行后处理.

2 研究方法

2.1 数据处理方法

在CMA-CMORPH融合数据中，选取东经120°30′～122°，北纬28°～29°30′，2008—2018年汛期4—10月逐小时降水格点数据.将数据进行6 h累计计算（北京时间8:00点作为起始），筛出数据中的异常值和缺失值，利用min-max标准化方法进行数据归一化，将6 h累计数据按照日期排序.将ECMWF数值预报降水数据转化成0～72 h预报时效的6 h累计降水预报.筛出各预报时效内数据的异常值和缺失值，利用min-max标准化对51个成员的各数据进行归一化处理，再将处理后的数据按照日期排序，以保证数据的一致性.单个成员无法很好地代表ECMWF数据的整体性，而集合平均（等权重算术平均）数据可以涵盖所有成员的预报能力，因此采用51个成员集合平均数据作为指标评判.将集合平均数据进行降尺度，使其分辨率为0.1°.利用实测降水资料的空间分布进行降尺度可以更好地捕捉降水的空间特征，但实际预报应用中不存在实测降水资料，因此本研究选用反距离权重插值法（IDW）进行数据降尺度，得到与CMA-CMORPH数据具有相同空间尺度的ECMWF数据用于后续评估.

2.2 结合CNN和LSTM的降水时空后处理模型

将CMA-CMORPH栅格数据和ECMWF降水预报数据作为输入数据，CNN提取分析空间特性，LSTM获取时间相关性[18]，CNN-LSTM进行预报降水数据的时空误差后处理.

2.2.1 模型基本原理 CNN基础构成包含3种类型的层级：卷积层、池化层以及全连接层.卷积层通过卷积运算来提取输入数据中的特征[18].在池化层中，滤波器将卷积层提取的特征信息进行压缩，并且保留提取到的最显著特征.CNN可以通过多层网络学习不同感知域对目标网格和周围区域的空间依赖性[19].LSTM是特殊的循环神经网络（RNN），它时序前馈连接的特性，使过去时刻的信息可以对当前时刻的输出产生影响.相较于传统ANN，LSTM具有4种结构：更新门、输入门、遗忘门和输出门.更新门用于储存过去信息的累计；遗忘门可以选择将过去信息中影响较小的因素“遗忘”，使模型累计的系统误差减少，得到比传统RNN更优的结果[19].

如图2所示为CNN-LSTM的模型结构，该模型是时间序列数据集从输入到结果输出的降水后处理过程.

图2 CNN-LSTM的模型结构Fig.2 Model structure of CNN-LSTM

2.2.2 模型结构 CNN-LSTM包括CNN空间特征提取、LSTM时间相关性分析2个模块.每个空间位置的降水不孤立存在，而是与周围空间相关联，利用CNN模块提取目标周围区域的降水信息，为当前目标点的修正提供参考.过去时间的降水信息会影响当前时刻降水，为此，利用LSTM模块时序前馈的特点提取降水前后方向的时间相关性.由此可知，CNN-LSTM能够兼顾时间、空间，可以运用于降水后处理中.

从CNN模块来看，由于预报数据的不确定性会导致单独成员的预报结果有较大的偏差，将51个预报成员的数据融合，生成数据集{51, 5, 5}，并将率定期所有数据组连接成连续时间序列数据集作为耦合模型的输入数据，即该数据集包括不同成员的预报信息要素，还能包含空间和时间信息.输入数据先通过二维卷积层Conv2D，再通过最大值池化层MaxPool提取最大值特征，利用Dropout方法防止模型过拟合.先将数据根据设定参数进行卷积、最大值池化，再向量平铺成一维向量.将CNN模块输出得到的向量与神经节点为256个的全连接层相连，再将全连接层转化成（16×16）的中间结果输入LSTM模块.从LSTM模块来看，16×16的二维中间向量，分别经过神经元结点为64、32、16的双向LSTM模型，每一步都通过Dropout方法，再由激活函数Relu输出16×16的二维向量与实测数据进行训练.

本研究模型的优化器采用Adam，训练迭代次数为1 000，每组训练的数据量为64，耦合模型的学习率最初为0.001，每迭代400次其学习率变为原来的0.1，模型损失函数为均方误差MSE、两幅图像相似度的指标SSIM两者加权的结果，权重设置分别为0.7、0.3.

式中：x、y均为对比图像，µx、µy分别为x、y像素的平均值，σx、σy分别为x、y像素方差，σxy为x、y的像素协方差，c1、c2均为稳定常数.

2.3 模型评估方法

以下4个性能指标用于定量评估CNN-LSTM的后处理能力：均方根误差RMSE、相对偏差RB、相关系数CC、平均绝对误差MAE.指标的计算式[20]分别为

式中：n为各个网格点观测值的数量，为网格预测降水，为网格实测观察降水，为网格预测降水的平均值，为网格实测观察降水的平均值.

为了进行分级检验评估，在ECMWF降水预报资料中提取24 h累计降水，计算永安溪、始丰溪、平原区域面雨量.根据《中短期天气预报质量检验办法》降水等级划分如表1所示，将降水预报和降水实况按照24 h累计降水量划分成对应等级.对各降水分级进行检验，各等级TS（threat score）评分为

表1 降水等级划分表Tab.1 Classification of precipitation levelsmm

式中：k为1～4级，分别代表小雨、中雨、大雨、暴雨及以上降水预报；NAk为k级降水预报正确次数；NBk为k级降水空报次数；NCk为k级降水漏报次数.ECMWF预报数据采用集合平均的数据进行TS评分计算.

3 结 果

评估验证期内不同预报时效的预报后处理数据，对比分析不同预报时效内永安溪、始丰溪、平原区域面雨量；选取验证期内的多场典型暴雨进行分析，着重分析永安溪、始丰溪、平原区域面雨量.

3.1 不同预报时效后处理效果评估

计算12个预报时效内每个网格点的数据，对比后处理前、后预报数据的指标来评估CNNLSTM.流域全部格点不同预报时效的各指标计算结果如图3(小提琴图)所示，其中虚线为25%、50%、75%分位线，t为预报时效.原始预报RMSE普遍超过3 mm，最大值超过6 mm；经CNNLSTM后处理，各个预报时效的RMSE均有所降低，基本不超过2 mm，最大值为3 mm.经CNNLSTM后处理，RB绝大部分落在-0.25附近区域，而原始预报的RB跨度大.可以看出，原始预报并不是预报时效越短预报技能越高，且预报技能也不随着预报时效增长呈基本的单调增减，这是气象预报领域较为常见的“预报不一致性”现象，即在连续多次预报中，前后两次预报的预报结果差异较大.该现象被认为主要与初始误差和模式误差有关，相关分析可参考文献[21].

图3 模型评估指标在所有格点不同预报时效后处理前、后的对比图Fig.3 Comparison of model assessment indexes before and after post-processing for all grids with different lead times

如图4所示为不同预报时效后处理前、后各指标的空间分布图.可以看出，不同预报时效后处理后的RMSE、RB、MAE、CC在全流域均有很大的提升.其中RMSE从2.6～6.0 mm下降为1.0～3.0 mm，后处理前、后RMSE下降的格点比例为100%；RB从-0.6～1.0稳定至-0.35～-0.15；MAE从0.7～1.8 mm下降为0.3～0.7 mm，CC从0.2～0.5上升至0.70～0.90；后处理前、后MAE、CC指标改善的格点比例为100%.综合所述，CNN-LSTM对于数值预报降水产品不同预报时效的预报精度均有显著提升，并且对于整个椒江流域空间精度提高具有普遍性.

图4 不同预报时效后处理前、后模型评估指标的空间分布图Fig.4 Spatial distribution of model assessment indexes before and after post-processing for different lead times

如图5所示，永安溪、始丰溪流域以及平原区域不同预见期面雨量结果相似.图中，下标数字0代表后处理前，下标数字1代表处理后.各个预报时效的RMSE均落在基准线下侧，说明后处理后数据的RMSE均小于后处理前数据的RMSE，且后处理前RMSE为2.2～4.0 mm，处理后RMSE为1.0～2.0 mm.各个预报时效的CC均落在基准线上方，后处理前的CC为0.0～0.4，处理后的CC＞0.8.各预报时效指标的MAE均落在基准线下方，后处理前的MAE为0.75～1.75 mm，处理后的MAE为0.25～0.75 mm.各个预报时效的RB分布在基准线两侧，其中后处理后的RB稳定在-0.25附近，后处理前的RB为-0.75～0.25，可见处理后的RB稳定在更合理的范围.综上所述，在永安溪、始丰溪流域，4个模型评估指标经过CNN-LSTM模型后处理后的预报面雨量数据较于后处理之前在数据精度方面有很大程度的提升.

图5 模型评估指标在永安溪、始丰溪流域面雨量不同预报时效后处理前、后的对比图Fig.5 Comparison of model assessment indexes of areal rainfall before and after post-processing over Yonganxi and Shifengxi watershed for different lead times

3.2 分等级降雨评估

如表2所示为计算所得永安溪、始丰溪、平原区域的各量级降水出现的累计频次N.如表3所示为ECMWF降水预报后处理前、后永安溪、始丰溪、平原区域的TS评分.由表可知，经过后处理3个区域的TS分别为0.900、0.906、0.894，均高于后处理之前的0.754、0.755、0.752，分别提升了19.4%、20%和18.9%.结果证实经CNN-LSTM后处理的雨季预报正确率有显著提升.

表2 永安溪、始丰溪、平原区域降水分级累计出现频次Tab.2 Cumulative occurrence frequency of precipitation classification in Yongan, Shifeng and plain regions

表3 ECMWF预报产品后处理前、后雨季3个区域的TS评分项Tab.3 TS scores items of ECMWF forecast products before and after post-processing at three regions during rainy season

如图6所示为各区域不同降水级别的TS评分图.后处理后的各等级累计降水TS评分较之后处理前的均有很大程度提升，尤其是小雨级别的TS，数值由0.77提升至0.91，中雨级别的TS从0.48提升至0.6，大雨级别的TS从0.3提升至0.4，暴雨级别的TS评分也有所提升.在后处理前、后的ECMWF预报产品在降水等级评估中，小雨的预报正确率最高，中雨次之，暴雨的预报正确率较低，这可能是暴雨事件较少导致的.

图6 ECMWF预报产品后处理前、后的等级累计降水TS评分对比图Fig.6 Comparison of TS scores of ECMWF forecast products before and after post-processing

3.3 典型暴雨后处理效果分析

选取2012年8月“海葵”台风、2015年8月“苏迪罗”台风作为典型暴雨，分析永安溪流域、始丰溪流域、平原区域未来3 d累计6 h面雨量.选取2012年8月6日8:00起未来72 h数据作为“海葵”台风的预报数据，并且利用分位线代表集合预报中95%置信区间.由图7可知(图中，p为降水量)，本场台风实际暴雨中心位于流域以北区域，流域内最大实测降水量为362.8 mm，原始预报降水中心偏离，最大雨量仅为226.6 mm，小于实测降水量；后处理后预报暴雨中心与实测降水一致，流域内累计降水量为332.6 mm.选取2015年8月7日8:00起未来72 h数据作为“苏迪罗”台风的预报数据.本场台风的暴雨中心落在椒江流域，实测累计降水的降水值为338.1 mm，原始预报降水最大雨量为308.6 mm，后处理后的预报降水最大值约为310.4 mm.由对比可知，在2场台风期间，后处理后的预报数据不管是降雨累计，还是暴雨中心都比处理前的更为精确.后处理后的预报数据对比原数据有相当程度的提升.

图7 2场台风72 h累计降水分布图Fig.7 Distribution of 72h cumulative precipitation for two typhoons

如图8所示为台风期间永安溪、始丰溪流域和平原区域逐6 h面雨量情况，为了进一步评估后处理前、后预报的差异，图中增加原始集合预报的95%的置信区间结果.在“海葵”台风期间，3个区域最大累计实测降水量分别为35、40、35 mm，预报后处理后累计降水分别为32、36、33 mm；三者预报后处理后的累计降水误差均不超过10%.在“苏迪罗”台风期间，3个区域最大累计实测降水量分别为50、60、70 mm，后处理后累计降水分别为45、54、64 mm.三者预报后处理后的累计降水误差均不超过10%，后处理前的累计降水误差较大.

图8 2场台风在3个区域的逐6 h面雨量对比图Fig.8 Comparison of 6 h areal rainfall of two typhoons at three regions

对2场台风预报进行评估，结果如表4所示.在“海葵”台风期间，3个区域后处理后的RMSE下降均超过90%，MAE在后处理后也有大幅度下降， CC在后处理后稳定在0.99，RB在后处理后提升效果均超过89%.在“苏迪罗”台风期间，3个区块后处理后的RMSE分别下降了91.3%、90.7%、91.6%；MAE分别下降了91.3%、90.2%、91.7%；后处理前的CC不足0.55，后处理后的均稳定在0.99；RB分别提升了91.3%、90.1%、91.5%.结果表明CNN-LSTM后处理的数据模拟典型台风的效果良好，面雨量模拟精度大幅度超过未经后处理的数据.结果也表明原始数值预报降水产品对台风降水预报能力较弱，通过耦合模型后处理后的降水产品对台风降水预报能力较强.

4 结 论

（1）相较于原始数值天气预报数据，后处理后的RMSE的最大值从6 mm降至3 mm；CC的最大值从0.6升至0.9，RMSE、CC的改善格点比例为100%.说明CNN-LSTM能够有效提高降水预报精度.

（2）后处理后永安溪、始丰溪流域和平原区域不同预报时效面雨量的指标有显著提升：RMSE全面提升，最大值从4 mm降至2 mm；后处理后的CC均超过0.8；MAE从0.75～1.75 mm降至0.25～0.75 mm；RB后处理后的值优于后处理前的，并稳定在-0.25.

（3）后处理后永安溪、始丰溪以及平原区域的雨季降水预报的TS评分均大于0.89，高于后处理前的TS评分；后处理后的各降水等级TS评分均好于后处理前的，其中小雨TS评分从不足0.8提升至0.91，中雨的TS评分从不足0.5提升至0.6，大雨TS评分从不足0.4提升至0.4.

（4）2场典型台风降水后处理后的累计降水量和暴雨中心评估准确性高.对于永安溪、始丰溪流域、平原区域逐6 h面雨量降水峰值的预报数据结果误差均不超过10%，各类指标评提升均超过85%，说明经过CNN-LSTM后处理的ECMWF降水产品可以用于台风降水预报.

（5）CNN-LSTM适用于数值降水预报后处理，能够提高定量降水预报的准确性，对防汛减灾都具有十分重要的意义.本研究仅针对南方湿润地区得出结论，后续将针对北方地区开展相关研究.