1 引　言

临近预报是指对0—2 h的短时间内天气现象变化的预报（Browning，1982）。天气雷达探测资料具有较高的时间和空间分辨率，成为临近预报的主要工具。传统的雷达回波外推方法主要有交叉相关算法（Tracking Radar Echoes by Correlation，TREC）（Li，et al，1995）、单体质心法（The Strom Cell Identification and Tracking，SCIT）（Witt，et al，1993）和光流法（Optical Flow）（Ayzel，et al，2019）。TREC是通过计算雷达回波资料在连续时次的空间最优相关，得到对流系统不同位置的移动矢量特征，并基于这些获得的移动矢量对雷达回波进行外推预报（陈明轩等，2007），该方法对变化较快的对流性降水的预报效果差，张亚萍等（2006）提出了一种基于差分图像的相关方法追踪雷达回波运动（Difference Image based Tracking Radar Echo by Correlations，DITREC），消除了TREC矢量场中由于回波快速变化而导致的无序矢量，在使用导出的DITRE场预报降水场时，其精度依赖于所采用的Z-I关系。SCIT主要用于雷暴单体的识别追踪，其在识别了三维雷暴单体后对单体质心路径进行追踪，适用于对强雷暴单体的追踪和临近预报（俞小鼎等，2012），但丢失了弱回波信息。光流法是通过计算雷达回波的光流场得到回波的运动矢量场，并基于运动矢量场对雷达回波进行外推。该方法与交叉相关算法不同之处在于光流法立足于变化，而不是选定不变特征（曹春燕等，2015），但是其在计算光流矢量和外推两步时存在累积误差。这些方法仅根据若干个时刻的雷达回波图像推测下一时刻的回波位置，且忽略了实际情况下雷达回波中的中小尺度大气系统的运动非线性，存在对历史雷达资料利用率不足和外推时效较短的局限。

针对上述传统方法存在的问题，引入神经网络进行雷达回波的外推，神经网络模型具有强大的非线性映射能力，且对数据的利用率大幅度提高（高隽，2003）。将神经网络应用于临近预报，目前还处于起步阶段，中外学者在这方面展开了一定的探索。冯汉中等（2004）运用支持向量机（Support Vector Machine，SVM）实现了降水量的预报。陈家慧等（2000）尝试将BP神经网络模型应用于雷达回波的临近预报中，BP模型较好的记忆联想功能展现了运用神经网络进行雷达资料外推的可行性。郭尚瓒等（2017）将多层感知器（Multilayer Perceptron，MLP）和光流法两个模式集成后实现了预测某区域36 min内的降雨概率。郭瀚阳等（2019）使用卷积GRU（Gated Recurrent Unit，门控循环单元）对雷达回波进行外推，并与传统雷达回波区域跟踪外推临近预报方法（TREC）进行效果对比，结果表明，深度学习模型更具有优势。Klein等（2015）在传统卷积神经网络（Convolutional Neural Networks，CNNs）结构上增加了动态卷积层，生成两个预测概率向量实现对降水回波的预测。施恩等（2018）在Klein等（2015）提出的动态卷积层的基础上加入循环神经网络（Recurrent Neural Network，RNN）构建了循环动态卷积神经网络（Recurrent Dynamic Convolutional Neural Networks，RDCNN），该网络在预报准确率和预报时效两方面均取得较好效果。传统的长、短时记忆单元（Long Short-Term Memory，LSTM）无法实现对空间特征的提取，针对这个问题，Shi等（2015）将输入到状态及状态到状态之间的变化改为卷积操作，提出了卷积LSTM（ConvLSTM）网络，为了适应实际情况中大部分运动的非时、空恒定，Shi等（2017）改进了模型，引入具有可学习卷积的TrajGRU（轨迹GRU）模型。Singh等（2017）也通过在循环神经网络的基础上增加卷积结构来适应雷达回波图像的时、空依赖性，实现了基于雷达回波图像序列的预测。韩丰等（2019）对采用神经网络和相关交叉算法进行雷达临近预报的长期和个例检验结果表明，神经网络在20 dBz和30 dBz检验项上的准确率显著提高。Wang等（2017）提出了基于LSTM单元改进的ST-LSTM单元（Spatiotemporal LSTM），并将其应用于一个新的端到端模型Predrnn。Agrawal等（2019）将降水预报看作由图片至图片的转换问题，并利用一种U-net结构的卷积神经网络来实现预报目的，是一种利用数据驱动、完全不使用大气物理模型来建立的短时临近降水预报模型。Sønderby等（2020）开发的模型MetNet可以预测未来8 h内的降水量。神经网络模型具有一定的通用性，雷达回波外推问题同时可以类比于视频帧画面内容的预测问题，Lotter等（2016）提出的采用CNN+LSTM结构设计的PredNet，可以很好地学习视频的特征、并可用于雷达回波的外推预报。Wang等（2019）将三维卷积和RNN集成为Eidetic 3D LSTM单元，使得该网络具有较好的对早期活动的识别能力。张德正等（2019）通过加深由卷积门控循环单元（GRU）构成的神经网络中的卷积层来提升对空间结构信息细节的表征能力，该网络不仅可以实现视频帧的预测还可以实现对雷达回波的外推预测。

本研究采用带有Causal-LSTM单元组成的神经网络模型（Wang，et al，2018）作为网络架构实现雷达回波的外推。为消除降水的不平衡和提高对强回波的预测准确率，使用带权重的损失函数训练网络模型，在训练过程中赋予强回波较大的权重来提升对强回波区域的预报能力，并通过测试集和4个个例对该预测神经网络和光流法进行对比评估。

2 方　法 2.1 循环神经网络

循环神经网络是一类以序列数据为输入，在序列的演进方向进行递归且所有节点（循环单元）按链式连接的递归神经网络（Goodfellow，et al，2016）。该网络包括前馈通路和反馈通路，由于反馈通路的存在，神经元的输出信号经过一个或几个时间步长之后能够重新作为其他神经元或自身输入，进而增强了网络处理时间序列的能力。图1为一个简单的循环神经网络结构。

图中， $x$ 是神经网络的输入， $y$ 是神经网络的输出， $h$ 是神经网络隐含层的输出，u是输出层到隐含层的权重矩阵，v是隐含层到输出层的权重矩阵，w是记忆单元到隐含层的权重矩阵。计算公式如下

${h_t} = f({\boldsymbol u} {\text •} {x_t} + {\boldsymbol w} {\text •}{h_{t - 1}})$

(1)

${y_t} = g({\boldsymbol v} {\text •} {h_t})$

(2)

式中，f、g表示激励函数， ${x_t}$ 是当前时刻的输入， ${h_{t - 1}}$ 是前一时刻隐含层的输出， ${h_t}$ 是当前时刻隐含层的输出， ${y_t}$ 是当前时刻网络的输出。

2.2 网络模型

本研究采用带有Causal LSTM单元组成的神经网络模型作为网络架构。Causal LSTM单元包含两个记忆状态：时间记忆单元（ $ {C}_{t}^{k} $ ）和空间记忆单元（ $ {M}_{t}^{k} $ ），Causal LSTM由三个级联部分组成，第1个级联用于更新时间记忆单元，第2个级联用于更新空间记忆单元，之后将更新后的 $ {C}_{t}^{k} $ 、 $ {M}_{t}^{k} $ 和输入 ${x_t}$ 共同传入第3个级联后输出该时刻的隐含层 $ {H}_{t}^{k} $ 。

与传统LSTM单元相比，其优点在于：①输入与状态以及状态与状态的转化改变为卷积操作，整个网络就可以在解决时间依赖问题的同时提取图像的空间特征；②三个级联的形式加深了网络的深度，增加更多的非线性操作使得特征放大，更有利于捕捉短期动态变化和突发情况。

除了短期的动态变化特性外，由于较长的传播，Causal LSTM会遭受梯度反向传播的困扰，时间记忆单元（ $ {C}_{t}^{k} $ ）可能会忘记过去的空间特征。为解决这一问题，网络模型中加入梯度高速公路单元（Gradient Highway Unit，GHU）解决梯度消失的问题。GHU的输入为目前的更低一层的输出和前一个时刻GHU的输入，即GHU连接了当前时刻以及前一个时刻的输入。GHU的加入使梯度不再是一条线传播，而是直接在第1层与第2层之间有个高速的传播。换言之，传播的距离缩短了。基于此，整个网络能够自适应地捕捉长期和短期的特征并且能够建立较深的网络进行学习。

2.3 损失函数

雷达回波外推的主要目的是实现降水回波特征变化的临近预报，利用深度学习进行雷达回波外推虽然取得了一些令人惊喜的结果，但仍旧处于起步阶段，行业内尚未达成对该种模型评判的统一标准。在现有的神经网络预测降水模型中，多使用单一阈值判别降雨与否（Shi，et al，2015）。不同强度的降水对人类活动产生的影响程度均不同，例如强降水出现的频率低但是影响较大，考虑到这种不平衡性，Shi等（2017）将雨强划分为6个等级，并给不同强度的降水赋予不同的权重。

在雷达定量测量降水的过程中，雷达反射率因子与雨强存在统计关系：Z=AI^b，其中I为雨强，单位mm/h，A、b为两个统计常数。Z-I关系不仅随地点、季节以及降水类型而变，即使在同一次降水过程中，A、b值也是随时间和空间而变化的。降雨预测的准确度依赖于雷达回波外推的准确度，即对降水预测模型的评估归结为对雷达反射率因子外推准确率的评估。在雷达气象学中，一般层状云降水的雷达反射率因子在20—30 dBz，对流云降水一般不低于40 dBz，因此，参考Shi等（2017）的雨强分级方案，对不同区间的雷达反射率因子赋予不同的权重，即根据每个像素点的雷达反射率因子大小分配权重。为消除降水的不平衡并提高对强回波的预报能力，对强回波赋予更大的权重以提高其在拟合过程中的影响程度，使其达到更好的收敛。

雷达图像每个像素点的回波强度（单位：dBz）的有效预测数值在[0，75 ]。将雷达图像有效预测数值划分为5个等级，第1级的取值区间为[0，20]（无雨或毛毛雨）；第2级的对应取值区间为（20，30]（小雨）；第3级的取值区间为（30，40]（中雨）；第4级的取值区间为（40，50]（大雨、暴雨）；第5级的取值区间为（50，75]（特大暴雨）。对不同级别的雷达反射率因子分配不同的权重。在此采用的权重为

$w(x) = \left\{ {\begin{array}{*{20}{l}} 1&{x \in [0,20]} \\ {10}&{x \in (20,30]} \\ {15}&{x \in (30,40]} \\ {30}&{x \in (40,50]} \\ {50}&{x \in (50,75]} \end{array}} \right.$

(3)

式中， $ w \left(x\right) $ 表示权重大小， $ x $ 表示像素点的回波数值。

本研究的预测神经网络采用的损失函数为带权重的均方误差（E_s）和带权重的平均绝对误差（E_a），表达式如下

${{{E}}_{\rm{s}}} = \frac{\,1\,}{\,n\,}\sum\limits^{N}_{n=1} {\sum\limits_{i = 1}^{400} {\sum\limits_{j = 1}^{400} {{w_{n,i,j}}} } } {({x_{n,i,j}} - {\hat x{_{n,i,j}}})^2}$

(4)

${E_{\rm a}} = \frac{\,1\,}{\,n\,}\sum\limits_{n=1}^N {\sum\limits_{i = 1}^{400} {\sum\limits_{j = 1}^{400} {{w_{n,i,j}}} } } \left| {{x_{n,i,j}} - {{\hat x}{_{n,i,j}}}} \right|$

(5)

式中，N为外推的总时间步长(每个步长是6 min)， $ {w }_{n,i,j} $ 表示第n时刻第（ $ {i},\;{j} $ ）像素点的权重， $ {x}_{n,i,j} $ 表示真实值， $ {\hat{x}}_{n,i,j} $ 表示预测值。

在实际试验中，采用的损失函数为E_s与E_a之和。

3 数据集构造 3.1 数据质量控制

试验使用福建、江苏、河南的新一代S波段多普勒天气雷达探测资料，该雷达有效探测距离为230 km，波束宽度为1.0°，反射率因子的距离库长为1 km。为确保样本数据的质量，首先对体扫反射率因子数据进行质量控制，质量控制包括滤除孤立噪声点与填补缺测值、剔除地物杂波和滤波平滑（庄薇等，2012；黄兴友等，2018；马雷，2018）三个方面。

3.2 雷达数据格点化

新一代天气雷达基数据采用极坐标的存储方式，空间分辨率不均匀，近距离的空间分辨率高，远距离的空间分辨率低。为了方便后续的研究，首先根据肖艳姣（2006）提到的方法将极坐标存储下的数据采用径向和方位上的最近邻居法与垂直线性内插法插值到三维笛卡尔坐标系下，网格水平分辨率为1 km×1 km，垂直分辨率也为1 km，生成的CAPPI图像分辨率为460×460。考虑到对流层中代表性平均引导气流的高度，其样本选用垂直高度为3 km处的CAPPI图像。

3.3 样本数据集的生成

训练时采用灰度图像可以提升训练的效率，因此将CAPPI图像处理为PNG格式存储的灰度图像，再对CAPPI灰度图像进行裁剪处理，保留图像中心的400×400区域。由于样本质量会对训练结果造成影响，面积过小的回波一般弱且演变快，对预报工作无太强指导意义，因此在样本数据集构造过程中，对个例进行筛选，将雷达有效探测范围内回波面积覆盖率大于1/10的个例选入样本数据集。样本中每个时次保存为一张PNG格式的CAPPI回波强度灰度图像，每组样本覆盖时长为4 h（输入时长2 h，预测时效2 h），每组样本中图像的时间间隔为6 min，即每组样本共有40幅雷达图像 $\{ {x_1},{x_2}, $ $ {x_3}, \cdot \cdot \cdot,{x_{20}},{y_1},{y_2},{y_3}, \cdot \cdot \cdot,{y_{20}}\} $ ，其中 $\{ {x_1},{x_2},{x_3}, \cdot \cdot \cdot,{x_{20}}\} $ 为输入图像序列， $\{ {y_1},{y_2},{y_3}, \cdot \cdot \cdot,{y_{20}}\} $ 为预测图像的标签。雷达回波的数据范围是0—75 dBz。图2为数据集构造流程。图3为数据集中一组样本可视化后的图像，x₁—x₂₀是输入图像，y₁—y₂₀是对照标签。

4 试验结果与分析 4.1 模型评价指标

对测试集的检验采用分阈值和预报时效的评估方法，使用20、30、40、50 dBz作为判别阈值，预报时长为120 min，每个预报时间步长为6 min，使用的评价指标为临界成功指数（CSI）、命中率（POD）、虚警率（FAR）。将预测值与观测值进行逐点对比，判别每个像素点是命中点（预测值与真值均大于判别阈值）、空报点（预测值大于判别阈值而真实值小于判别阈值）或漏报点（预测值小于判别阈值而真实值大于判别阈值）。a、b、c分别表示预测图像中的命中数、空报数和漏报数，三个评估指标的计算公式如下

$\begin{split}\\ \\{\rm{CSI}} = \frac{a}{{a + b + c}}\end{split}$

(6)

${\rm{POD}} = \frac{a}{{a + c}}$

(7)

${\rm{FAR}} = \frac{b}{{a + b}}$

(8)

在个例检验中采用结构相似性（SSIM）（Wang，et al，2004）来衡量两幅图像的相似度，SSIM从灰度、对比度和结构三个方面来估计降质图像的视觉质量。SSIM的计算公式如下

${\rm{SSIM}}(x,y) = \frac{{(2{\mu _x}{\mu _y} + {C_1})(2{\sigma _{xy}} + {C_2})}}{{(\mu _x^2 + \mu _y^2 + {C_1})(\sigma _x^2 + \sigma _y^2 + {C_2})}}$

(9)

式中，x、y为两张图像的像素值， $ {\mu }_{x} $ 为 $ x $ 的均值， $ {\mu }_{y} $ 为 $ y $ 的均值， $ {\sigma }_{x}^{2} $ 为 $ x $ 的方差， $ {\sigma }_{y}^{2} $ 为 $ y $ 的方差， ${\sigma _{xy}}$ 为 $ x $ 和 $ y $ 的协方差， $ {C}_{1}={\left({K}_{1}L\right)}^{2}{\text{、}}{C}_{2}={\left({K}_{2}L\right)}^{2},L $ 是像素值的动态范围，文中的变化范围是0—75， $ {K}_{1}=0.01 $ ， $ {K}_{2}=0.03 $ 。结构相似性的范围为0到1，SSIM值越大，表明两幅图越相似，当两张图一模一样时，SSIM=1。

4.2 测试集检验

试验的数据集由训练集和测试集组成，采用的试验数据集为福建2018年6—9月、2019年6—9月，江苏2018年6—9月、2019年6—9月，河南2017年5—9月、2018年6—9月、2020年6—9月SA雷达体扫资料。表1给出了光流法和预测神经网络在测试集上的指标评估结果。两个网络分别在30、60、90、120 min外推时效下对20、30、40、50 dBz这4个反射率阈值进行检验。由表1可知，在相同预报时效的情况下，4个反射率检验阈值的临界成功指数、命中率均表现为：预测神经网络优于光流法，虚警率表现为预测神经网络优于光流法。随着预报时效的延长，两个模型在相同检验阈值的情况下，均表现为临界成功指数、命中率随着时间延长而下降，虚警率随着时间延长而上升，这表明随着时间的推移，两个模型的预报准确率都在下降。从表1中也可以看出，同一模型在相同预报时效的情况下，对越强回波的预报准确率越低，但是使用带权重的损失函数训练的预测神经网络模型在40 dBz和50 dBz检验阈值项上，临界成功指数、命中率较光流法均有所提高，虚警率有所下降，表明预测神经网络对强回波的预报能力高于光流法。

4.3 个例分析 4.3.1 层状云降水过程1

个例为河南2017年5月3日发生的一次层状云降水过程，运用前2 h的雷达图像预报2 h降水情况。图4为预报时效为6、30、60、90、120 min的地面实况、光流法预测结果和预测神经网络的外推结果，从地面实况（图4a—e）可以看出，随着时间的推移，降水过程的30 dBz回波面积略有减小，光流法（图4f—j）在预报过程中，回波形态基本没有变化，而预测神经网络的回波（图4k—o）体现出了30 dBz回波略有减小的趋势。从两个模型预测图像与地面实况图像的结构相似性对比结果（图5）可以看出，在最初的约12 min内光流法的预报效果优于预测神经网络，但是其预测准确率下降较快，而预测神经网络在预测准确率上一直保持着相对较高的水平。

4.3.2 层状云降水过程2

个例为河南2018年6月18日发生的一次层状云降水过程，运用前2 h的雷达图像预报2 h降水情况。图6给出预报时效为6、30、60、90、120 min的地面实况、光流法预测结果和预测神经网络的外推结果，由于层状云降水回波变化较为缓慢，在预报30 min时间内，两个模型的预报准确率均较高。可以从地面实况（图6a—e）中看出左下角的小面积回波随着时间推移存在分裂现象，光流法的预报结果（图6f—j）并未体现这种变化的趋势，但是预测神经网络（图6k—o）成功预报出了这种回波的变化趋势。从两个模型预测图像与地面实况图像的结构相似性对比结果（图7）可以看出，在最初24 min预报时效内，光流法的预报结果要优于预测神经网络，但是在预报时效30 min后准确率下降幅度较大，而预测神经网络的预报准确率虽然随着预报时效的推移有所下降，但SSIM一直保持在0.8左右。从预报时效上看，预测神经网络较光流法更具有优越性。

4.3.3 对流云降水过程1

选取河南2020年6月24日发生的一次飑线天气过程，运用前2 h的雷达图像预报2 h降水情况。图8为预报时效为6、30、60、90、120 min的地面实况、光流法预测结果和预测神经网络的外推结果，通过对比预测图像与地面实况（图8a—e），预测神经网络（图8k—o）与实况相态上的相似度较高，光流法（图8f—j）在半小时后形变较为严重，且该方法预测的回波位置较实际偏右，而预测神经网络在预报位置上与地面实况对应较好。从两个模型预测图像与地面实况图像的结构相似性对比结果（图9）可以看出，两个模型在最初的18 min内都具有较高的预测准确率，但随着预测时间的推移，预测神经网络的SSIM始终高于光流法。

4.3.4 对流云降水过程2

选取2018年8月18日发生在福建的一次对流云降水过程，运用前2 h的雷达图像预报后2 h的降水情况。图10给出预报时效为6、30、60、90、120 min的地面实况、光流法预测结果和预测神经网络的外推结果，本次的对流云降水雷达回波强度最强超过50 dBz，且持续时间较长。从地面实况（图10a—e）可以看出，强回波为块状且面积一直保持一定的大小。光流法预测的图像（图10f—j）中，在预报时效60 min时，就出现了较大的形变，呈现出明显的分支状，但此时预测神经网络外推的雷达回波（图10k—o）中强回波为块状且保持一定的面积大小。但是预测神经网络在90、120 min预报时间段内强回波预测的位置较实际偏右，说明该预测神经网络虽然在强回波的预测上较光流法有所提升，但预报准确率的提升仍有一定的限制。从两个模型预测图像与地面实况图像的结构相似性对比结果（图11）可以看出，预测神经网络的预测图像与地面实况的相似度总体高于光流法。由于强回波变化较快，两个模型对其预测结果的准确率均低于层状云降水预测结果。

5 总结与讨论

文中介绍了神经网络在雷达回波外推上的应用，并使用带有Causal-LSTM单元组成的网络模型作为架构，在此基础上使用带权重的损失函数训练网络模型以消除降水的不平衡性和提高强回波预报的准确度，并与光流法的预报结果进行比较，得出以下结论：

（1）雷达回波外推的本质是时、空序列的预测问题，神经网络具有较强的非线性映射能力，可以充分利用大量的历史数据对雷达回波进行非线性外推。文中使用的预测网络模型具有提取空间特征的卷积模块和解决时间记忆问题的改进后的LSTM结构单元。结果表明预测神经网络对雷达回波外推的效果优于光流法。

（2）考虑到不同降水强度出现的频率不同，强降水对人类活动影响较大，网络在训练过程中赋予强回波较大的权重以提高对强回波的预报能力。测试集检验结果表明，预测神经网络2 h的回波外推在40和50 dBz检验阈值项上，临界成功指数、命中率较光流法均有所提高，虚警率有所下降，表明该网络对强回波的预报能力高于光流法。

（3）4个降水过程个例表明，光流法和预测神经网络对层状云降水过程的预报准确率均高于对流云降水，主要是由于层状云降水回波变化缓慢，对流云降水强回波变化较快。尽管两种临近预报方案的预报效果随预报时间的延长而下降，但是预测神经网络比光流法仍具有优越性。

（4）基于深度学习外推的4个降水个例的预测雷达回波在后期出现模糊现象，可能的原因是雷达回波的预测具有不确定性，随着时间的推移，对细节的预测难度更大，因此网络会更注重回波整体形态上的演变，而丢失了细节性的纹理。

现在使用的预测神经网络模型的输入仅为多个时刻雷达图像的反射率因子信息，可以考虑加入风场信息及数值预报环境场，使得网络可以学习到更多的物理变化机制，对天气系统的演变可能有更好的认识，进一步提高预报的准确度。同时，类LSTM单元或卷积模块也可以进一步优化，可以通过加深卷积等来获取更丰富抽象的变化特征。