基于集成学习的巢湖面雨量计算方法

王 杰，徐 祥，罗晓丹，张 萌，黄 澈，洪冠中，汪 翔

（1.安徽省农村综合经济信息中心，安徽 合肥 230031；2.广东省清远市气象局，广东 清远 511515；3.安徽省公共气象服务中心，安徽 合肥 230031）

0 引 言

巢湖是中国五大淡水湖之一，地处中纬度地带，由于其位于南北气候过渡带，冷暖气团常在此交汇而造成暴雨等恶劣天气频发［1］。建国以来，巢湖发生较大洪涝灾害的年份共计18 年［2］，尤其2020 年最为严重，与巢湖相关的主要河湖洪水位均突破历史极值，严重威胁人民群众的生命财产安全。巢湖各条支流、水库、堤防的抗洪能力都与流域和湖面面雨量息息相关，因此开展面雨量计算和监测可为应急管理部门提供数据信息支撑，为防汛抗洪的科学决策和合理部署提供重要依据。

目前，面雨量的计算方法主要有等值线法、数值法、算术平均法和网格点平均法等［3］。其中，等值线法计算较为复杂，难以满足实际应用要求。数值法是指泰森多边形法或三角形法，该方法通过确定区域内各气象观测站的权重系数计算面雨量，适用于气象观测站或降雨量分布不均匀的区域。算术平均法是对区域内各站测得的同期雨量求算术平均从而计算面雨量，仅适用于区域面积小、地形起伏不大且观测站多而分布又较均匀的场景［4］。网格点平均法［5］是利用GIS 技术将区域按设定的经纬度分辨率进行网格划分，通过计算网格点上的降雨量实现面雨量计算，其中网格点降雨量多基于天气雷达、气象卫星和地面观测等资料的某种或多种的融合进行估算，从时空分辨率精度及空间覆盖范围方面考虑，采用雷达数据估测网格点降水具有明显优势。

传统的雷达估测降水常采用Z-R 关系法构建雷达降水估测模型［6］，通过雷达反射率因子Z 直接反演降水R。Z-R 关系随地区、季节、降水类型而变化，甚至同一次降水过程中亦可能发生变化，这种变化直接影响雷达降水估测精度［7］。前人研究多基于计算复杂、专业性较强的数值模式拟合Z-R 关系，而以集成学习为代表的人工智能技术在雷达估测降水中的应用鲜有报道。近年来，随着人工智能技术的发展，集成学习因其优越的泛化性能，被学者们应用于气象方面研究。张钧民等［8］利用XGBoost 算法开展多源降水数据融合方法研究；于霞等［9］以雷雨大风的46 个相关特征属性数据为输入，降雨和大风的数值为输出，利用XGBoost 算法建立雷雨大风预测模型；陈志月等［10］选用GPR、CatBoost、XGBoost 这3 种集成学习算法预测江西地区水面蒸发量；陈双等［11］基于云顶温度、中层融化参数等构建的决策树判别模型，可较好地提升临界气温下雨、雪的判别准确率。集成学习虽在气象领域不同场景开展应用研究，但此类算法在雷达降水估测的适用性还需进一步探究。

本文以巢湖为研究对象，提出一种基于集成学习算法和网格点平均法相结合的面雨量计算方法，通过对巢湖流域边界内的雷达三维拼图数据和气象观测站降水数据进行预处理，构建以雷达多层拼图数据为特征、气象观测站小时降水量为标签的数据集，利用所构建数据集对模型进行训练迭代，结合GIS 技术，将巢湖区域按经纬度网格划分并与雷达拼图在空间上进行叠加，根据QPE 模型计算各网格点降雨量，从而实现巢湖面雨量的自动估算，以期为提升巢湖防灾减灾能力提供重要信息支撑。

1 资料与方法

1.1 研究区域概况

巢湖流域地处安徽省中部，位于东经116°24'～118°0'，北纬30°58'～32°06'，面积为13486 km2［12-13］。巢湖作为巢湖流域的中心，从南、西、北三面汇入，西面、北面均为山区，是巢湖重要支流的发源地，东面经巢湖闸由裕溪河及牛屯河注入长江（见图1）。巢湖及流域地势轮廓都是东西长南北窄，流域内地形复杂多样，有山地、平原、丘陵等多种地形类型，且西高东低，中间地势低洼平坦，区域地形地貌极易形成山洪冲击、江洪顶托等不利因素，从而使得巢湖防汛抗洪工作面临严峻挑战［13］。

图1 巢湖流域及气象观测站点分布图

1.2 数据来源

本 文以2019 年6 月5 日01 时—6 日01 时 的一次降水过程为例，开展巢湖面雨量计算方法研究。数据包括雷达三维拼图数据、气象观测站降水数据、巢湖流域DEM 高程数据和基于气象陆面数据同化系统（CMA Land Data Assimilation System，CLDAS）的降水数据［14］。其中，雷达三维拼图数据来自灾害天气短时临近预报系统（Severe Weather Automatic Nowcast System，SWAN）［15］，主要包含经纬度、时间、层数、站点信息等29 个字段信息，拼图按高度0.5 km～19 km 排列共21 层，每层包含800×900 个网格，空间分辨率0.01°×0.01°，时间分辨率为6 min。图2为3 km 高度巢湖流域的雷达拼图，6 月5 日14 时雷达回波强度普遍小于21时，这表明6月5日21时降水概率和强度普遍大于14时。降水数据和CLDAS数据均来自气象大数据云平台“天擎”，降水数据包括巢湖流域范围内375 个气象观测站，其中国家站10 个，区域站365 个（见图1）；CLDAS 数据为空间分辨率为0.0625°×0.0625°，时间分辨率为1 h 的降水资料。本文雷达三维拼图数据和气象观测站降水数据作为构建雷达降水估测模型的数据集，基于CLDAS 得到的面雨量数据用于检验本文模型的精度。

图2 3 km雷达拼图回波变化

1.3 面雨量计算方法

1.3.1 基于集成学习算法构建雷达QPE模型

鉴于集成学习算法能较好地反应地面降水与雷达拼图各层降水数据之间的非线性关系，为提高雷达估测降水的精度，本文选用3 种经典集成学习算法RF、XGBoost、LightGBM 开展雷达估测降水方法研究，以期获得性能较好的雷达QPE模型。

RF 属于Bagging 类算法，通过自助法（Bootstrap）对训练集样本和特征有放回地随机采样，构建许多相互不关联的回归决策树［16］。RF 并行训练所有决策树，每棵决策树得到一个预测结果，将所有树的预测结果结合取平均值，从而得到整个模型的预测结果。XGBoost、LightGBM 属于Boosting 类算法，两者是梯度提 升 决 策 树（Gradient Boosting Decision Tree，GBDT）的进化版本［17］，其中，XGBoost 是在决策树构建过程中加入了正则项，还使用一阶、二阶梯度函数来优化损失函数，这样做是为了避免过拟合和降低模型复杂度［17-18］；而LightGBM 是通过采用直方图算法（Histogram）、基于梯度的单边采样算法（GOSS）、互斥特征捆绑算法（EFB）分别对分裂点数、样本量和特征数量3个变量进行优化［19］，并在叶子生长策略上采用带深度限制的按叶子生长（leaf-wise）算法和对缓存命中率进行优化，在提高预测效率和准确率的同时防止模型过拟合。引入3 种集成学习算法建立的雷达QPE模型如式（1）：

其中，j表示RF、XGBoost、LightGBM表示采用集成学习算法获得的降水预测数据；fEM表示基于集成学习算法构建的非线性回归预测模型；Radarlayeri表示雷达拼图的第i层降水数据。本文选择雷达拼图的5 个低仰角数据作为训练集的特征输入，这样做一方面是由于雷达降水信息主要反映在雷达体扫5 个低仰角上，另一方面可以减少模型训练数据量。

1.3.2 基于网格点降雨量计算面雨量

巢湖面雨量是整个湖面区域单位面积上的平均降水量［4］。基于网格点降雨量的计算方法可以利用区域内一切可获得的降水信息，能较精确地计算面雨量，故本文采用基于GIS 技术的网格点平均法计算巢湖面雨量。其具体实现步骤如下：

1）通过GIS 技术，将巢湖区域按0.01°×0.01°经纬度网格划分并与雷达拼图在空间上进行叠加，建立雷达拼图经纬度网格与巢湖区域的映射关系，并对湖面边界内网格点进行标记（见图3）。

图3 雷达拼图与巢湖边界的映射关系

2）根据雷达QPE模型计算各网格点的降雨量。

3）基于步骤2 的计算结果，对湖面区域所有网格点降雨量进行算术平均，以实现湖面面雨量估算。具体计算公式如式（2）所示，其中n为网格数为雷达QPE模型预测结果。

1.4 模型超参数调优

不同超参数会对模型性能造成较大差异［20］，鉴于超参数调优需要花费大量的计算时间和资源，本文选择网格搜索法提高超参数搜索效率，通过在Python语言环境下加载sklearn（Version: 0.23.2）开源库的GridSearchCV函数实现［21-22］。具体步骤如下：

1）根据RF、XGBoost、LightGBM 这3 种算法模型特点和超参数的原理选择需要调整的超参数组合，并将每一个超参数值设置一个比较宽松的调参区间［23］。

2）采用网格搜索法进行训练迭代，每次训练结束后根据模型评价指标和基于CLDAS 得到的面雨量数据评估调参效果，逐渐缩小调参范围，直至得到最优调优结果。

3）根据步骤2 中调优结果，在超参数组合中选择某一超参数为变量，其他超参数值按步骤2 中结果保持不变。改变此超参数的值开展训练迭代验证对模型性能的影响，若影响程度较小则可将此超参数剔出最佳参数组合，后续训练中使用超参数默认值。反之，此超参数则是影响模型性能关键超参数。其他超参数按此方法依次进行验证。

1.5 模型评价指标

本文以决定系数（R2）、平均绝对误差（MAE）和均方根误差（RMSE）作为雷达QPE 模型的评价指标［24］，计算公式为：

2 实验与结果分析

2.1 实验环境

实验硬件环境为16 GB NVIDIA Tesla V100 GPU，操作系统为CentOS Linux7.5，选用的机器学习开源架构为Tensorflow（Version: 2.3.0），通过Tensorflow调用GPU训练［25］雷达QPE回归预测模型。

2.2 数据集制作

表1 给出了气象观测站降水数据预处理部分数据实例。

表1 气象观测站降水数据预处理部分数据实例

1）气象观测站降水数据预处理：流域内气象观测站降水原始数据有9377 条，为了使训练的结果更加真实有效，剔除降水缺测记录、异常记录、没有通过质量控制的数据和不完整观测记录。同时，由于0.1 mm 的降水数据为非有效数据，会影响训练模型的评价指标［26］，因此删除降水量为0.1 mm 的观测记录。经过上述操作后最终筛选出7376条有效数据。

2）雷达三维拼图数据预处理：提取流域内气象观测站空间经纬度对应的雷达21 个仰角反射率，将6 min/次的拼图在特定时次内求算术平均得到逐小时平均反射率。表2 中列举了处理后的逐小时平均反射率（layer0～layer4）部分数据实例。其中，dBZ 表示雷达反射率因子的单位。

表2 模型数据集部分实例

3）制作以小时降水量为标签的数据集：提取气象观测站各时次降水量，按站点和时间，标注雷达数据，形成以雷达5 个低仰角逐小时平均反射率（layer 0～layer4）为特征，以小时降水量（pre）为标签的数据集，建立气象观测站小时降水量与雷达逐小时平均反射率映射。

2.3 模型训练

数据集分为训练集和验证集，训练集和验证集按8:2 划分。模型训练流程如图4 所示，采用RF、XGBoost、LightGBM 这3 种集成学习算法模型分别开展训练迭代，并基于模型评价指标进行不同集成算法下模型性能对比。设定不同超参数开展训练对比实验，训练完成记录每次训练时长、GPU 内存消耗、评价指标值、输出2019 年6 月5 日01 时—6 月6 日00 时巢湖实况面雨量并保存雷达QPE模型。

图4 基于集成学习的巢湖面雨量训练流程图

2.4 超参数调优结果

将每次训练结束后超参数的调优结果和效果进行对比，得到RF、XGBoost、LightGBM 模型最佳超参数组合和调优结果。表3以数据形式直观给出了3种算法模型最佳超参数组合、获得调优结果前确定的调参范围和得到调优结果的训练时间。其中，n_estimators 表示回归树数量；max_depth 表示树的深度；min_samples_split 表示拆分内部节点所需最少样本数；min_samples_leaf 表示在叶节点处需要最小样本数；min_child_weight 表示子节点所需的实例权重的最小总和；learning_rate 表示学习率；num_leaves 表示每棵回归树上的叶子数。调参前超参数数值代表含义，以n_estimators 为例，（100，500，10）表示调参范围为100～500，步长为10。

表3 模型超参数组合和训练时间

本文以RF 模型的n_estimators 为例，验证表3 调优结果的可信度和超参数是否为模型关键超参数。图5 为n_estimators 值改变时，MAE、RMSE 和R2这3个评价指标的变化情况。由图4可见，当n_estimators变化时，对模型评价指标影响明显，结合n_estimators=480 时RF 模型面雨量预测结果（见表4），可以判定n_estimators 是RF 模型的关键超参数，表3 中RF 模型的n_estimators调优结果可信。

表4 集成学习算法+网格点平均法与CLDAS面雨量计算结果比较 单位：mm

图5 不同超参数值对模型性能的影响

2.5 集成学习算法模型评估

经过多次训练迭代，采用RF、XGBoost 和Light-GBM 这3 种集成学习算法模型预测降水的评价指标结果如表5所示。

表5 不同算法模型评价指标结果

由表5 可知，3 种算法模型训练时长和内存开销差别很小。在训练期最优模型是RF 预测模型（MAE=0.492 mm/h，RMSE=1.254 mm/h，R2=0.933）；在验证期XGBoost 的R2值（R2=0.697）虽然略高，但MAE 和RMSE 表现不如其他2 种模型，综合分析各项评价指标，RF 预测模型（MAE=1.342 mm/h，RMSE=3.435 mm/h，R2=0.689）优 于XGBoost 和LightGBM 模 型。综上，基于RF 算法的雷达QPE 估测模型在本次降水过程中表现最佳。

2.6 实验结果验证

输出RF、XGBoost、LightGBM 这3 种模型2019 年6 月5 日01—6 月6 日00 时巢湖面雨量的计算结果，并与基于CLDAS 得到的面雨量结果对比，以此验证基于集成学习算法和网格点平均法相结合的方法在实际场景中应用效果。

由表4 可知，在2019 年6 月5 日01 时—12 时，基于RF 的雷达QPE 模型和网格点平均法相结合的面雨量计算结果与CLDAS 中面雨量计算结果均为0 mm；在13 时—21 时和23 时，两者结果相差在1 mm以内，但在22 时和2019 年6 月6 日00 时两者计算结果相差均大于1 mm，这点在图6中体现得更为直观。图6是4种模型面雨量计算结果趋势图，由图可见，基于RF 的雷达QPE 模型和网格点平均法相结合的面雨量计算方法与基于CLDAS 的计算结果更加接近，但在后面3 个时刻虽然4 种计算方法结果趋势一致，但离散程度增大，XGBoost偏离程度最大。

图6 4种巢湖面雨量计算方法结果趋势图

3 结束语

本文提出了一种基于集成学习算法与网格点平均法相结合的巢湖面雨量计算方法，通过对巢湖流域内雷达三维拼图数据和气象观测站降水数据进行预处理，建立了模型训练数据集，构建了3 种基于集成学习算法的雷达QPE 模型，使用GridSearchCV 函数对模型超参数进行了调优，评价了3 种集成学习算法在雷达降水预测中的表现，最后引入基于CLDAS 的面雨量计算方法对本文面雨量计算方法进行了验证。结果表明：

1）RF、XGBoost 和LightGBM 这3 种网络模型中，RF模型在验证期表现最佳，评价指标MAE、RMSE、R2分别为1.342 mm/h、3.435 mm/h、0.689，为适用于本次降水过程的预测模型。

2）对比分析了本文与CLDAS 这2 种方法得到的湖面面雨量计算结果，计算结果数值虽有差异，但总体趋势一致，说明使用本文提出的方法计算湖面面雨量可行。

本文为进一步开展巢湖及其流域防汛抗洪方法研究、建立科学合理的湖面面雨量阈值、增强防汛抗洪监测预警能力［27］提供了新的方法，当前可作为数值模式面雨量计算的有力补充。但本文方法还存在一些不足，对比本文方法得到的面雨量结果与数值模式计算结果在后面3 个时刻发现，计算结果与实际值偏离程度增大，说明有些因素没有考虑到。下一步将研究地形、降水相态、气象观测站与雷达拼图在空间不同对应关系等因素对模型的影响，以期进一步提高巢湖面雨量计算精度。