北京地区汛期降水时空分布的统计降尺度研究*

郝 翠 于 波 戴 翼 智协飞 张迎新

1 北京市气象台,北京 100097 2 南京信息工程大学,气象灾害教育部重点实验室/气象灾害预报预警与评估协同创新中心,南京 210044

提 要：降水的时空降尺度方案一直以来是智能网格预报业务的难点,通过对比多种时间和空间的降尺度方法,凝练出适用于北京地区的最优方案。空间降尺度方面,对比反距离权重法、普通克里金法、最近邻法、双线性插值法、三维普通克里金法等5种方法的空间分布特征表明,双线性插值法在北京地区的应用效果最好,误差最小且ETS评分最高;时间降尺度方面,采用基于区域数值模式(睿图、CMA-MESO)产品的逐时分配和平均分配两种方案,其中睿图逐时分配、CMA-MESO逐时分配和平均分配在RMSE和MAE的误差表现差距不大,但睿图逐时分配在ETS上的效果最显著,且在强降水时段表现也较优,说明从预报准确率角度采用睿图逐时分配的优势更为明显。双线性插值和睿图逐时分配作为北京地区客观降水预报方法的时空降尺度方案,能够支撑智能网格业务提供精细化的预报产品,其成果可为相关业务研究提供借鉴。

引 言

数值预报后处理产品已成为数值产品释用的重要参考依据(赵声蓉等,2012;唐健等,2018;郝翠等,2019;张宏芳等,2019;金巍等,2020),但目前数值预报后处理产品的时空分辨率仍不能满足当前精细化预报服务的需求,因此发展适用于本地区的时空降尺度技术成为预报业务的关键问题。

降水具有时空不连续性、非正态分布等特点(王海霞和智协飞,2015;董旭光等,2018),因此时空降尺度一直以来是降水精细化预报业务的难点,不同降尺度方法所得到的插值结果不尽相同。在进行统计降尺度之前,必须先将基于站点或粗分辨率的网格产品插值到更细分辨率网格上(王亚男和智协飞,2012;智协飞等,2016;刘绿柳等,2018;吴遥等,2020)。对于空间降尺度,由于降水的局地性很大,尤其是在地形复杂区域,地理位置虽接近但降水量可能会差别很大,不同的空间插值方法产生的降水分布差异很大(邬伦等,2010),因此需要遴选出本地适用性最好的插值方法;目前广泛应用且效果较好的空间插值方法包括最近邻(Nearest)插值、双线性(Bilinear)插值(陈胜东等,2015)、反距离权重(inverse distance weighting,IDW)插值(包慧濛等,2015)、普通克里金(ordinary Kriging,OK)插值(陈胜东等,2015)等;而基于地形信息的三维插值方法(吴链,2017)进行空间降尺度,应用效果也较好。根据中国气象局的业务要求,本文将基于观测的站点资料插值成空间分辨率为5 km的格点产品,采用IDW、OK、Nearest、Bilinear、三维普通克里金(three-dimensional ordinary Kriging,TOK)5种插值方法对北京地区的降水进行空间降尺度试验,该结果可以为粗分辨率数值模式及数值预报后处理产品提供参考。

时间降尺度的研究相对较少。刘凑华等(2021)通过模拟雨带位置和形态来实现时间降尺度,其利用数值模式逐3 h预报进行位置订正、动态重构和频率匹配,将逐日降水客观预报降尺度为逐1 h。由于目前区域数值模式已有逐1 h预报产品,本文将观测的日降水资料基于区域数值模式逐1 h的预报产品进行时间降尺度,并与平均分配的结果进行对比分析,相关结论可为粗时间分辨率的数值模式释用及数值模式预报后处理产品提供参考。

本文针对降水客观预报产品进行业务化应用时遇到的时空降尺度问题,若直接利用降水客观预报产品或数值模式的日降水预报进行时空插值,则难以对比时空插值的真实效果。为最大程度反映时空插值的效果,利用观测的日降水资料和中尺度区域数值模式预报产品,进行时空降尺度试验,探讨北京地区最优的降水时空降尺度技术方法,为智能网格预报业务建设提供参考和依据。

1 资料和方法

考虑到观测数据的完整性和准确性,观测资料采用2019年4—9月北京地区共364个自动气象观测站(以下简称自动站,2015年自动站总数)逐时降水质量控制产品,并将逐时产品处理成逐日(前日21时至当日20时)降水产品。数值预报产品选取对北京预报性能较好的北京市气象局睿图模式(RMAPS)和中国气象局CMA-MESO的逐小时降水预报产品,预报时长为24 h。其中,睿图模式为北京市城市气象研究院于2002年开始研发的快速更新循环数值预报系统,至今已升级更新了3个版本,对2007年以来北京地区重大天气过程无一漏报(陈敏等,2014;仲跻芹等,2010;2017),截至目前可提供空间分辨率为3 km的地面要素及降水的集成预报(刘梦娟和陈敏,2014;闵晶晶,2014);CMA-MESO为中国气象局自主研发的区域数值天气预报系统,2006年7月投入业务应用,2014年起提供空间分辨率为10 km的地面降水预报产品,在天气预报业务中发挥了重要作用(沈学顺等,2020) 。本文所用时间,除特殊说明外,均为北京时。

1.1 空间降尺度方案

利用15 m空间分辨率的数字高程资料实现三维空间的降尺度处理。对比验证的观测资料为北京地区512个自动站(2019年自动站总数)降水实况(图1)。近年来北京地区的观测站点逐年增加,而数值模式后处理技术需要相对稳定的多年观测资料,因此进行时空降尺度未采用较新的512个自动站的实况信息,而是仍采用364个站(2015年自动站总数)的实况信息。

注:蓝点:2015年,364个;圆点:2019年,512个。图1 北京地区自动站分布图Fig.1 Distribution of automatic weather stations (dots) across the Beijing Region

常用的5种空间降尺度方案具体为:(1)IDW:以插值点与样本点间距离的倒数为权重进行加权平均(Bartier and Keller,1996)。(2)OK:以空间统计学作为理论基础,对空间分布数据进行全局变异函数计算,获取最优、线性、无偏的内插估计;其本质是由可用样本数据的线性组合来获得待插值点的降水量,其中加权值不仅与观测点和待插值点的距离有关,同样与观测点的位置、空间结构相关(Pardo-Igúzquiza and Dowd,1998)。(3)Nearest:以距离网格点最近的站点值为此网格点值来进行插值(陈欢欢等,2007)。(4)Bilinear:采用2个变量的插值函数的线性插值扩展,利用网格点附近的4个站点值的水平、垂直方向上的线性内插获取网格点的值(叶金印等,2013)。(5)TOK:即在二维插值的基础上加入高程信息,将海拔高度作为经度、纬度之外的第三维度进行插值(曾怀恩和黄声享,2007)。

利用上述方法将实况自动站观测资料插值至5 km×5 km网格上。针对空间降尺度,本文采用交叉验证法(Seaman,1983) 评价各方法的效果,即将逐日364个站的降水资料按8∶2的比例分成2组,用291个站资料进行插值,73个站的资料进行检验,为保证结果稳定,每日资料随机分组10次(图2)。

图2 北京地区时空降尺度方案技术路线图Fig.2 The technology roadmap of spatio-temporal downscaling methods for the Beijing Region

1.2 时间降尺度方案

采用的数值模式资料为2019年5月1日至9月30日的北京睿图模式和CMA-MESO的逐时降水预报,以及对应的364个站降水实况(逐日和逐时)。本文设计2种方案(图2):(1)当数值模式(睿图逐时预报、 CMA-MESO逐时预报)预报有降水时,按照当日数值模式逐1 h的降水比例分配;当数值模式预报无降水时,平均分配;当数值模式预报有降水、观测无降水时,按观测无降水分配。(2)对逐日观测降水进行逐时平均分配(平均分配)。

检验方法为计算各站点在整个降水过程中各时次观测值与估计值间的均方根误差(root mean square error,RMSE)、平均绝对误差(mean absolute error,MAE)、 公平风险(equitable threat score,ETS)评分。其中,MAE反映估计值的误差范围,取值越接近0表明误差越小;RMSE反映估计值的灵敏度和极值情况,取值越接近0表明误差越小(高华喜和殷坤龙,2007;宋丽琼等,2008);ETS反映预报准确率,是对TS评分的改进,对空报或漏报进行惩罚,使评分更公平,取值范围为-1～1,越接近1表明预报效果越好。具体计算为:

(1)

(2)

(3)

2 结果分析

2.1 空间降尺度

对比5种降尺度方案(IDW、OK、Nearest、Bilinear、TOK)的空间分布特征,通过计算各站点在整个降水过程中逐时次观测值与预报值间的误差,检验了不同插值方法误差的空间分布特征(表1)。结果表明,Bilinear的效果最好,RMSE和MAE最小,ETS评分最高;对10 mm以上降水的插值效果评估也表明Bilinear插值效果最好,ETS(>10 mm)可达0.65;OK的效果次之,Nearest、IDW的结果较差,这与高歌等(2007)对逐日降水的空间插值的结论一致,OK的结果相对于IDW较好。加入高程信息后的TOK效果最差,误差反而增大,因此将地形高度与降水量做相关分析,分析表明二者相关系数仅为-0.01,几乎无相关性。北京地区的复杂地形对降水的影响主要体现在对降水的触发上,即降水发生频次,而对于降水量来说,二者并无相关性,这与袁宇锋(2017)的研究结论一致。因此,对北京地区而言,在进行空间插值时,不考虑次级信息(高程)的二维效果反而比三维效果要好。

表1 不同空间降尺度方案的检验结果Table 1 The test results of different spatial downscaling schemes

由于北京地区地形复杂,西北为山地,东南为平原,按海拔高度将北京地区分为山区(海拔高度200 m及以上)和平原(海拔高度200 m以下)两部分,分别占总面积的59.29%和40.71%。对山区和平原5种降尺度方案的空间分布特征分析(表2)表明:(1)对山区而言,Bilinear效果最好,RMSE和MAE最小,分别为5.30 mm和2.81 mm,其次是OK和Nearest;(2)平原地区,IDW效果为最好,其RMSE和MAE分别为5.57 mm和3.25 mm,其次是Bilinear,RMSE和MAE分别为5.76 mm和3.28 mm。综合山区和平原的效果分析,Bilinear、Nearest、IDW、OK差距不大,但TOK与其他4种方法有明显差距。

表2 山区(海拔200 m及以上)和平原(海拔200 m以下)不同空间降尺度方案的检验结果Table 2 The test results of different spatial downscaling schemes for mountains (above sea level≥200 m) and plains (above sea level<200 m)

对北京地区3个主要的降水分布型(全市降雨、北多南少、东多西少)的检验表明(图3), OK虽然在降水量大的区域插值效果较好,但会增加降水分布范围,Nearest和IDW在刻画局部点的特征时容易失真造成误差较大,TOK效果最差,Bilinear在各种降水分布类型中都与观测实况最为接近(图3～图5),这与税军峰等(2019)的研究结论一致。从ETS评分的检验效果看,Bilinear的评分也是最高的。张永顺(2020)在对比海温插值方法时,同样认为Bilinear优于Nearest和三次插值。

图3 2019年6月3日北京地区小雨型降水不同空间降尺度方案的分布(a)实况,(b)Nearest,(c)IDW,(d)Bilinear,(e)OK,(f)TOKFig.3 Distribution of different spatial downscaling schemes for light rain across Beijing Region on 3 June 2019 (a) observation, (b) Nearest, (c) IDW, (d) Bilinear, (e) OK, (f) TOK

从降水类型看,2019年6月3日受高空槽和低层切变线影响,此次过程属于对流性降水,Bilinear方法的分布与实况最为接近(图3)。同样作为对流性降水,但对流性更强的2019年8月2日受高空槽和低层偏东风造成的降水,IDW、Bilinear、Nearest的分布与观测实况都很接近,其中IDW在降水量较大的分布区域与观测实况最为接近(图4)。2019年7月23日受副热带高压外围偏南暖湿气流和高空槽系统影响,以及大尺度稳定性降水和对流性降水的共同作用,此次过程属于混合性降水,IDW、Bilinear、Nearest的分布与观测实况都很接近,其中IDW在降水量较大的分布区域与观测实况最为接近(图5)。由于华北地区汛期降水以混合性和对流性为主,大尺度稳定性降水较少,且常发生在春秋两季,为了分析对大尺度稳定性降水的空间降尺度效果,本文回算了2018年4月21日的降水个例,此次降水过程整体雨势平稳,属于典型的稳定性降水,Nearest、Bilinear、IDW的分布与观测实况都很接近,其中Nearest和Bilinear在降水量较大的分布区域与观测实况最为接近,Bilinear整体分布较Nearest 更平滑(图6)。总体来看,在不同类型的降水个例中,Bilinear不一定是表现最好的插值方法,但它却是一个表现效果最为稳定的方法。

图4 2019年8月2日北京地区东多西少型降水不同空间降尺度方案的分布(a)实况,(b)Nearest,(c)IDW,(d)Bilinear,(e)OK,(f)TOKFig.4 Distribution of different spatial downscaling schemes for precipitation more in the east and less in the west of Beijing Region on 2 August 2019 (a) observation, (b) Nearest, (c) IDW, (d) Bilinear, (e) OK, (f) TOK

图5 2019年7月23日北京地区北多南少型降水不同空间降尺度方案的分布(a)实况,(b)Nearest,(c)IDW,(d)Bilinear,(e)OK,(f)TOKFig.5 Distribution of different spatial downscaling schemes for precipitation more in the north and less in the south of Beijing Region on 23 July 2019 (a) observation, (b) Nearest, (c) IDW, (d) Bilinear, (e) OK, (f) TOK

图6 2018年4月21日北京地区大尺度稳定性降水不同空间降尺度方案的分布(a)实况,(b)Nearest,(c)IDW,(d)Bilinear,(e)OK,(f)TOKFig.6 Distribution of different spatial downscaling schemes for large-scale stable precipitation in Beijing Region on 21 April 2018 (a) observation, (b) Nearest, (c) IDW, (d) Bilinear, (e) OK, (f) TOK

2.2 时间降尺度

对2019年4—9月逐时的统计表明(图7),基于睿图的逐时分配、基于CMA-MESO的逐时分配和平均分配在RMSE上的表现差距不大,表明三者在与观测值偏差上的表现基本一致;MAE的结果表明,睿图逐时分配和CMA-MESO逐时分配在大部分时段比平均分配的误差小,睿图逐时分配在误差上的表现相对于CMA-MESO更有优势;而从ETS评分看,不论是0.1 mm还是10 mm以上的降水,基于睿图模式的逐时分配都是3种方案中最高的,平均分别为0.22和0.05,说明从预报准确率看睿图逐时分配的优势更明显;CMA-MESO在0.1 mm以上降水的ETS评分最差,为0.08,但在10 mm以上ETS评分中较平均分配略高,为0.003。CMA-MESO逐时分配的效果与睿图逐时分配的效果差距较大,原因可能与2个数值模式的空间分辨率差异较大有关,CMA-MESO的空间分辨率为10 km,插值到自动站可能也会产生较大误差,未来利用自2020年起业务化运行的CMA-MESO 3 km预报产品时(黄丽萍等,2022),可能会有较好的结果,并且由于睿图模式融合了地基GPS和京津冀29部雷达的实时资料,对降水临近预报更接近实况,因此睿图逐时分配的方案ETS评分较高。这与刘凑华等(2021)的结果一致,其通过位置订正和动态重构方法对逐日降水预报做逐时分配,ETS评分也有显著提高。

图7 2019年4—9月不同时间分配方案的逐时降水检验评分(a)RMSE,(b)MAE,(c)ETS(>0.1 mm),(d)ETS(>10 mm)Fig.7 Verification scores of different temporal allocation schemes for hourly precipitation from RMAPS, CMA-MESO and observation between April and September 2019 (a) RMSE, (b) MAE, (c) ETS (>0.1 mm), (d) ETS (>10 mm)

以2019年8月6—7日强降水时段(6日23时至7日02时)为例(图8),基于睿图逐时分配的空间效果与实况空间分布一致性较高,CMA-MESO没有预报出该时段内的强降水分布特征,平均分配在强降水上没有任何体现。

图8 2019年8月(a～d)6日23时,(e～h)7日00时,(i～l)7日01时,(m～p)7日02时不同时间分配方案北京地区强降水时段的分布(a,e,i,m)实况,(b,f,j,n)睿图,(c,g,k,o)CMA-MESO,(d,h,l,p)平均分配Fig.8 Spatial distribution of different temporal allocation schemes for heavy precipitation cases in Beijing Region at (a-d) 23:00 BT 6, (e-h) 00:00 BT 7, (i-l) 01:00 BT 7 and (m-p) 02:00 BT 7 August 2019(a, e, i, m) observation, (b, f, j, n) RMAPS, (c, g, k, o) CMA-MESO, (d, h, l, p) average

3 时空降尺度方案在降水客观方法预报中的应用

目前,北京市气象台降水客观技术方法预报基于睿图数值模式预报产品,采用频率匹配和偏差订正方法,获取对模式的逐日预报订正结果(以下简称降水客观方法预报),最大程度地减小预报误差,自动得到各时段最优的降水预报结果。降水客观方法预报采用的是上文的时空降尺度方案,与睿图区域模式预报相比,2019年6—9月北京地区0～3 d的平均晴雨预报准确率提高了4.1%;0～3 d的小雨TS评分提高了2.0%(图略)。

2019年8月6—7日降水个例(图9)分析表明,逐日的预报相对观测实况较大,降水客观方法预报在睿图模式的基础上订正后虽然更接近实况的降水分布,但降水量仍然偏大。降水客观方法预报的逐日预报根据睿图逐时分配的结果(图9f,9i,9l,9o)也较为接近实况,在强降水中心的预报上表现较优。因此,降水客观方法预报采用Bilinear空间插值和睿图逐时分配应用效果良好,可以为业务提供精细化的预报结果,可满足业务需求。

图9 2019年8月(a～c)6日21时至7日20时,(d～f)6日23时,(g～i)7日00时,(j～l)7日01时,(m～o)7日02时不同时间分配方案北京地区强降水时段的分布(a,d,g,j,m)实况,(b,e,h,k,n)睿图预报,(c,f,i,l,o)降水客观方法预报Fig.9 Spatial distribution of different temporal allocation schemes for heavy precipitation cases in Beijing Region (a-c) from 21:00 BT 6 to 20:00 BT 7, (d-f) at 23:00 BT 6, (g-j) at 00:00 BT 7, (j-l) at 01:00 BT 7 and (m-o) at 02:00 BT 7 August 2019 (a, d, g, j, m) observation, (b, e, h, k, n) RMAPS forecasts, (c, f, i, l, o) forecasts by objective prediction method

4 结论与讨论

本文讨论的时空降尺度方法,都是基于数值模式后处理的降水客观预报产品,数值预报模式和降水客观预报产品的时空分辨率精度不一,在业务应用中需要考虑时空降尺度问题,主要是为了加强数值模式及数值预报后处理产品的精细化业务应用,在不影响原效果的情况下满足预报产品精细化和智能化的业务需求。

对空间降尺度的研究表明,Bilinear插值效果最好,其RMSE和MAE为最小,ETS评分最高;加入高程信息的TOK误差较大,可以认为北京地理高程对降水的影响主要与降水发生频次有关,与降水量关系不大。对不同的降水分布型的空间分布结果表明,Bilinear的效果都是与实况最接近的,OK在降水量大的区域插值效果较好,但会增加降水分布范围,Nearest和IDW在刻画局部点的特征时容易失真造成误差较大,TOK效果最差。对山区和平原分区的分析统计表明,山区Bilinear效果最好,而平原IDW效果最好,从山区和平原综合的效果分析,Bilinear、Nearest、IDW、OK差距不大,但TOK与其他四种方法差距较大。对不同降水类型的个例分析表明,Bilinear、Nearest、IDW 效果最好且差距不大,其中Bilinear不一定是表现最好的插值方法,但它是一个表现最为稳定的方法。

对时间降尺度的研究表明,睿图逐时分配、CMA-MESO逐时分配和平均分配在RMSE上的表现差距不大,表明二者在与观测值的偏差上的表现基本一致;而从ETS评分看,不论是0.1 mm还是10 mm以上的降水,睿图逐时分配比其他方法要显著,说明从预报准确率角度睿图逐时分配的优势更明显。

本研究的时空降尺度方案均未考虑降水的发生频率、空间分布特点、时间移动规律,但本研究的结论可为未来继续深入进行降水预报的时空降尺度研究打下较好基础。

上述时空降尺度方案已应用于北京市气象台的降水客观方法预报,并具有较好的业务应用效果。此外,相关时空降尺度方案还可用于北京地区的智能网格业务系统中,提高灾害天气的网格预报准确率,为业务人员提供更为丰富、更加精细的参考产品。