2 国家气象中心, 北京 100081
3 北京城市气象研究所, 北京 100081
2 National Meteorological Center, Beijing 100081
3 Institute of Urban Meteorology, China Meteorological Administration, Beijing 100081
大气是一个混沌系统,初始条件或预报模式的微小误差均会导致预报技巧丧失,这使得确定性数值预报不可避免地具有误差或者说不确定性。集合预报是解决"单一"确定性数值预报"不确定性"问题的新一代随机动力预报理论和方法(Leith, 1974),通过积分略有差异的多个初值或采用多个模式来估计预报误差概率分布,提供更多的确定性预报信息(Steven Tracton and Kalnay, 1993;Molteni et al., 1996;陈静等,2003;马旭林等, 2008, 2015;Saito et al., 2012)。集合预报的关键问题是如何合理描述数值预报初值误差源和分布特征,即采用何种初值扰动方法,能够正确描述初值误差三维结构,并正确给出初值误差(或者说大气可预报性)随时间的演变。理想的初值扰动方法不仅能使扰动场特征与实际分析资料中可能的误差分布一致,保证每个初值场都有可能代表大气的实际状态,而且能使每个扰动初值在模式中的演变方向尽可能发散,保持合理的初值扰动增长率,故合理的初值扰动结构及扰动增长特征对集合预报是至关重要的(Molteni et al., 1996)。
自1992年欧洲中期天气预报中心(ECMWF)和美国环境预报中心(NCEP)首次建立业务集合数值天气预报系统以来,目前常用的业务集合数值天气预报初值扰动方法主要划分为三类。第一类方法是产生沿模式相空间最不稳定方向发展的扰动场,如奇异向量法(Singular Vectors,SVs)(Buizza and Palmer, 1995; Buizza et al., 2005;Kim and Morgan, 2002)及增长模繁殖法(Breeding Growth Method, BGM)(Toth and Kalnay, 1993, 1997;肖玉华等,2011),此类方法的优点是可以获得增长最快的扰动结构,缺点是不能较好的代表资料同化分析过程中观测资料的不确定性。第二类方法是考虑同化分析中观测资料的不确定性,如CMC的观测资料扰动法(Houtekamer et al., 1996)、NCEP的集合卡尔曼滤波同化分析扰动法(Torn, 2010)、ECMWF的SVs奇异向量和集合资料同化(Ensemble of Data Assimilation,EDA)扰动组合法,即热带地区初值扰动场由EDA同化系统产生,中高纬地区初值扰动场则继续采用SVs方法获得(Buizza et al., 2008)。第三类方法既考虑扰动场沿着模式相空间最不稳定方向发展,同时又考虑观测资料分布的不确定性,如集合变换(Wei et al., 2008)和集合变换卡尔曼滤波(Ensemble Transform Kalman Filter,ETKF)初值扰动方法(Wang et al., 2004;黄红艳等,2016;庄潇然等,2017)。英国气象局和韩国气象局建立了基于ETKF初值扰动方法的全球集合预报业务系统(Bowler et al., 2009;Kay and Kim, 2014)。
ETKF是基于集合变换方法(Ensemble Transform,ET)和卡尔曼滤波理论(Kalman Filter,KF)发展而来的,是一种次优的卡尔曼滤波方案(Wang and Bishop, 2003),由(Bishop et al., 2001)首先针对适应性观测提出,并用于集合预报初值扰动的生成,通过集合变换和无量纲化求解与观测分布有关的误差协方差矩阵,可以快速估计出不同观测所造成的分析误差,反映观测资料对初值不确定性的影响,且能通过繁殖循环来使初值扰动识别快速增长的误差方向。Wang and Bishop(2003)将ETKF与增长模繁殖法进行对比,结果显示,ETKF初值扰动方法产生的分析扰动能够反映观测资料的空间分布,且扰动相互正交,有效地解决了BGM中分析误差方差固定不变等缺陷。为了使集合扰动相对于集合平均中心化,Wang et al.(2004)引入了两种集合中心化方案,结果表明球面单型中心化方案优于常用的集合扰动正负对称法。Wei et al.(2006)将ETKF初值扰动方法应用在NCEP(National Centers for Environmental Prediction)业务环境预报中,结果表明ETKF产生的扰动具有较大的自由度。Bowler et al.(2008)和Kay and Kim(2014)研究了ETKF初值扰动特征,结果显示,基于ETKF方法的初始扰动集合离散度分布与观测密度及大气动力属性有关。
自2008年起,中国气象局数值预报中心建立了基于ETKF初值扰动方法的GRAPES-REPS(Global and Regional Assimilation and Prediction Enhanced System-Regional Ensemble Prediction System)区域集合预报系统。王太微(2008)采用BGM和ETKF两种初值扰动方法进行了GRAPES区域集合预报系统对比试验,结果表明两种方法均具有一定的概率预报技巧,但ETKF方法对降水落区的预报效果略优于BGM方法。龙柯吉等(2011)对ETKF初值扰动方法的观测误差方差进行修订开展GRAPES区域集合预报试验,结果显示ETKF方案优于随机扰动方法。2014年8月,基于ETKF初值扰动方法的GRAPES区域集合预报系统实现业务化运行(张涵斌等,2014b)。2016年3月,在ETKF初值扰动循环基础上,初值扰动方法更新为多尺度混合方法(Multi-Scale Blending, MSB),即将T639全球集合预报大尺度扰动叠加在ETKF初值扰动循环产生的天气尺度及中小尺度扰动中,以增加集合离散度(Zhang et al., 2015),系统核心技术依然是ETKF初值扰动方法,故深入评估ETKF方法,是认识并提高GRAPES区域集合预报系统对预报不确定性描述能力的关键。但过去多从系统运行角度上评估ETKF方法,对ETKF初值扰动结构及扰动增长特征的分析诊断较少。为进一步认识ETKF初值扰动在中国区域的结构和演变特征,改进GRAPES区域集合预报效果,将利用15 km和10 km两种水平分辨率的GRAPES-REPS区域集合预报模式,对2015年6月1~15日进行集合预报初值扰动模拟试验,深入分析ETKF初值扰动结构及扰动增长特征,为发展更优的ETKF初值扰动方法提供依据。
2 GRAPES-REPS模式ETKF初值扰动方法简介 2.1 ETKF初值扰动数学公式ETKF初值扰动方法是基于卡尔曼滤波理论发展而来,由集合预报扰动方差和观测误差方差快速估计分析误差的一种初值扰动技术。假设集合预报成员数为K,第i个集合预报成员的预报扰动和分析扰动定义为
$ {\mathit{\boldsymbol{Z}}^{\rm{f}}} = \frac{1}{{\sqrt {K - 1} }}\left({\mathit{\boldsymbol{z}}_{\rm{1}}^{\rm{f}}, \mathit{\boldsymbol{z}}_2^{\rm{f}}, \cdots \mathit{\boldsymbol{z}}_K^{\rm{f}}} \right), {\mathit{\boldsymbol{P}}^{\rm{f}}} = {\mathit{\boldsymbol{Z}}^{\rm{f}}}{\left({{\mathit{\boldsymbol{Z}}^{\rm{f}}}} \right)^{\rm{T}}}, $ | (1) |
$ {\mathit{\boldsymbol{Z}}^{\rm{a}}} = \frac{1}{{\sqrt {K - 1} }}\left({\mathit{\boldsymbol{z}}_{\rm{1}}^{\rm{a}}, \mathit{\boldsymbol{z}}_2^{\rm{a}}, \cdots \mathit{\boldsymbol{z}}_K^{\rm{a}}} \right), {\mathit{\boldsymbol{P}}^{\rm{a}}} = {\mathit{\boldsymbol{Z}}^{\rm{a}}}{\left({{\mathit{\boldsymbol{Z}}^{\rm{a}}}} \right)^{\rm{T}}}. $ | (2) |
将预报扰动向量Zf更新为分析扰动Za是利用变换矩阵T来实现的,即通过T矩阵对当前时刻预报扰动向量进行线性组合产生当前时刻的分析扰动(张涵斌等,2014a):
$ {\mathit{\boldsymbol{Z}}^{\rm{a}}} = {\mathit{\boldsymbol{Z}}^{\rm{f}}}\mathit{\boldsymbol{T}}, $ | (3) |
变换矩阵T的计算公式为
$ \mathit{\boldsymbol{T}} = \mathit{\boldsymbol{C}}{\left({\mathit{\boldsymbol{ \boldsymbol{\varGamma} }} + \mathit{\boldsymbol{I}}} \right)^{ - 1/2}}{\mathit{\boldsymbol{C}}^{\rm{T}}}, $ | (4) |
其中,C是扰动向量Zf投影到观测空间中的协方差矩阵(以下简称投影矩阵E)的特征向量,Γ则是包含了所有特征值的(K-1)×(K-1)对角矩阵(Wang et al., 2004)。
$ \mathit{\boldsymbol{E}} = {\left({{\mathit{\boldsymbol{R}}^{ - 1/2}}\mathit{\boldsymbol{H}}{\mathit{\boldsymbol{Z}}^{\rm{f}}}} \right)^{\rm{T}}}\left({{\mathit{\boldsymbol{R}}^{ - 1/2}}\mathit{\boldsymbol{H}}{\mathit{\boldsymbol{Z}}^{\rm{f}}}} \right), $ | (5) |
其中,R为观测误差协方差矩阵,H为线性化观测算子,HZf为预报扰动向量Zf在观测空间H的投影矩阵。由式(5)可见,采用了观测误差协方差R的平方根对HZf进行归一化处理,因而归一化观测空间预报扰动向量是正交的,由此保证了分析扰动Za的正交性,(Γ+I)-1/2是特征向量矩阵C的尺度化因子。如果集合预报平均是真实大气的最佳估计,则集合预报初始扰动场之和应该为0,Wang et al.(2004)通过采用球面单型中心化方法,在转置矩阵T后面乘上特征向量矩阵C的转置矩阵,使得集合初始扰动场变换为无偏的,实现ETKF产生的初值扰动量之和等于0。
由于集合成员数K远远小于模式预报相空间的自由度,故公式(3)得到的分析扰动方差小于真实的分析误差方差,为解决此问题,在公式中引入放大因子(Wang and Bishop, 2003),故n时刻的分析扰动为
$ \mathit{\boldsymbol{Z}}_n^{\rm{a}} = \mathit{\boldsymbol{Z}}_n^{\rm{f}}{\mathit{\boldsymbol{T}}_n}{\mathit{\Pi }_n}, $ | (6) |
放大因子Πn的计算公式为
$ {\mathit{\Pi }_n} = {\mathit{\Pi }_{n - 1}}\sqrt {{\alpha _n}}, $ | (7) |
其中,参数α具体公式为
$ {\alpha _n} = \frac{{\mathit{\boldsymbol{\tilde d}}_n^{\rm{T}}{{\mathit{\boldsymbol{\tilde d}}}_n} - N}}{{{\rm{trace}}\left({\mathit{\boldsymbol{\tilde HP}}_n^{\rm{f}}{{\mathit{\boldsymbol{\tilde H}}}^{\rm{T}}}} \right)}} = \frac{{\mathit{\boldsymbol{\tilde d}}_n^{\rm{T}}{{\mathit{\boldsymbol{\tilde d}}}_n} - N}}{{\sum\limits_{i = 1}^{K - 1} {{\lambda _i}} }}, $ | (8) |
式中,
控制预报模式版本采用GRAPES-MESO V4.0,垂直坐标采用地形追随高度坐标,垂直层次为50层,集合预报成员数为14个集合扰动成员,加上控制预报共15个集合成员,集合预报成员物理过程略有差异(张涵斌等,2014b),控制预报初值采用GRAPES-3DVAR系统生成的分析场,背景场和侧边界条件由数值预报中心T639全球集合预报系统提供,初值扰动采用ETKF方法。
GRAPES-REPS区域集合预报的ETKF初值扰动方法采用12 h集合扰动循环计算方案(即将集合预报系统12 h短期预报场构成的预报扰动更新为当前时刻的分析扰动),计算投影矩阵E采用的变量为GRAPES-MESO模式面纬向风U和经向风V,采用的垂直层次为11层。观测资料采用模拟观测资料,即将T639全球模式分析场插值到观测站点上,模拟观测站点选择参考真实观测站分布(如图 1所示),观测站数为1100个左右,扰动分析场由ETKF方法计算的14个扰动场与控制预报初值相加获得。
为确保所分析的GRAPES区域集合预报初值扰动结构和扰动增长特征具有较好的可信度,分别在10 km和15 km两种水平分辨率下进行ETKF初值扰动集合预报试验,预报时效为72 h,预报范围覆盖整个中国区域(15°~55°N,70°~140°E)。试验时段为2015年6月1日12:00(协调世界时,下同)至6月15日12:00,其中6月1日至10日为ETKF分析扰动α参数计算适应调整期,不参与评估检验。
3 ETKF分量特征由2.1节可知,ETKF方法产生分析扰动向量Za的数学计算公式有两个分量,一是由公式(4)计算变换矩阵T,二是由公式(7)计算放大因子。变换矩阵T是由12 h预报扰动向量计算产生,因此首先分析ETKF产生的分析扰动特征向量矩阵是否具有正交性,下面通过对10 km与15 km水平分辨率下集合分析扰动特征值谱的分析,来说明GRAPES-REPS集合预报系统ETKF扰动的正交性特征。
3.1 沿正交方向的误差方差分布图 2给出了2015年6月11日00:00至15日12:00(每日两次,共计10次)平均的12 h预报扰动与分析扰动协方差矩阵的特征值谱分布(由第一个特征值进行归一化处理得到)。归一化特征值谱分布越平缓说明大量的特征向量参与描述预报不确定性,反之,特征值谱分布越陡峭说明只有少量特征向量参与描述预报的不确定性。从图中可看出,10 km(图 2a)与15 km(图 2b)试验的模式12 h预报扰动特征值谱分布较陡峭,说明特征值相差较大,而分析扰动协方差矩阵的特征值分布较为平缓,特征值大小基本相当。故通过集合卡尔曼滤波后,集合方差可以平均分配到集合子空间的全部特征向量上,保证每个扰动成员均参与描述误差方差,ETKF初值扰动方案产生的分析扰动在标准化观测空间中相互正交且不相关,15 km和10 km试验中集合预报均具有这一特性。
放大因子Πn是由不同循环时间αn参数确定的,以确保集合离散度与集合平均均方根误差量级相当。在ETKF集合扰动计算初期,由于初始集合离散度与预报均方根误差量级不协调,造成αn参数值很不稳定;经过约3~4天的扰动循环计算后,初始集合离散度与预报均方根误差量级逐渐趋于协调,αn参数逐渐稳定于0.8~1.2左右,放大因子的分布也趋于稳定。图 3给出了10 km与15 km试验中,αn参数值(图 3a)及放大因子Πn(图 3b)随时间的演变,初始时刻的αn参数值及放大因子Πn均为1。由图可知,在2015年6月1日12:00至4日12:00试验时段(即初始扰动适应调整期),αn参数值很不稳定,且10 km试验的参数值大于15 km试验的,这是因为在此过程中,与15 km试验相比,10 km试验的集合平均均方根误差显著增长,而集合离散度增长幅度较小,故由公式(8)可知,10 km试验的均方根误差与离散度比值(α参数值)大于15 km试验的,故10 km试验稳定后的放大因子较大。之后,15 km和10 km分辨率的αn参数值均维持在0.8~1.2之间,平均值接近于1,放大因子Πn也趋于稳定,15 km试验的放大因子维持在30~40之间,10 km试验的放大因子维持在40~50之间,故不同分辨率的ETKF初值方案α参数及放大因子均具有较好的稳定性。
合理的初值扰动结构及扰动增长特征对预报不确定性的准确描述和产生较好的概率预报技巧有重要的影响(Palmer, 1993;Magnusson et al., 2008),以下通过对2015年6月11至15日10 km及15 km水平分辨率下集合预报扰动方差准确率、动能谱增长特征、扰动能量演变、日变化特征及扰动一致性演变的分析,来揭示GRAPES-REPS集合预报系统ETKF初值扰动增长特征。
4.1 集合预报扰动方差准确率集合预报扰动方差准确率是集合预报方差与控制预报误差方差的比值,集合预报方差指的是集合成员与集合平均差异的平方,控制预报误差方差是控制预报与相应分辨率的GRAPES模式格点分析场差异的平方。集合预报扰动方差准确率描述了集合扰动对预报误差的捕捉能力(Wang and Bishop, 2003),其具体计算方法如下:首先,计算某一变量在某一层次上所有格点的控制预报误差方差及集合预报方差,然后,将所有格点的集合预报方差值由小到大排序,分为50个等级,对每一等级的控制预报误差方差和集合预报方差分别取平均构成散点图,散点图越接近对角线,集合扰动越能更好捕捉预报误差的结构。
图 4给出了2015年6月11日00:00至15日12:00(每日两次,共计10次)平均的850 hPa温度与速度集合预报扰动方差准确率散点图。集合扰动方差准确率散点图的斜率大于0,小的集合预报方差对应小的控制预报误差方差,反之亦然,说明集合扰动可以有效地捕捉预报误差的结构,且对于6 h(图 4a)和60 h(图 4b)预报的850 hPa温度而言,10 km试验的温度扰动方差准确率较15 km试验的更接近于对角线,表明10 km试验的温度扰动方差具有更高的准确率。对850 hPa速度而言,尽管10 km试验与15 km试验在6 h预报的扰动方差准确率相当(图 4c),但随预报时效的延长,在60 h预报时效内(图 4d),10 km试验的集合扰动方差准确率略优于15 km试验。850 hPa温度变量的控制预报误差方差在预报初期远大于集合预报方差(图 4a),随着预报时效延长,控制预报误差方差及集合预报方差均呈增加趋势(图 4b),但集合扰动方差准确率散点图斜率更大,说明集合预报方差小于控制预报误差方差在预报后期更加明显,即集合预报方差发散度不足。850 hPa风场的控制预报误差方差虽然在预报初期小于集合预报方差(图 4c),但预报后期也明显大于集合预报方差(图 4d),其他层次及变量的集合扰动方差准确率散点图分布类似,均存在集合预报方差(即集合离散度平方)随预报时效延长发散度不足的问题,不再赘述。
图 5是2015年6月11日00:00开始预报(间隔12 h)的850 hPa温度预报扰动方差准确率散点图,图 5a、c、e分别是12 h、36 h和60 h的预报扰动方差准确率,均是夜间12:00的预报,而图 5b、d、f则对应清晨00:00的预报。从图可见,10 km分辨率的温度扰动方差准确率较15 km的获得明显的改善,特别是夜晚12:00的预报改进更为显著,而对清晨00:00的影响相对较小。
总体而言,GRAPES-REPS集合扰动可以有效地捕捉预报误差的结构,但普遍存在集合离散度不足的问题。此外,提高模式水平分辨率可以改善温度场和风场的集合扰动准确率,且夜晚低层温度扰动方差准确率的改进略优于清晨。
4.2 动能谱增长特征动能谱分析可以有效地描述数值模式所能识别的尺度信息(郑永骏等,2008)。本文采用二维离散余弦变换(2-Dimensional Discrete Cosine Transform,2D-DCT)(Denis et al., 2002)对不同水平分辨率模式的扰动场进行谱分解。图 6给出了2015年6月11日00:00至15日12:00(每日两次,共计10次)平均的200 hPa、500 hPa和850 hPa动能谱分布。从图中可看出,无论是10 km试验还是15 km试验,集合扰动波谱能量峰值区均位于2000~5000 km的大尺度扰动范围内,表明ETKF初值扰动方法生成的扰动场以大尺度扰动为主,中小尺度扰动能量相对较小。随着预报时效增加,扰动波谱能量呈增加趋势,如以15 km试验为例,在0 h、24 h、48 h预报时效内,850 hPa的800 km波长对应的扰动波谱能量分别为128.96 K2、301.5 K2、502.33 K2,850 hPa的150 km波长对应的扰动波谱能量分别为9.68 K2、22.32 K2、29.74 K2,其他层次及波长的扰动波谱能量变化趋势类似,说明GRAPES-REPS可以较好的捕捉到不同尺度波动预报不确定性。另外,进一步对比10 km(红色)与15 km(蓝色)模式在不同层次的动能谱分布特征,可以看出,10 km模式在中α尺度及大尺度波长范围(400~10000 km)内的扰动波谱能量明显高于15 km模式,尤以中高层大尺度扰动波谱能量增加最为显著,说明水平分辨率提高可以增加中高层大尺度扰动波谱能量,对中尺度略有改进,但由于850 hPa的主要天气系统为低涡、低空急流等中尺度系统,水平分辨率提高对850 hPa预报后期大尺度扰动波谱能量提高不显著。
综上所述,GRAPES-REPS的ETKF初值扰动方法生成的扰动场以大尺度扰动为主,水平分辨率提高有助于增加中高层大尺度扰动的波谱能量,对中尺度略有改进,且随预报时效延长,扰动波谱能量逐渐发展,能够较好的捕捉到不同尺度波动预报不确定性。
4.3 扰动能量演变扰动能量演变可以较好地体现出预报误差随预报时效的变化。为说明集合扰动总能量的演变特征,引用Palmer et al.(1998)提出的适用于资料同化及天气预报研究的集合扰动总能量进行分析,定义为
$ P\left({i, j, k} \right) = \frac{1}{2}\left[ {{{u'}^2}\left({i, j, k} \right) + {v^2}\left({i, j, k} \right)} \right] + \frac{{{c_p}}}{{{T_r}}}{T'^2}\left({i, j, k} \right), $ | (9) |
其中,u'、v'、T'分别代表纬向风U、经向风V、温度T的扰动,扰动值是各个集合成员的预报与集合平均的差值,Tr是参考温度,cp是干空气的定压比热,i、j、k分别为水平及垂直格点数。
通过计算所有扰动成员在各个垂直层次上的区域平均扰动总能量,再对集合成员求平均便可得到10 km及15 km模式下不同预报时效的扰动总能量垂直分布图。图 7给出了2015年6月11日00:00至15日12:00(每日两次,共计10个时次)平均的扰动总能量垂直分布。从图中可看出,10 km与15 km模式各层扰动总能量随预报时效的延长均呈发展趋势,能够表示预报误差的增长,以10 km模式为例,在250 hPa高度层上,0 h预报时效的扰动总能量近似为20 J kg-1,相应的72 h预报时效的扰动总能量为38 J kg-1,能量增幅为18 J kg-1。在800 hPa高度层上,0 h预报时效的扰动总能量近似为9 J kg-1,相应的72 h预报时效的扰动总能量为18 J kg-1,能量增幅为9 J kg-1,其他层次扰动总能量变化趋势类似。其次还可看到,10 km与15 km模式扰动总能量在垂直方向分布不均匀,扰动总能量最大值位于250 hPa,次大值位于低层925 hPa,650 hPa的扰动总能量小于其余层次。最后我们可以注意到,在同样的预报时效内,10 km模式各垂直层次的扰动总能量均大于15 km模式相应层次的扰动总能量,说明10 km模式的集合成员离散度更大,更有可能包含真实大气状态。
为了进一步说明GRAPES-REPS集合扰动总能量在低层和高层分量组成,图 8给出了2015年6月11日00:00至15日12:00(每日两次,共计10个时次)平均的动能扰动能量(虚线)与动能和内能扰动总能量(实线)在0 h(图 8a)和72 h(图 8b)预报时效内的垂直分布。从图中可知,10 km与15 km模式的扰动总能量最大值均位于250 hPa,即急流轴附近,扰动总能量以动能为主,与Bowler et al.(2008)的研究结果相吻合;扰动总能量的次大值位于近地面,在此高度上,内能扰动能量占据主导地位,这可能与地面长波辐射、感热、潜热加热引起温度变化明显有关。故扰动总能量在低层包含的内能成分较多,在高层包含的动能成分较多。
下面我们分析扰动总能量水平分布特征。图 9给出了2015年6月11日00:00的集合成员平均500 hPa扰动总能量在0 h和72 h预报时效的水平分布。从图中可以看出,两种水平分辨率模式下,在动力不稳定区域,如0 h预报时效对应的内蒙古高原大槽及72 h预报时效对应的贝加尔湖南侧及青藏高原东北侧的脊区,扰动总能量量级较大,而在等位势线分布较为平直的地方,扰动总能量量级相对较小,说明500 hPa扰动总能量具有明显的流型依赖特性,其他层次的扰动总能量分布类似,在这里不再赘述。另外,10 km模式的扰动总能量(图 9a、b)明显大于15 km模式相应时效的扰动总能量(图 9c、d),故水平分辨率提高可以增加不稳定区域的扰动总能量,对天气系统的响应更明显,动力学特征更显著。
为了进一步说明扰动结构的水平分布特征,图 10给出了2015年6月11日00:00的10 km模式第4个集合成员500 hPa扰动量u'、v'及T'在0 h和72 h预报时效的水平分布。从图中可见,无论是风场还是温度场的扰动,扰动量大值区(深色阴影区域)均位于槽脊区(即动力不稳定区),而在等位势线分布较为平直的地方,扰动量级较小(浅色阴影区域),说明扰动结构具有明显的随流型依赖特征,与上述扰动总能量水平分布特征相吻合。
4.3节得出5天(10个时次)平均各层次扰动总能量随预报时效延长呈发展趋势,那么,某一时次(如2015年6月11日00:00)各层扰动总能量随预报时效变化规律均是一致的吗?下面就此问题展开进一步讨论。图 11给出了2015年6月11日00:00开始预报的纬向平均第3层模式面上集合成员平均扰动总能量的经度—时间剖面。从图中可见,两种水平分辨率模式的低层扰动总能量随预报时效延长整体呈发展趋势,但依然存在明显的日变化特征,尤以我国西部地区(青藏高原)日变化特征最为显著,在12 h、36 h和60 h[对应12:00(青藏高原当地时间17:30PM)的预报]均对应着扰动总能量大值区,如36 h和60 h的扰动总能量最大值超过18 J kg-1,而在0 h、24 h、48 h、72 h[对应00:00(青藏高原当地时间5:30AM)的预报]扰动总能量较小,可见低层扰动总能量具有下午增大、清晨减小的日变化特征,而中高层扰动总能量不存在这种日变化特征(图略)。图 8显示低层扰动总能量主要是由扰动内能组成,故低层扰动总能量日变化主要是由扰动内能日变化引起的,由于下午对流活动较旺盛,集合成员相对集合平均低层温度分布较为发散,温度扰动值较大,故12:00的扰动内能较大,而在后半夜,大气中对流相对下午偏弱,集合成员相对集合平均低层温度分布较为集中,温度扰动值较小,故00:00的扰动内能较小,这种低层扰动内能日变化特征与青藏高原本身温度日变化特征相吻合。
集合预报平均是集合预报成员的数学平均,代表了集合预报的中心分布;均方根误差(RMSE)是指集合预报与相应分辨率的GRAPES模式格点分析场之间的差异;集合离散度(ensemble spread)是衡量预报成员发散度的指标,即集合成员相对集合平均的变化幅度。好的集合预报系统,集合平均均方根误差与离散度的比值(即一致性)要近似等于1,这样才能确保预报集合最大可能的涵盖大气真实状态。图 12是2015年6月11日00:00至15日12:00(每日两次,共计10个时次)平均的200 hPa、500 hPa和850 hPa温度T与纬向风U的控制预报均方根误差、集合平均均方根误差和集合离散度随预报时效的演变。从图 12中可看出,在各个预报时效内,集合平均的均方根误差均小于控制预报均方根误差,即集合预报对确定性预报具有一定的改善作用。另外,随着预报时效的延长,10 km与15 km模式等压面要素的集合离散度及均方根误差均呈增长趋势,如对10 km模式的850 hPa温度而言,在6 h到72 h预报时效内,集合平均均方根误差从1.25 K增加到2.2 K,集合离散度从1.1 K增加到1.6 K,这种增长趋势能在一定程度上体现预报的不确定性,但我们可以明显地看出,集合离散度增长率小于均方根误差增长率,故集合离散度随预报时效的演变存在发散度不足的问题,这与4.1节结论相一致。
进一步对比10 km与15 km模式不同高度层等压面要素的均方根误差和离散度差异可发现,对低层等压面气象要素850 hPa温度(图 12c)、850 hPa纬向风U(图 12f)而言,10 km模式的控制预报与集合平均均方根误差小于15 km模式,其中以850 hPa温度变化最为显著,与麻巨慧(2012)的结论相符,与此同时,10 km模式的集合离散度大于15 km模式,故10 km GRAPES-REPS区域集合预报低层要素的集合预报一致性更接近于1,集合预报包含实况的可能性更高。而对中高层等压面要素,如200 hPa温度(图 12a)、200 hPa纬向风U(图 12d)、500 hPa温度(图 12b)、500 hPa纬向风U(图 12e)而言,在各个预报时效内,水平分辨率提高使得集合离散度均得到较大的增长,而相应的集合平均均方根误差基本保持不变或略有增长,集合预报一致性更接近于1,故水平分辨率提高可以改进温度和风场的预报性能,其他变量的变化类似,在这里不再赘述。
4.5.2 近地面要素图 13是2015年6月11日00:00至15日12:00(每日两次,共计10个时次)平均的2 m温度与10 m风的控制预报、集合平均均方根误差、集合离散度随预报时效的演变。从图中可看出,10 km与15 km模式的集合离散度与均方根误差随预报时效延长均呈增长趋势,能够表示预报的不确定性,同时集合离散度增长率小于均方根误差增长率,即集合离散度不足。另外,与15 km模式相比,10 km模式的集合预报一致性更接近于1,故水平分辨率提高对近地面温度和风场具有明显的改进,同时,陆面是大气热量、动量、水汽的源汇项,近地面要素的改善可以通过陆气相互作用影响大气,进而影响模式的预报效果。
综上所述,GRAPES-REPS的集合均方根误差及集合离散度增长趋势可以表示预报的不确定性,但存在集合离散度增长率不足的问题。另外,水平分辨率提高可以明显地改进等压面及近地面温度场及风场的集合预报效果。
5 结论与讨论中国气象局数值预报中心建立了GRAPES区域集合预报业务系统,初值扰动采用ETKF方法,为进一步认识ETKF初值扰动在GRAPES区域集合预报系统中的结构和演变特征,改进GRAPES区域集合预报效果,本文利用15 km和10 km两种水平分辨率的GRAPES-REPS区域集合预报模式,对2015年6月1~15日进行集合预报初值扰动模拟试验,通过对ETKF初值扰动分量、集合预报扰动方差准确率、动能谱分布、扰动能量演变、日变化、集合离散度及均方根误差等特征进行分析,以揭示GRAPES区域集合预报ETKF初值扰动分量、扰动结构与扰动增长特征,得出如下结论:
(1) ETKF初值扰动方法产生的分析扰动在标准化观测空间中相互正交且不相关,且初值扰动方案α参数值及放大因子具有较好的稳定性,初值扰动分量特征较为合理。
(2) ETKF初值扰动方法能够较好的捕捉预报误差的结构,反应不同尺度波动预报不确定性,初值扰动结构以大尺度扰动为主,且低层以内能扰动为主,高层以动能扰动为主,并随预报时效延长,具有明显的流依赖变化动力学特征,扰动能量逐渐发展。
(3) ETKF初值扰动总能量和集合离散度随预报时效的演变呈发展趋势,能够代表预报误差的增长,但离散度增长率小于均方根误差增长率。
(4) 提高模式水平分辨率可以增加中高层大尺度扰动波谱能量,提高等压面及近地面温度场及风场的集合概率预报技巧。
值得关注的是,尽管GRAPES-REPS ETKF初值扰动结构与扰动增长特征较为合理,但仍然存在离散度不足的问题,且在青藏高原地区低层内能扰动总能量存在着明显的日变化特征,需要进一步研究青藏高原初值扰动结构的合理性。
Bishop C H, Etherton B J, Majumdar S J. 2001. Adaptive sampling with the ensemble transform Kalman filter. Part Ⅰ:Theoretical aspects [J]. Mon. Wea. Rev., 129(3): 420-436. DOI:10.1175/1520-0493(2001)129<0420:ASWTET>2.0.CO;2
|
Bowler N E, Arribas A, Mylne K R, et al. 2008. The MOGREPS short-range ensemble prediction system [J]. Quart. J. Roy. Meteor. Soc., 134(632): 703-722. DOI:10.1002/qj.234
|
Bowler N E, Arribas A, Beare S E, et al. 2009. The local ETKF and SKEB:Upgrades to the MOGREPS short-range ensemble prediction system [J]. Quart. J. Roy. Meteor. Soc., 135(640): 767-776. DOI:10.1002/qj.394
|
Buizza R, Palmer T N. 1995. The singular-vector structure of the atmospheric global circulation [J]. J. Atmos. Sci., 52(9): 1434-1456. DOI:10.1175/1520-0469(1995)052<1434:TSVSOT>2.0.CO;2
|
Buizza R, Houtekamer P L, Pellerin G, et al. 2005. A comparison of the ECMWF, MSC, and NCEP global ensemble prediction systems [J]. Mon. Wea. Rev., 133(5): 1076-1097. DOI:10.1175/MWR2905.1
|
Buizza R, Leutbecher M, Isaksen L. 2008. Potential use of an ensemble of analyses in the ECMWF ensemble prediction system [J]. Quart. J. Roy. Meteor. Soc., 134(637): 2051-2066. DOI:10.1002/qj.346
|
陈静, 薛纪善, 颜宏. 2003. 物理过程参数化方案对中尺度暴雨数值模拟影响的研究[J]. 气象学报, 61(2): 203-218. Chen Jing, Xue Jishan, Yan Hong. 2003. The impact of physics parameterization schemes on mesoscale heavy rainfall simulation (in Chinese)[J]. Acta Meteorologica Sinica, 61(2): 203-218. DOI:10.11676/qxxb2003.019
|
Denis B, Côté J, Laprise R. 2002. Spectral decomposition of two-dimensional atmospheric fields on limited-area domains using the discrete cosine transform (DCT) [J]. Mon. Wea. Rev., 130(7): 1812-1829. DOI:10.1175/1520-0493(2002)130<1812:SDOTDA>2.0.CO;2
|
Houtekamer P L, Lefaivre L, Derome J, et al. 1996. A system simulation approach to ensemble prediction [J]. Mon. Wea. Rev., 124(6): 1225-1242. DOI:10.1175/1520-0493(1996)124<1225:ASSATE>2.0.CO;2
|
黄红艳, 齐琳琳, 刘健文, 等. 2016. 多物理ETKF在暴雨集合预报中的初步应用[J]. 大气科学, 40(4): 657-668. Huang Hongyan, Qi Linlin, Liu Jianwen, et al. 2016. Preliminary application of a multi-physical ensemble transform Kalman filter in precipitation ensemble prediction (in Chinese)[J]. Chinese Journal of Atmospheric Sciences, 40(4): 657-668. DOI:10.3878/j.issn.1006-9895.1508.14308
|
Kay J K, Kim H M. 2014. Characteristics of initial perturbations in the ensemble prediction system of the Korea Meteorological Administration [J]. Wea. Forecasting., 29(3): 563-581. DOI:10.1175/WAF-D-13-00097.1
|
Kim H M, Morgan M C. 2002. Dependence of singular vector structure and evolution on the choice of norm [J]. J. Atmos. Sci., 59(21): 3099-3116. DOI:10.1175/1520-0469(2002)059<3099:DOSVSA>2.0.CO;2
|
Leith C E. 1974. Theoretical skill of Monte Carlo forecasts [J]. Mon. Wea. Rev., 102(6): 409-418. DOI:10.1175/1520-0493(1974)102<0409:TSOMCF>2.0.CO;2
|
龙柯吉, 陈静, 马旭林, 等. 2011. 基于集合卡尔曼变换的区域集合预报初步研究[J]. 成都信息工程学院学报, 26(1): 37-46. Long Keji, Chen Jing, Ma Xulin, et al. 2011. The preliminary study on ensemble prediction of GRAPES-meso based on ETKF (in Chinese)[J]. Journal of Chengdu University of Information Technology, 26(1): 37-46. DOI:10.3969/j.issn.1671-1742.2011.01.008
|
麻巨慧. 2012. NCEP全球集合预报系统框架的优化研究[D]. 南京信息工程大学博士学位论文. Ma Juhui. 2012. Research on optimizing the configuration of NCEP global ensemble forecast system[D]. Ph. D. dissertation (in Chinese), Nanjing University of Information Science and Technology.
|
马旭林, 薛纪善, 陆维松. 2008. GRAPES全球集合预报的集合卡尔曼变换初始扰动方案初步研究[J]. 气象学报, 66(4): 526-536. Ma Xulin, Xue Jishan, Lu Weisong. 2008. Preliminary study on ensemble transform Kalman filter-based initial perturbation scheme in GRAPES global ensemble prediction (in Chinese)[J]. Acta Meteorologica Sinica, 66(4): 526-536. DOI:10.11676/qxxb2008.050
|
马旭林, 时洋, 和杰, 等. 2015. 基于卡尔曼滤波递减平均算法的集合预报综合偏差订正[J]. 气象学报, 73(5): 952-964. Ma Xulin, Shi Yang, He Jie, et al. 2015. The combined descending averaging bias correction based on the Kalman filter for ensemble forecast (in Chinese)[J]. Acta Meteorologica Sinica, 73(5): 952-964. DOI:10.11676/qxxb2015.062
|
Magnusson L, Källén E, Nycander J. 2008. Initial state perturbations in ensemble forecasting [J]. Nonlinear Processes in Geophysics, 15(5): 751-759. DOI:10.5194/npg-15-751-2008
|
Molteni F, Buizza R, Palmer T N, et al. 1996. The ECMWF ensemble prediction system:Methodology and validation [J]. Quart. J. Roy. Meteor. Soc., 122(529): 73-119. DOI:10.1002/qj.49712252905
|
Palmer T N. 1993. Extended-range atmospheric prediction and the Lorenz model [J]. Bull. Amer. Meteor. Soc., 74(1): 49-66. DOI:10.1175/1520-0477(1993)074<0049:ERAPAT>2.0.CO;2
|
Palmer T N, Gelaro R, Barkmeijer J, et al. 1998. Singular vectors, metrics, and adaptive observations [J]. J. Atmos. Sci., 55(4): 633-653. DOI:10.1175/1520-0469(1998)055<0633:SVMAAO>2.0.CO;2
|
Saito K, Seko H, Kunii M, et al. 2012. Effect of lateral boundary perturbations on the breeding method and the local ensemble transform Kalman filter for mesoscale ensemble prediction [J]. Tellus A:Dynamic Meteorology and Oceanography, 64(1): 11594. DOI:10.3402/tellusa.v64i0.11594
|
Steven Tracton M, Kalnay E. 1993. Operational ensemble prediction at the national meteorological center:Practical aspects [J]. Wea. Forecasting., 8(3): 379-400. DOI:10.1175/1520-0434(1993)008<0379:OEPATN>2.0.CO;2
|
Torn R D. 2010. Performance of a mesoscale ensemble Kalman filter (EnKF) during the NOAA high-resolution hurricane test [J]. Mon. Wea. Rev., 138(12): 4375-4392. DOI:10.1175/2010MWR3361.1
|
Toth Z, Kalnay E. 1993. Ensemble forecasting at NMC:The generation of perturbations [J]. Bull. Amer. Meteor. Soc., 74(12): 2317-2330. DOI:10.1175/1520-0477(1993)074<2317:EFANTG>2.0.CO;2
|
Toth Z, Kalnay E. 1997. Ensemble forecasting at NCEP and the breeding method [J]. Mon. Wea. Rev., 125(12): 3297-3319. DOI:10.1175/1520-0493(1997)125<3297:EFANAT>2.0.CO;2
|
王太微. 2008. 中尺度模式不确定性与初值扰动试验研究[D]. 中国气象科学研究院硕士学位论文. Wang Taiwei. 2008. Impact of initial perturbation on meso-scale model uncertainty[D]. M. S. thesis (in Chinese), Chinese Academy of Meteorological Sciences.
|
Wang X, Bishop C H. 2003. A comparison of breeding and ensemble transform Kalman filter ensemble forecast schemes [J]. J. Atmos. Sci., 60(9): 1140-1158. DOI:10.1175/1520-0469(2003)060<1140:ACOBAE>2.0.CO;2
|
Wang X, Bishop C H, Julier S J. 2004. Which is better, an ensemble of positive negative pairs or a centered spherical simplex ensemble? [J]. Wea. Rev., 132(7): 1590-1605. DOI:10.1175/1520-0493(2004)132<1590:WIBAEO>2.0.CO;2
|
Wei M Z, Toth Z, Wobus R, et al. 2006. Ensemble transform Kalman filter-based ensemble perturbations in an operational global prediction system at NCEP [J]. Tellus A, 58(1): 28-44. DOI:10.1111/j.1600-0870.2006.00159.x
|
Wei M Z, Toth Z, Wobus R, et al. 2008. Initial perturbations based on the ensemble transform (ET) technique in the NCEP global operational forecast system [J]. Tellus A, 60(1): 62-79. DOI:10.1111/j.1600-0870.2007.00273.x
|
肖玉华, 何光碧, 陈静, 等. 2011. 区域集合预报增长模繁殖扰动方法研究[J]. 高原气象, 30(1): 94-102. Xiao Yuhua, He Guangbi, Chen Jing, et al. 2011. Study on perturbing method in regional BGM ensemble prediction system (in Chinese)[J]. Plateau Meteorology, 30(1): 94-102.
|
张涵斌, 陈静, 智协飞, 等. 2014a. 基于GRAPES_Meso的集合预报扰动方案设计与比较[J]. 大气科学学报, 37(3): 276-284. Zhang Hanbin, Chen Jing, Zhi Xiefei, et al. 2014a. Design and comparison of perturbation schemes for GRAPES_Meso based ensemble forecast (in Chinese)[J]. Transactions of Atmospheric Sciences, 37(3): 276-284. DOI:10.3969/j.issn.1674-7097.2014.03.003
|
张涵斌, 陈静, 智协飞, 等. 2014b. GRAPES区域集合预报系统应用研究[J]. 气象, 40(9): 1076-1087. Zhang Hanbin, Chen Jing, Zhi Xiefei, et al. 2014b. Study on the application of GRAPES regional ensemble prediction system (in Chinese)[J]. Meteor. Mon., 40(9): 1076-1087. DOI:10.7519/j.issn.1000-0526.2014.09.005
|
Zhang Hanbin, Chen Jing, Zhi Xiefei, et al. 2015. Study on multi-scale blending initial condition perturbations for a regional ensemble prediction system [J]. Advances in Atmospheric Sciences, 32(8): 1143-1155. DOI:10.1007/s00376-015-4232-6
|
郑永骏, 金之雁, 陈德辉. 2008. 半隐式半拉格朗日动力框架的动能谱分析[J]. 气象学报, 66(2): 143-157. Zheng Yongjun, Jin Zhiyan, Chen Dehui. 2008. Kinetic energy spectrum analysis in a semi implicit semi Lagrangian dynamical framework (in Chinese)[J]. Acta Meteorologica Sinica, 66(2): 143-157. DOI:10.11676/qxxb2008.015
|
庄潇然, 闵锦忠, 王世璋, 等. 2017. 风暴尺度集合预报中的混合初始扰动方法及其在北京2012年"7.21"暴雨预报中的应用[J]. 大气科学, 41(1): 30-42. Zhuang Xiaoran, Min Jinzhong, Wang Shizhang, et al. 2017. A blending method for storm-scale ensemble forecast and its application to Beijing extreme precipitation event on July 21, 2012 (in Chinese)[J]. Chinese Journal of Atmospheric Sciences, 41(1): 30-42. DOI:10.3878/j.issn.1006-9895.1605.15233
|