大气科学  2015, Vol. 39 Issue (6): 1225-1236   PDF    
热带风压场平衡特征及其对GRAPES系统中同化预报的影响研究Ⅱ:动力与统计混合平衡约束方案的应用
王瑞春1,2, 龚建东2 , 张林2, 陆慧娟2    
1 南京信息工程大学大气科学学院, 南京 210044;
2 中国气象局数值预报中心, 北京 100081
摘要: 研究I的结果表明:线性平衡方程(LBE)在热带地区不适用,而进一步改进方向是削弱LBE在该区域的约束程度。本文以此为基础,在GRAPES(global/regional assimilation and prediction system)全球变分同化系统中引入动力与统计混合平衡约束方案。新方案在逐层求解LBE的基础上增加垂直方向的线性回归,回归系数随纬度和高度变化。针对背景误差协方差的分析表明,新方案可以更好的保证独立分析变量间预报误差不相关的基本要求,并大幅度减小热带地区平衡气压预报误差方差的量值和占总方差的比例。单点试验结果表明,与LBE方案相比,新方案对中、高纬影响很小,但在热带地区成功实现了风、压场分析的解耦,两者分析更为独立。并且,虽未考虑具体波动模态,但新方案给出的风、压场协相关结构与研究I的理论分析结果相近。一个月的同化循环与预报结果表明,引入新方案后,赤道外地区的同化预报效果为中性偏正,而热带地区风场的同化预报效果显著提高,LBE方案中平流层低层的风场同化预报异常被基本消除。
关键词: 变分资料同化     平衡约束     线性平衡方程     线性回归     GRAPES(global/regional assimilation and prediction system)    
Tropical Balance Characteristics between Mass and Wind Fields and Their Impact on Analyses and Forecasts in GRAPES System. Part Ⅱ:Application of Linear Balance Equation-Regression Hybrid Constraint Scheme
WANG Ruichun1,2, GONG Jiandong2 , ZHANG Lin2, LU Huijuan2    
1 School of Atmospheric Science, Nanjing University of Information Science & Technology, Nanjing 210044;
2 Numerical Weather Prediction Center of China Meteorological Administration, Beijing 100081
Abstract: The results of Part I of this study show that the imposition of the Linear Balance Equation(LBE) in the tropics is not correct and additional steps should be taken to decouple the analysis of mass and wind fields. Based on the results in Part I, this paper developed a LBE-regression hybrid balance constraint in GRAPES(global/regional assimilation and prediction system) global variational data assimilation system. In the new scheme, after the calculation of LBE on each level, a vertical regression whose coefficients could vary at different latitudes and model levels was introduced. By analyzing and comparing the implied background error covariance of different schemes, the authors found that the new scheme could efficiently reduce the unreasonable correlations of control variables. Additionally, the forecast error variance of the balanced pressure and its ratio to the whole pressure were also both dramatically reduced in the tropics. Results of the single-observation experiments indicated that, the new scheme had little impacts on the mid- and high-latitude compared to the LBE scheme, but it did successfully decouple the mass/wind analyses in the tropics. Although equatorial wave modes were not explicitly considered in the new scheme, the structure of covariance was consist with the theory analysis results based on those modes in Part I. The results of cycle analysis and forecast experiments in one month showed that, the new scheme could bring slightly positive results in the extratropics and significantly improve the wind analysis and forecast accuracy in the tropics, the abnormally large tropical wind errors in the LBE scheme were dramatically suppressed.
Key words: Variational data assimilation     Balance constraint     Linear balance equation     Linear regression     GRAPES    
1 引言

研究I的结果表明,GRAPES(global/regional assimilation and prediction system)全球变分资料 同化系统(以下简称GRAPES-VAR)采用的风、压场平衡约束——线性平衡方程(Linear Balance Equation,以下简称LBE)在热带地区并不适用,会造成虚假平衡(王瑞春等,2015)。LBE主要表达了罗斯贝波模态下的风、压场配置,但该模态在热带区域的短期预报误差中并不占主导,补充考虑其他赤道波动的影响就需要削弱LBE对热带风、压场的约束程度,使得两者的分析变得更加独立。作为研究的第II部分,本文致力于在GRAPES-VAR中引入更加合理的平衡约束方案,对热带虚假平衡问题做针对性修正,以提高该区域的风场分析效果。

热带地区缺少类似中、高纬准地转这样的主导机制(Holton,1992),基于自身动力学特征构造平衡约束的研究进展十分缓慢。目前,这方面研究的一个主要思路是采用相互正交的赤道波动作为特征分量构造预报(背景)误差协方差矩阵(B矩阵)(Daley,1993; Žagar et al.,2004; Körnich and Källén,2008),也即与研究I中构造风、压场协相关的方案类似。然而,由于在波动权重的确定、误差垂直结构的设定以及如何与中、高纬衔接等问题上仍不十分清晰,这些研究仅局限于在正压浅水模型中做一些理想试验,距离业务应用有很大距离。另一个努力方向是采用比LBE更加复杂的平衡方程,例如非线性平衡方程。由于考虑了流场曲率的作用,非线性平衡方程在处理强旋转系统(例如台风)时要更加精确(Fisher,2003; 庄照荣等,2006; 万齐林和薛纪善,2007),但由于其仍是基于中、高纬大气运动的尺度特征简化得到,因而在热带也是不适用的(Daley,1991; Barker et al.,2004)。

在业务同化系统中,为减小热带风、压场的平衡约束,一个可行方法是采用统计方案。该方案采用线性回归直接统计不同变量间气候态的平衡约束,其引入的最初目标是解决Lorenz垂直离散方案下求解温度的欠定问题(Parrish et al.,1997)。Derber and Bouttier(1999)将该方法应用于欧洲中期数值预报中心(ECMWF)变分同化系统时发现,风、压场分析在中、高纬与LBE相近,而热带地区两者分析接近独立。由于构造和计算简单,统计方案在业务和研究中得以广泛应用(Berre,2000; Wu et al.,2002; Huang et al.,2009; Kleist et al.,2009a)。王瑞春等(20122014)在将该方案应用于GRAPES-VAR的研究中也发现,统计得到的平衡约束在热带地区要远小于LBE,可以帮助减小虚假平衡问题。

然而,Kleist et al.(2009b)基于NCEP(National Centers for Environmental Prediction)变分同化系统的研究指出,统计得到的平衡约束在中、高纬并不如LBE稳健。此外,统计方案只能考虑线性约束,无法进一步向非线性扩展。因此,英国气象局和ECMWF在其变分同化系统中将平衡方程和统计方案结合到一起,构造了动力与统计混合平衡约束方案(以下简称“混合方案”)(Lorenc et al.,2000Fisher,2003)。混合方案先逐层求解线性或非线性平衡方程,并在此基础上引入垂直方向的回归统计。平衡方程的使用保证了约束在中、高纬的稳健性,并且还可以向非线性平衡方程扩展;而统计模块的引入可以减小平衡方程不适用地区的虚假平衡(Barker et al.,2004)。鉴于混合方案具有取长补短的优势,并有进一步升级的潜力,本文针对GRAPES-VAR的改进也选用该方案。

虽然统计方案与混合方案已有不少业务应用先例,但之前研究主要集中于方案的具体构造以及业务性能评估上,而对方案减小热带虚假平衡问题的内在机制缺乏详细讨论。变分同化中,平衡约束主要是在针对B矩阵的物理变换,即独立分析变量的构造中引入的(Bannister,2008)。LBE应用于热带时,其虚假平衡问题对独立分析变量以及B矩阵的构造究竟会产生什么样的不利影响,引入统计模块后能否避免,目前国内外研究中鲜有这方面的详细分析。此外,虽然统计方案与混合方案得到的热带风、压场约束远小于LBE给定的值,但统计过程中并未细化考虑不同赤道波动的影响,那么其结果能否与基于赤道波动模态的理论分析相符,也值得详细对比与分析。针对该问题,Žagar et al.(2004)做过单个个例的对比,其结果表明ECMWF变分系统的混合方案在对流层高层与波动理论分析结果较为一致。为科学认识引入统计模块对构造热带平衡约束的帮助,本文在引入混合方案后,将基于GRAPES-VAR对上述两个问题做进一步详细分析。

本文第2节给出混合方案的基本框架以及统计模块实施的主要技术细节;第3节基于GRAPES短期预报误差样本,分析了混合方案对独立分析变量以及B矩阵构造的影响;第4节则利用单点理想观测试验考察了混合方案给定的热带风、压场平衡特征,并与研究I的理论分析结果作对比;第5节通过接近业务实际的同化循环与预报试验评估了混合方案对GRAPES同化预报性能的影响;最后第6节给出结论和讨论。

2 混合方案的引入 2.1 方案设计

GRAPES-VAR采用增量形式的目标函数,变量间的平衡约束在物理变换部分引入(庄世宇等,2005薛纪善等,2008)。物理变换是针对B矩阵的预条件变换之一,通过抽取分析变量间的预报误差协相关获取一组误差不相关的独立分析变量(Bannister,2008)。

目前,GRAPES-VAR的分析变量包括水平风场(uv),质量场$\pi $(Exner函数,也称无量纲气压)以及比湿q。由于纬向风u与径向风v的预报误差存在高度相关,利用赫姆霍兹速度分解定理将它们转换为旋转风和散度风,分别用流函数$\psi $和势函数$\chi $表示。此外,本文关于平衡约束的讨论不涉及水汽,为简化叙述将GRAPES-VAR的(动力学)分析变量记为$x = {\left[ {\psi ,\chi ,\pi } \right]^T}$。为构造独立分析变量,GRAPES-VAR采用旋转风$\psi $的全量表征大气运动中的平衡部分,将质量场$\pi $拆分为与$\psi $相平衡的${\pi _{\rm{b}}}$以及非平衡的${\pi _{\rm{u}}}$,并暂不考虑散度风和其它变量间的相关,也即:

$\left\{ \matrix{ \psi \;\, = \psi , \hfill \cr \chi \;\,{\kern 1pt} {\kern 1pt} = \chi , \hfill \cr {{\rm{\pi }}_{\rm{u}}} = \pi - {\pi _{\rm{b}}} = \pi - N\psi , \hfill \cr} \right.$ (1)
式中,平衡算子N表达了旋转风与质量场之间的平衡约束,也即抽取了$\psi $与$\pi $间的预报误差协相关。新的独立分析变量可以表示为${x_{\rm{u}}} = {\left[ {\psi ,\chi ,{\pi _{\rm{u}}}} \right]^T}$,同化框架中假设它们之间不再相关。上述物理变换方案与英国气象局变分系统一致(Lorenc et al.,2000),与ECMWF变分系统相比则缺少了散度风与旋转风、质量场之间的约束(Derber and Bouttier,1999),不过上述两项均是小量,它们与本研究的关系将在讨论部分给出。

在现有框架中,GRAPES-VAR采用线性平衡方程LBE求解算子N(下文将该方案称作“LBE方案”)。目前,LBE方案在模式面(地形高度追随坐标面)上直接求解LBE,表达式(具体推导参见薛纪善等,2012)为

${\pi _{{\rm{b}}1}} = {1 \over {{c_p}{{\bar \theta }_{\rm{B}}}}}{\nabla ^{ - 2}}\left[ {\nabla \cdot \left( {f\,\nabla \psi } \right)} \right]{\rm{ }},$ (2)
式中,${c_p}$是摩尔定压热容,${\bar \theta _{\rm{B}}}$是整层平均的背景场位温,f是科氏参数。为下文比较的方便,这里将LBE方案计算得到的$\pi $的平衡部分记为了${\pi _{{\rm{b}}1}}$,相应的非平衡部分记为${\pi _{{\rm{u}}1}}$。

研究I的理论分析表明,式(2)在热带地区并不适用,与实际情形相比其约束过强,因而会造成虚假平衡。为克服该问题,并保留LBE在中、高纬稳健的优点,本文引入动力与统计混合平衡约束方案。混合方案将N算子的求解分为两步实施:第一步是在模式面上利用原有的动力学推导得到的平衡方程,目前就是式(2)给出的LBE,逐层计算一个初步的平衡气压${{\rm{\pi }}_{{\rm{b}}1}}$;第二步是利用线性回归统计得到的系数R(统计方法在下文给出)对${{\rm{\pi }}_{{\rm{b}}1}}$的垂直廓线作加权求和处理,得到最终的平衡气压${{\rm{\pi }}_{{\rm{b}}2}}$,其在第p层上的值为

${\pi _{{\rm{b}}2}}\left( p \right) = \sum\limits_{k = 1}^K {R\left( {k,p} \right){\pi _{{\rm{b}}1}}} (k),$ (3)
式中,K为总的模式层数。从式(3)的表达可以看出,混合方案计算每一层上${\pi _{{\rm{b}}2}}$时用到了所有层次上${\pi _{{\rm{b}}1}}$的信息,但权重较大的主要为相邻层次上的信息(见下文)。混合方案引入后,$\pi $的非平衡部分也会发生相应改变,被记为${\pi _{{\rm{u}}2}}$。

式(3)与式(2)相结合,就构成了变分框架中混合方案的主要计算公式,其结构十分简明。不过,混合方案顺利实施的关键在于统计得到一个实用、稳健的系数R

2.2 系数R的统计

统计R的值需要获取一组短期预报误差样 本,与研究I一样仍采用NMC方法,也即采用模式预报到同一时刻,不同预报时效(48 h和24 h)的预报场的差值作为短期预报误差的替代样本(Parrish and Derber,1992)。样本获取的具体设置与研究I相同,模式水平分辨率为1°×1°,垂直层次为36层。不过与研究I的理论分析工作不同,这里的系数R会在业务系统中被反复调用,需尽可能减小统计噪音的影响。为此,本文选取一年的NMC样本,具体时间段为2009年12月至2010年11月。

NMC方法首先获取的是模式预报变量uv、$\pi $的短期预报误差样本,之后将uv转换为$\psi $,并进而采用式(2)计算得到${\pi _{{\rm{b}}1}}$的误差样本。在剔除$\pi $和${\pi _{{\rm{b}}1}}$样本中的时间平均后,逐层、逐纬圈建立如下的多元线性回归模型:

$\pi {\rm{(}}p{\rm{)}} = \sum\limits_{k = 1}^K {R(k,p){\pi _{{\rm{b}}1}}} (k) + {\pi _{{\rm{u}}2}}(p).\;$ (4)
逐层统计是考虑了大气垂直分层对平衡特征的影响,例如研究I中不同层次上赤道波动模态的比例有很大差异;逐纬圈统计是考虑科氏参数随纬度变化对平衡约束的影响,这对于区分中、高纬和热带不同的大气运动特征十分关键。

注意到,式(4)的多元回归模型中,预报因子是K个层次上的${\pi _{{\rm{b}}1}}\left( k \right)$。由于大气具有三维连续性,临近层次上${\pi _{{\rm{b}}1}}$的值一般均比较接近,也即相邻预报因子的相关系数很大。此时,统计求解会面临所谓的多重共线性问题,统计结果严重依赖具体样本,稳健性差(Wilks,2006)。为改善统计结果,这里与王瑞春等(2014)关于纯粹统计方案的研究一样,采用Lorenc et al.(2000)使用的岭回归方法(Ridge Regression method)替代经典的多元线性回归求解式(4)。关于该方法的详细的讨论可参考上述两者的工作,这里不再展开。下面我们主要分析统计得到的系数R的基本特征。

图 1a给出了在第13层上(约为532 hPa)统计得到的回归系数随纬度变化的情况。从图中可以看出,回归系数在南北半球具有较好的对称性。以北半球为例分析,系数R的一个主要特征是大值区主要集中在与13层相毗邻的几个层次上,也即混合方案计算得到的${\pi _{{\rm{b}}2}}$主要由相邻层次上的${\pi _{{\rm{b}}1}}$加权求和得到,其他层次统计得到的情形与之类似。另一个主要特征是,从中、高纬向热带地区推进时,R的值逐步减小。中、高纬地区R的数值较大,特别是在40°N北的地区;而从40°N逐渐向热带靠近的过程中,R的数值不断减小,在赤道上已几乎完全等于零。在该情形下,根据式(3),混合方案在热带计算得到的${\pi _{{\rm{b}}2}}$将远小于LBE方案中的${\pi _{{\rm{b}}1}}$。

图 1 模式(a)13层(约为532 hPa)和(b)28层(约为71 hPa)上根据式(4)统计得到的系数R的结构(等值线间距0.02) Fig. 1 Structures of R in equation (4) at model levels (a) 13 (approx. 532 hPa) and (b) 28 (approx. 71 hPa), contour intervals are 0.02

对流层大部分区域的统计结果均与图 1a相近,但到平流层低层之后,热带地区的R值出现较大变化。图 1b给出了第28层(约为71 hPa)的统计结果,其与图 1a的主要差异是赤道附近的R出现明显的负值。根据研究I,热带地区70 hPa上,赤道罗斯贝波比例最小,Kelvin波比例大大增加,风、压场平衡特征由Kelvin波主导。LBE主要表达赤道罗斯贝模态下的风、压场配置,与Kelvin波主导的情形相反,因而统计中系数R出现明显负值。这种情形在平流层低层的其他层次上也有所表现,与研究I图 3b中赤道罗斯贝波比例大幅减小的层次相对应。

另外,对比图 1两幅图可以发现,在中、高 纬区域,28层统计得到R值要小于13层的相应 值,北半球表现得更为明显。说明即使在中、高纬区域,平流层LBE的适用性也要小于对流层中层。这一现象在多个业务中心关于B矩阵的统计分析中均有表现(Ingleby,2001; Polavarapu et al.,2005),其具体原因仍有待进一步研究。

至此,我们已经给出了混合方案的基本公式和系数R的离线统计方案。根据R随纬度变化的特征来看,热带地区LBE方案求得的${\pi _{{\rm{b}}1}}$经式(3)之后会被大幅减小。那么这对独立变量以及B矩阵的构造会产生什么样的影响,下一节将具体分析。

3 混合方案对B矩阵构造的影响

本文2.1节中已经指出,风、压场的平衡约束是通过物理变换引入的,该变换将分析变量间的协相关关系(或称平衡约束)以平衡算子的形式抽取 出来,将它们转换为独立分析变量。理想的独立分析变量的预报误差间完全独立,不存在交叉协相关,但这难以实现。在业务系统中,物理变换的目标是寻找一组预报误差协相关尽可能小的分析变量,并忽略该组变量间剩余的相关,也即假设它们之间相互独立(下文称为“独立性假设”)(Parrish and Derber,1992; Derber and Bouttier,1999; 朱宗申和胡铭,2002)。这样的假设显然会带来误差,而如果物理变换抽取的平衡约束关系越合理、全面,其设计的独立分析变量间的误差协相关越小,那么假设带来的误差也就越小。

基于独立性假设,GRAPES-VAR中分析变量$x = {\left[ {\psi ,\chi ,\pi } \right]^T}$对应的B矩阵可以用独立分析变量${x_{\rm{u}}} = {\left[ {\psi ,\chi ,{\pi _{\rm{u}}}} \right]^T}$以及平衡算子N来表示:

$B = \left[ {\matrix{ {C\left( \psi \right)} & 0 & {C\left( \psi \right){N^T}} \cr 0 & {C\left( \chi \right)} & 0 \cr {NC\left( \psi \right)} & 0 & {NC\left( \psi \right){N^T} + C\left( {{\pi _u}} \right)} \cr } } \right]{\rm{,}}$ (5)
式中,$C\left( \right)$表示各独立分析变量的自协方差矩阵,$NC\left( \psi \right)$和$C\left( \psi \right){N^T}$表示旋转风和质量场间的交叉协方差矩阵,而$NC\left( \psi \right){N^T}$表示由旋转风和给定平衡约束导出的${\pi _{\rm{b}}}$的自协方差矩阵(记为$C\left( {{\pi _{\rm{b}}}} \right)$),其与$C\left( {{\pi _{\rm{u}}}} \right)$累加后构成分析变量${\rm{\pi }}$的自协方差矩阵(记为$C\left( \pi \right)$)。而根据统计学知识,如果任意拆分$\pi $为${\pi _{\rm{b}}}$与${\pi _{\rm{u}}}$两部分,那么$C\left( \pi \right)$应该由三部分组成:$C\left( {{\pi _{\rm{b}}}} \right)$、$C\left( {{\pi _{\rm{u}}}} \right)$以及${\pi _{\rm{b}}}$与${\pi _{\rm{u}}}$的交叉协方差$C\left( {{\pi _{\rm{b}}},{\pi _{\rm{u}}}} \right)$(黄嘉佑,2010)。之所以略去$C\left( {{\pi _{\rm{b}}},{\pi _{\rm{u}}}} \right)$,其依据就是独立性假设忽略了$\psi $与${\pi _{\rm{u}}}$间的相关,那么由$\psi $线性导出的${\pi _{\rm{b}}}$与${\pi _{\rm{u}}}$间也就相互独立,$C\left( {{\pi _{\rm{b}}},{\pi _{\rm{u}}}} \right)$自然为零。通过这里的分析可以看出独 立性假设对B矩阵构造十分重要,因而可以通过比较该假设成立的好坏分析不同方案的性能。

我们以对流层中层为例,图 2给出了短期预报误差样本中不同变量在北半球13层上的相关系数,均为纬圈平均值。其中,$r\left( {\psi ,\pi } \right)$给出的是分析变量$\psi $和$\pi $间的相关系数,可以看到其数值很大,需要采用物理变换将它们转换为独立分析变量。LBE方案获得的$r\left( {\psi ,{\pi _{{\rm{u1}}}}} \right)$与$r\left( {\psi ,\pi } \right)$相比,中、高纬地区的相关被成功压制,除极地以外的大部分地区已接近于零。这说明基于LBE的物理变换在中、高纬是成功的,独立性假设能较好成立。然而,随着纬度的减小,$r\left( {\psi ,{\pi _{{\rm{u}}1}}} \right)$的绝对值逐渐增大,在热带地区已超过$r\left( {\psi ,\pi } \right)$相应值,也即获取的独立分析变量间的相关反而超过了分析变量间的相关。注意到,$r\left( {\psi ,{\pi _{{\rm{u}}1}}} \right)$的值在赤道上再次接近零值,这并非说明LBE方案在赤道附近变好,而是由于$r\left( {\psi ,\pi } \right)$的值在南北半球符号相反,赤道附近的值本身为零。为进一步说明问题,图 2b给出了$r\left( {{\pi _{{\rm{b}}1}},{\pi _{{\rm{u1}}}}} \right)$的值,它在热带地区一直保持较大的负相关,赤道附近也是如此。上述情形说明,LBE方案中独立性假设在热带地区不成立。

图 2 模式13层的GRAPES短期预报误差样本中不同变量间的相关系数(r)(纬向平均值) Fig. 2 Zonally averaged correlation coefficients (r) between different variables of short-range forecast errors in GRAPES system at model level 13. $\psi $ and $\pi $ represent original stream function and Exner pressure, respectively; ${{\rm{\pi }}_{{\rm{b1}}}}$ and ${{\rm{\pi }}_{{\rm{u1}}}}$ represent balanced and unbalanced Exner pressure in LBE (linear balance equation) scheme, respectively; ${\pi _{{\rm{b2}}}}$ and ${\pi _{{\rm{u2}}}}$ represent balanced and unbalanced Exner pressure in hybrid scheme, respectively

进一步的,考察混合方案中$r\left( {\psi ,{\pi _{{\rm{u2}}}}} \right)$与$r\left( {{\pi _{{\rm{b}}2}},{\pi _{{\rm{u2}}}}} \right)$值。在中、高纬地区,上述值与LBE方案表现相当(北极附近的绝对值有所减小);而在热带地区以及副热带地区,它们的绝对值与LBE方案相比大幅减小,更加接近于零。其原因在于,在式(4)的多元回归模型中,最小二乘理论要求预报变量$\pi $的方差要等于回归估计值${\pi _{{\rm{b2}}}}$的方差与残差${\pi _{{\rm{u2}}}}$的方差之和,这也意味着${\pi _{{\rm{b2}}}}$与${\pi _{{\rm{u2}}}}$间的交叉协方差必须为零(黄嘉佑,2010)。注意到,图 2b的值并不完全等于零,这是由于本文计算中采用岭回归替代了经典回归,方差最小原则有所放松。但总体而言,混合方案的确能更好地保证独立性假设的成立,热带地区尤为如此。其他层次的情形与上述分析类似。

那么独立性假设成立好坏对B矩阵量值的影响如何,我们通过对比两个方案统计得到的预报误差方差予以说明。由于两个方案在中、高纬差异较小,这里主要关注热带地区。首先对比两个方案统计得到的$\pi $平衡部分的预报误差方差(图 3a),它们是根据同样的$\psi $样本经不同的平衡约束获取。如图,混合方案${\pi _{{\rm{b2}}}}$的方差要明显小于LBE方案${\pi _{{\rm{b1}}}}$的值,这可与图 1中系数R的值在热带地区接近零的情形相互印证。注意到,${\pi _{{\rm{b1}}}}$的方差在许多层次上已超过了$\pi $样本本身的方差值(黑色空心圈),而如果将${\pi _{{\rm{b1}}}}$与${\pi _{{\rm{u1}}}}$的方差的方差相加,得到LBE方案系统中实际使用的$\pi $的误差方差(见式(5)),其值更要远远大于NMC样本值。这种不合理情形出现的原因在于,${\pi _{{\rm{b1}}}}$和${\pi _{{\rm{u1}}}}$在热带具有过高的负相关(见图 2b),两者的交叉协方差(负值)已大到不可忽略。而混合方案由于引入统计模块,${\pi _{{\rm{b2}}}}$和${\pi _{{\rm{u2}}}}$几乎不相关,因而两者方差的和恰好等于NMC样本中$\pi $的方差。进一步的,图 3b给出了$\pi $平衡部分方差占同化系统中实际使用的$\pi $总方差的比例,混合方案也要明显小于LBE方案。

图 3 GRAPES热带预报误差样本中(a) $\pi $ (空心圈)和 $\pi $ 平衡部分(黑色线)的预报误差方差(单位:10−7),以及 $\pi $ 平衡和非平衡部分的方差和(灰色线)(单位:10−7),(b) $\pi $ 平衡部分方差占同化系统中 $\pi $ 总方差的比例。上述结果均为热带地区(20°S~20°N)的平均值,其中实线是LBE方案的结果,虚线是混合方案的结果 Fig. 3 (a) Averaged forecast error variances of different variables in GRAPES system in the tropics (20°S–20°N) (units: 10−7) and (b) ratio of variances (variance of balanced $\pi $ / variance of total $\pi $ used in the assimilation system)

综上,在热带地区,LBE方案由于存在虚假平衡问题,独立性假设不成立;混合方案通过引入统计模块保证了该假设的成立,系统中实际使用的$\pi $的总方差更真实反映了NMC样本应有值,平衡部分方差及其占总方差的比例均明显减小。那么这些变化对同化系统中的信息传递会产生什么样的影响,下一节将通过单点理想观测试验做进一步考察。

4 单点理想观测试验

单点理想观测试验是资料同化框架设计中考察、评估B矩阵结构的简单、有效方案。热带地区观测以质量场观测为主,风场直接观测很少,这里主要考察给定理想气压观测(均低于背景场1 hPa),LBE方案和混合方案强迫出的u风分析增量(分析场减去背景场)的差异,并将这里的结果与研究I的理论分析作对比。

图 4中理想观测均位于模式13层的180°(经度),而纬度分别设定为45.5°N、20.5°N和0.5°N。可以看到,在中、高纬45.5°N处,由于LBE本身适用性很好,混合方案加入统计模块后对原有的平衡约束影响很小,两个方案强迫出的u风几乎完全一致。而到了20.5°N,由于LBE适用性已不如中、高纬,混合方案对原有LBE做了一定程度削弱,强迫出的u风比LBE方案已有一定程度的减小。再进一步到赤道上,混合方案强迫出的u风已远远小于LBE方案,该区域上的风、压场分析接近独立。这样的结果与第2.2节和第3节对混合方案性能的分析是相对应的,热带地区原有LBE给定的虚假平衡约束被大幅削弱。两个方案中,由气压增量导出的v风间的差异与u风情形相一致。

图 4 理想气压观测(低于背景场1 hPa)位于模式13层的180°(经度)时强迫出的同一层上的u风分析增量(单位:m s−1):(a)、(c)、(e)LBE方案;(b)、(d)、(f)混合方案。纬度:(a、b)45.5°N;(c、d)20.5°N;(e、f)0.5°N Fig. 4 Analysis increments (units: m s−1) of zonal wind at model level 13 generated by a simulated pressure observation that is 1 hPa lower than the background at the same level at 180°, and at (a, b) 45.5°N, (c, d) 20.5°N, (e, f) 0.5°N. (a), (c), (e) The results of the LBE scheme; (b), (d), (f) the results of the hybrid scheme

图 4中不同纬度上u风分析增量的差异与Daley(1996)采用奇异向量方案削弱LBE虚假平衡的理论研究结果相一致(见其文中Fig. 7)。进一步的,图 4中单点位于赤道的情形与研究I在该区域的理论分析结果也具有较好的对应关系。LBE方案图 4e的结果与单纯考虑赤道罗斯贝波情景下风、压场协相关特征(研究I的图 4a)相近,气压(或位势高度h)与u风间存在显著负相关,且纬圈方向的相关尺度很大。而混合方案图 4f的结果虽未如研究I考虑所有赤道波动模态后的情形那样,hu风的相关转为很弱的正值(研究I图的5b),但协相关相对于单纯罗斯贝波情形大幅减小至接近于零的特征是相符的。

图 4e、4f一样,图 5中单点观测也位于0.5°N,但其垂直层次移到了第28层上,也即由对流层中层移至了平流层低层。在图 5中,混合方案强迫出的u风也是远小于LBE方案。然而,与图 4f不同的是,图 5b中强迫出了负的u风增量。上述情形说明,在赤道平流层低层,混合方案中气压与u风的协相关为正。这与图 1b中该区域负的R值相对应,也与研究I中70 hPa上hu正的协相关特征相一致(见研究I的图 6)。不过,与研究I的结果相比,这里给出的正相关偏小,结合图 4f未能表现弱的正相关,说明不区分具体波动的混合方案可能只部分反映了短期预报误差样本中Kelvin波的特征。

图 5 同图4e、4f,但垂直层次移到了模式28层 Fig. 5 As Figs. 4e and 4f, but for the results at model level 28

图 6 2013年5月,热带地区(20°S~20°N)u风分析相对于FNL再分析资料的(a)均方根误差(RMSE)与(b)偏差(Bias) Fig. 6 Analysis verification scores of u wind in the tropics (20°S–20°N) for the period of May 2013: (a) RMSE; (b) bias. Analysis results for each experiment were verified against FNL (final) data

综上,混合方案与LBE方案相比,质量观测信息向风场传递的基本特征在中、高纬变化很小,但热带地区的传递范围和量值均大幅度减小,风、压场分析接近独立。上述情形与研究I的理论分析结果具有较好的对应关系,其对GRAPES的同化预报的影响将在下一节中给出。

5 同化循环与预报试验 5.1 同化结果检验

为综合考察混合方案性能,特别是其对热带风场同化预报的影响,本文设计了接近实际业务的同化循环与预报试验。同化循环时段选取在2013年5月1日06时至31日18时,同化分析场经数字滤波后提供给模式作为预报初始场,模式的6 h预报场又提供给同化系统作为下一时次分析的背景场,如此往复循环。预报模式以及同化内、外循环分辨率均为1°×1°,垂直层次为36层。同化中使用的观测资料包括了目前GRAPES-VAR能业务使用的大部分常规和非常规观测资料,具体情形见表 1

表 1 同化循环试验中使用的观测资料种类 Table 1 Types of observations assimilated in the analysis cycle experiments

首先考察本文期望改善的热带地区的风场分析效果。图 6给出了试验时段内,热带地区(20°S~20°N)u风分析的均方根误差和偏差的垂直廓线,均相对FNL(final)再分析资料检验得到。从图中可以看出,LBE方案的u风分析在平流层低层(50 hPa~70 hPa)出现了异常大的分析误差和偏差,这与研究I给出的GRAPES短期误差样本中赤道罗斯贝波占比最小的区域刚好保持一致,LBE在该区域最不适用。而混合方案通过大幅削弱LBE在该区域的虚假平衡,显著改进了u风分析效果。与之对应,混合方案对该区域v风分析也有明显改进。

LBE方案之所以在热带平流层低层出现异常大的风场分析误差,一方面当然与其自身在该区域的虚假平衡有关。但另一个重要原因也必须指出,也即GRAPES系统自身的气压预报在该区域存有一定的正偏差。图 7给出了第28层(约为71 hPa)上经时间和纬圈平均的气压分析增量(∆p)和u风分析增量(∆u)。可以看到,该区域∆p为一致的负值,根据图 5a给出的LBE方案的信息传递机制,会在赤道附近强迫出很强的正的∆u,使得u风分析的正偏差不断增加。而混合方案中∆p虽然仍为负值,但根据图 5b显示的信息传递机制,其强迫出的∆u要小得多,且方向还发生了改变,变为负的增量,这有助于控制该区域风场的正偏差。

图 7 2013年5月,模式28层上的(a)气压分析增量(∆p)和(b)u风分析增量(∆u)的纬向、时间平均值 Fig. 7 Zonal and time averages of (a) pressure and (b) u wind analysis increments at model level 28 for the period of May 2013

热带地区平流低层的模式偏差是GRAPES发展过程中的现实问题,而同化系统中不合理的平衡约束LBE使得该问题在同化循环中被不断放大。混合方案通过给定合理的平衡约束减小了模式偏差的影响,使得整个GRAPES同化预报系统变得更为稳健、可靠。混合方案在显著改善热带地区风场分析效果的同时,也部分提高了该区域内高度场和温度场的分析效果,这可以从图 8中看出。

图 8 2013年5月,热带地区(20°S~20°N)位势高度分析(a)和温度分析(b)相对于FNL再分析资料的均方根误差(RMSE) Fig. 8 Analysis (a) geopotential height and (b) temperature RMSE in the tropics (20°S–20°N) for the period of May 2013. Analysis results for each experiment were verified against FNL (final) data

与LBE方案相比,混合方案对赤道外地区的同化分析呈中性偏正的效果,这可从下文针对预报效果的检验中一并看出。另外,值得补充说明的是,混合方案中的系数R是事先离线统计好的,同化系统直接使用,增加统计模块对同化分析的耗时影响很小,在多核并行条件下几乎可以忽略。

5.1 预报结果检验

利用上述同化循环得到的分析场,在每天12时(UTC)做72小时预报,进一步检验混合方案对GRAPES预报性能的影响。这里预报效果的检验也是相对于FNL再分析资料进行的,具体指标采用中国气象局数值预报中心的标准化预报检验工具GETv1.01计算得到。

首先检验赤道外地区(20°N~90°N,20°S~90°S)的情形,表 2给出了两个方案的72 h预报场在对流层中层(500 hPa)和高层(100 hPa)常用检验指标上的对比。为方便比较,表 2给出的是混合方案的相对提高率,对于距平相关(Anomaly Correlation,AC)而言,计算方法是混合方案的评分减去LBE方案的评分,结果除以LBE方案的评分;而对于均方根误差(RMSE)而言,计算方法是LBE方案的值减去混合方案的值,结果除以LBE方案的值。根据上述计算方式,表中正值表示引入混合方案后有正效果,负值表示有负效果。从表 2中可以看出,混合方案在赤道外区域与LBE方案差异很小,表现为略微偏正的效果,这与前文讨论的LBE本身在中、高纬适用性较好,引入统计模块后对其修正很小的结论相一致。

表 2 赤道外地区(20°N90°N20°S90°S)72 h预报场中,混合方案与LBE方案相比,距平相关(AC)的相对提高率以及RMSE的相对减小率 Table 2 Rate of AC(anomaly correlation)scores increase and RMSE decrease in hybrid scheme compared to LBE scheme for 72-h forecasts in the extratropics(20°N-90°N,20°S-90°S)

对于热带地区风场而言,经72 h预报之后,混合方案较LBE方案仍有明显改进,这可以从表 3看出。与图 6针对分析结果的检验相一致,混合方案对高层风场的改进大于低层,平流层低层的改进尤为明显。在出现分析异常的50 hPa上,混合方案中u风72 h预报场的均方根误差要比LBE方案小一半以上,v风相应值也大幅减小,LBE方案中虚假平衡造成的同化预报异常被基本消除。

表 3 赤道地区(20°S~20°N)72 h预报场中,与LBE方案相比,混合方案风场RMSE的相对减小率 Table 3 Rate of RMSE decrease in hybrid scheme compared to LBE scheme for 72-h forecasts in the tropics(20°S-20°N)
6 结论与讨论

针对线性平衡方程LBE作为风、压场平衡约束在热带应用时出现的虚假平衡问题,本文在GRAPES-VAR中引入了动力与统计混合平衡约束方案,以改善热带地区风场的分析效果。混合方案分为两步来实施:第一步利用LBE根据流函数$\psi $逐层计算初步的平衡气压${\pi _{b1}}$;第二步利用逐层、逐纬圈统计得到的系数R对${\pi _{b1}}$的垂直廓线做加权处理,以得到最终平衡气压${\pi _{b2}}$。统计样本采用NMC方法获得,样本长度达到一年,回归方案采用岭回归技术。通过样本统计与数值试验,主要得到以下结论:

(1)根据统计得到的系数R,混合方案求得的${\pi _{b2}}$主要由与其相邻几层的${\pi _{b1}}$加权求和得到,相距较远层次上的权重很小。而从R在不同纬度上的变化来看,热带地区的值与中、高纬相比大幅度减小,在对流层中层其数值接近于零,而在平流层低层转为负值。上述情形表明,在热带地区,混合方案第一步采用LBE引入的虚假平衡会在第二步的统计模块被大大削弱。

(2)变分同化中B矩阵构造的一个基本假设是:独立分析变量间的预报误差不相关。LBE方案在中、高纬能较好的保证假设的成立,但其在热带获得的独立分析变量高度相关,假设完全不成立。混合方案通过引入统计模块,减小了LBE在热带使用时的虚假平衡,保证了该区域内假设的成立。这样的差异使得,混合方案在热带地区实际使用的$\pi $的总方差更真实的反映了NMC样本应有情形,$\pi $平衡部分方差及其占总方差的比例都要明显减小。

(3)单点理想观测试验的结果表明,与LBE方案相比,混合方案对中、高纬风、压场协相关修正很小,但大幅减小了热带地区风、压场的约束程度,两者分析接近独立。这一结果与研究I中基于赤道波动模态的理论分析结果相匹配。尤其是在平流层低层,新方案能部分反映该区域风、压场协相关由Kelvin波主导的特征。

(4)2013年5月的同化循环与72 h预报试验结果表明,混合方案在赤道外地区呈中性偏正效果;而在热带地区,混合方案通过减小虚假平衡,避免了质量观测信息向风场分析的不合理传播,提高了同化系统的稳健程度,进而显著改善了该区域内风场的同化预报效果,平流层低层的效果最为显著,LBE方案在该区域的同化预报异常被基本消除。

在本文撰写期间,混合方案由于整体性能较优,已被GRAPES平行试验系统选为平衡约束的首选方案。不过本文研究主要针对现有N算子进行了优化,未涉及散度风与旋转风以及散度风和质量场间的平衡约束算子(分别记为ML)。关于M算子对GRAPES-VAR的影响,我们已做初步研究,它的作用主要集中于中、高纬边界层上,对热带地区的分析影响很小(王瑞春等,2012)。而对于L算子,Derber and Bouttier(1999)研究表明它对热带地区同化分析具有一定影响,但Chen et al.(2013)的工作表明,引入L算子对WRF在热带地区的预报效果影响很小。后续工作中将研究和评估L算子对GRAPES-VAR同化预报的影响。

最后,从本研究工作中可以看出,热带地区风、压场间的联系要远小于中、高纬。这也意味着,热带地区的风场分析难以从大量卫星质量观测中获取有用信息。为根本改善该区域风场分析就需要进一步挖掘云导风的潜力,并高度关注有关国际组织的星载多普勒雷达直接测风计划(Riishojgaard et al.,2012)。

致谢 感谢中国气象局数值预报中心各位专家老师给予本研究的指导和帮助,感谢两位匿名审稿专家和编辑老师对文章提出的宝贵意见。

参考文献
[1] Bannister R N. 2008. A review of forecast error covariance statistics in atmospheric variational data assimilation. Ⅱ:Modelling the forecast error covariance statistics[J]. Quart. J. Roy. Meteor. Soc., 134(637):1971-1996.
[2] Barker D M, Huang W, Guo Y R, et al. 2004. A three-dimensional variational data assimilation system for MM5:Implementation and initial results[J]. Mon. Wea. Rev., 132(4):897-914.
[3] Berre L. 2000. Estimation of synoptic and mesoscale forecast error covariances in a limited-area model[J]. Mon. Wea. Rev., 128(3):644-667.
[4] Chen Y D, Rizvi S R H, Huang X Y, et al. 2013. Balance characteristics of multivariate background error covariances and their impact on analyses and forecasts in tropical and Arctic regions[J]. Meteor. Atmos. Phys., 121(1-2):79-98.
[5] Daley R. 1991. Atmospheric Data Analysis[M]. Cambridge:Cambridge University Press, 457pp.
[6] Daley R. 1993. Atmospheric data assimilation on the equatorial beta plane[J]. Atmos.-Ocean, 31(4):421-450.
[7] Daley R. 1996. Generation of global multivariate error covariances by singular-value decomposition of the linear balance equation[J]. Mon. Wea. Rev., 124(11):2574-2587.
[8] Derber J, Bouttier F. 1999. A reformulation of the background error covariance in the ECMWF global data assimilation system[J]. Tellus, 51A(2):195-221.
[9] Fisher M. 2003. Background error covariance modeling[C].//ECMWF seminar on recent developments in data assimilation for atmosphere and ocean, 8-12 September 2003, Reading:ECMWF, 45-63.
[10] Holton J R. 1992. An Introduction to Dynamic Meteorology[M]. 3rd ed. New York:Academic Press, 507pp.
[11] 黄嘉佑. 2010. 气象统计分析与预报方法[M]. 第3版. 北京:气象出版社, 298pp. Huang Jiayou. 2010. Statistical Analyses and Prediction Methods in Atmospheric Science(in Chinese)[M]. 3rd ed. Beijing:China Meteorological Press, 298pp.
[12] Huang X, Xiao Q, Barker D M, et al. 2009. Four-dimensional variational data assimilation for WRF:Formulation and preliminary results[J]. Mon. Wea. Rev., 137(1):299-314.
[13] Ingleby N B. 2001. The statistical structure of forecast errors and its representation in the Met. office global 3-D variational data assimilation scheme[J]. Quart. J. Roy. Meteor. Soc., 127(571):209-231.
[14] Kleist D T, Parrish D F, Derber J C, et al. 2009a. Introduction of the GSI into the NCEP global data assimilation system[J]. Wea. Forecasting, 24(6):1691-1705.
[15] Kleist D T, Parrish D F, Derber J C, et al. 2009b. Improving incremental balance in the GSI 3DVAR analysis system[J]. Mon. Wea. Rew., 137(3):1046-1060.
[16] Körnich H, Källén E. 2008. Combining the mid-latitudinal and equatorial mass/wind balance relationships in global data assimilation[J]. Tellus, 60(2):261-272.
[17] Lorenc A C, Ballard S P, Bell R S, et al. 2000. The Met. Office global three-dimensional variational data assimilation scheme[J]. Quart. J. Roy. Meteor. Soc., 126(570):2991-3012.
[18] Parrish D F, Derber J C. 1992. The national meteorological center's spectral statistical interpolation analysis system[J]. Mon. Wea. Rev., 120(8):1747-1763.
[19] Parrish D F, Derber J C, Puser R J, et al. 1997. The NCEP global analysis system:Recent improvements and future plans[J]. J. Meteor. Soc. Japan, 75(1B):359-365.
[20] Polavarapu S, Ren S Z, Rochon Y, et al. 2005. Data assimilation with the Canadian middle atmosphere model[J]. Atmos.-Ocean, 43(1):77-100.
[21] Riishojgaard L P, Ma Z Z, Masutani M, et al. 2012. Observation system simulation experiments for a global wind observing sounder[J]. Geophys. Res. Lett., 39(17), doi:10.1029/2012GL051814.
[22] 万齐林, 薛纪善. 2007. 曲率修正线性平衡方程及其在变分同化风压约束中的应用[J]. 热带气象学报, 23(5):417-423. Wan Qilin, Xue Jishan. 2007. The curvature-modification linear balance equation and its application to the balance between wind and pressure in variational data assimilation[J]. Journal of Tropical Meteorology(in Chinese), 23(5):417-423.
[23] 王瑞春, 龚建东, 张林. 2012. GRAPES变分同化系统中动力平衡约束的统计求解[J]. 应用气象学报, 23(2):129-138. Wang Ruichun, Gong Jiandong, Zhang Lin. 2012. Statistical estimation of dynamic balance constraints in GRAPES variational data assimilation system[J]. Journal of Applied Meteorological Science(in Chinese), 23(2):129-138.
[24] 王瑞春, 龚建东, 张林, 等. 2014. 利用整层模式大气统计求解GRAPES-3DVAR动力平衡约束的数值试验[J]. 热带气象学报, 30(4):633-642. Wang Ruichun, Gong Jiandong, Zhang Lin, et al. 2014a. Numerical experiments on statistical estimation of dynamic balance constraints in GRAPES-3DVR with whole layers of model atmosphere[J]. Journal of Tropical Meteorology(in Chinese), 30(4):633-642.
[25] 王瑞春, 龚建东, 张林, 等. 2015. 热带风压场平衡特征及其对GRAPES系统中同化预报的影响研究 I:平衡特征分析[J]. 大气科学, 39(5):953-966. Wang Ruichun, Gong Jiandong, Zhang Lin, et al. 2015. Tropical balance characteristics between mass and wind fields and their impact on analyses and forecasts of GRAPES system. Part I:Balance characteristics[J]. Chinese Journal of Atmospheric Sciences(in Chinese), 39(5):953-966, doi:10.3878/j.issn.1006-9895.1412.14233.
[26] Wilks D S. 2006. Statistical Methods in the Atmospheric Sciences[M]. 2nd ed. Burlington:Academic Press, 627pp.
[27] Wu W, Purser R J, Parrish D F. 2002. Three-dimensional variational analysis with spatially inhomogeneous covariances[J]. Mon. Wea. Rew., 130(12):2905-2916.
[28] 薛纪善, 庄世宇, 朱国富, 等. 2008. GRAPES新一代全球/区域变分同化系统研究[J]. 科学通报, 53(20):2408-2417. Xue Jishan, Zhuang Shiyu, Zhu Guofu, et al. 2008. A study of GRAPES new generation global/regional variational data assimilation[J]. Chinese Science Bulletin(in Chinese), 53(20):2408-2417.
[29] 薛纪善, 刘艳, 张林, 等. 2012. GRAPES全球三维变分同化系统模式变量分析版[R]. 北京:中国气象局数值预报中心, 105pp. Xue Jishan, Liu Yan, Zhang Lin, et al. 2012. GRAPES-3DVAR Version-GM(in Chinese)[R]. Beijing:Numerical Weather Prediction Center of China Meteorological Administration, 105pp.
[30] Žagar N, Gustafsson N, Källén E. 2004. Variational data assimilation in the tropics:The impact of a background-error constraint[J]. Quart. J. Roy. Meteor. Soc., 130(596):103-125.
[31] 朱宗申, 胡铭. 2002. 一种区域格点三维变分分析方案——基本框架和初步试验[J]. 大气科学, 26(5):684-694. Zhu Zongshen, Hu Ming. 2002. A regional grid three-dimensional variational analysis scheme-Basic construction and preliminary tests[J]. Chinese Journal of Atmospheric Sciences(in Chinese), 26(5):684-694.
[32] 庄世宇, 薛纪善, 朱国富, 等. 2005. GRAPES全球三维变分同化系统——基本设计方案与理想实验[J]. 大气科学, 29(6):872-884. Zhuang Shiyu, Xue Jishan, Zhu Guofu, et al. 2005. GRAPES global 3D-Var system-Basic scheme design and single observation test[J]. Chinese Journal of Atmospheric Sciences(in Chinese), 29(6):872-884.
[33] 庄照荣, 薛纪善, 朱宗申, 等. 2006. 非线性平衡方案在三维变分同化系统中的应用[J]. 气象学报, 64(2):137-148. Zhuang Zhaorong, Xue Jishan, Zhu Zongshen, et al. 2006. Application of nonlinear balance scheme in three-dimensional variational data assimilation[J]. Acta Meteorologica Sinica(in Chinese), 62(2):137-148.