引言

延伸期天气预报是近20多年来的国际大气科学的研究热点，很多气象专家从不同的角度对延伸期预报方法进行了研究(Waliser et al, 2003)。目前，预报员利用数值天气中期预报模式输出结果，再结合经验修正，使天气预报的可预报时效达10 d左右。而针对10 d以上的延伸期天气预报，模式输出结果的参考意义十分有限。正因为如此，延伸期过程预报方法的研发成为目前气象学界的热点也是难点问题。延伸期预报困难的原因在于其预报时效超越了确定性预报的理论上限(两周左右)，而预报对象的时间尺度又小于气候预测的月、季时间尺度。随着卫星遥感等多种观测手段进步, 近10年来气候观测获得的各种数据迅速增加, 这些科学大数据反映和表征着复杂的自然现象与关系。从全球海量数据中提取出部分有效数据, 可以获得比过去更全面的大尺度低频变化信息表征区域极端天气变化信号, 为建立预报时效较长的简化的延伸期预测模型提供了很好的基础。所以，进一步采用不同的途径和方法，进行预报研究是非常重要的。所以, 采用不同的途径和方法, 进行预报研究是非常重要的。

总体而言，延伸期天气预报的基本方法主要有动力模式(丑纪范等，2010；Vitart et al, 2010; 陈官军等2012；郑志海等，2013；Fu et al, 2013)和统计方法两类(von Storch et al, 1990; 杨秋明, 1998; Waliser et al, 1999; Mo, 2001; Wheeler et al, 2004; Jiang et al, 2008; 琚建华等, 2010; 贾小龙等, 2012; 谷德军等，2013; 牛法宝等, 2013); 另外, 还有低频天气图(孙国武等, 2013; 陈青等，2014；信飞等，2014)、神经网络模型(Love et al, 2009)、灰参数微分方程(杨秋明, 1990)、大气物理量分解(钱维宏, 2012; 钱维宏等, 2013a；2013b)以及各种方法的综合集成等方法也取得了一定的预报效果。在这些方法中, 其主要途径之一是基于大气季节内振荡(interseasonal oscillation, ISO), 利用其动力机制分析可传播信号的放大和对预报区域的影响，建立合适的统计预报方程, 可以作出远高于随机判断和单纯依靠统计计算的有效的10~30 d延伸期天气预报。

目前统计模式对ISO的预报能力要比动力模式好, 在动力模式对ISO的模拟和预测进一步得到改进之前, 它是一种研究ISO可预报性的更加有效的途径。1990年以来, 基于ISO的统计方法及其对其可预报性的研究和预报试验逐渐增多, 预报时效也逐步延长。近5年来, 这些试验研究迅速增加, 部分方法已用于实时预测。这些方法主要有时滞回归分析、非整数波功率谱分析、奇异谱分析(SSA)、主振荡型分析(POP)、自回归模型(AR)等时间序列分析等方法, 一般平均预报时效可以达到20 d左右。由于ISO主要由大气内部动力过程决定, 它的强度存在显著的年际变化, 不同年份的ISO型变化特性及其与影响区域极端天气的联系和可预报性也有明显差异。所以, 需要对不同年份的ISO演变过程更详细的描述和分析, 才能提高预报稳定性和进一步延长预报时效。因此, 对引起区域极端天气(如强降水过程等)的ISO进行更多的预测试验, 从不同角度改进统计方程的建模方案对于提高预报精度是十分重要的。

大气各种ISO型对中国不同区域的降水过程形成和变化有不同的影响，过去国内已取得了许多有意义的成果(冯俊阳等，2012；张耀华等，2012；牛法宝等, 2013；刘炜等，2014)，这些研究对于降水10~30 d延伸期预报方法的改进有重要作用。对于东亚夏季风区降水延伸期预报，过去的预测研究主要基于东亚和南亚热带和副热带地区以及欧亚中高纬度地区ISO的变化(梁萍等, 2012; 甘晶晶等, 2013)，南半球环流ISO的影响考虑很少，不能反映南北半球环流之间的相互作用对东亚季风降水ISO的影响。因此，在延伸期预测研究中，考虑南半球低频环流的显著的作用，将是改进预报效果的可能途径之一。本文将用南半球热带外地区环流低频主成分, 对2013年初夏长江下游地区降水20~30 d低频分量(与强降水过程密切相关)(杨秋明, 2009), 由动态数据来驱动复杂低频变化过程与系统的构建，建立多变量时滞回归(multivariable lagged regression, MLR)预测模型, 进行独立预测试验, 延长预报时效; 同时, 也对2001—2012年资料进行历史回报试验研究, 讨论了可预报性的年际变化及其影响因子。

1 资料

850 hPa经向风场选用NCEP/NCAR逐日再分析风场资料(Kalnay et al, 1996)全球范围内的2.5°×2.5°的格点, 时间是2001—2013年3—8月; 采用一阶Butterworth滤波器进行20~30 d带通滤波, 得到全球850 hPa低频经向风场资料阵。逐日长江下游地区降水量采用区域30.5°~32.0°N、118.0°~122.5°E中5个代表站平均值(MICAPS逐日站点降水资料), 时间是3月15日至8月31日。长江下游地区低频降水预测试验的时间是6月1日至7月30日。

2 方法

本文采用基于主成分的MLR模型:

设某区域D中的M个观测样本和格点为N的滤波资料阵_MS_N=(s_ij), i=1, 2, …, M; j=1, 2, …, N, 主成分分析(principal component analysis, PCA)后[_MS_N=_MT_LLV_N^T, T=(t_j1(i))是时间系数矩阵]的前L个主要分量是t_j1(i), i=1, 2, …, M; j₁=1, 2, …, L; 长江下游地区低频降水(r_lcj)与时间滞后为τ的主要分量t_j1满足线性方程:

$ {{\mathit{r}}_{\mathit{lcj}}}\left(\mathit{i} \right)={{a}_{0}}\left(\tau \right)+\sum\limits_{{{j}_{1}}=1}^{\mathit{L}}{{{a}_{{{j}_{1}}}}\left(\tau \right)}{{\mathit{t}}_{{{j}_{1}}}}\left(i-\tau \right) $

(1)

系数a_j1, j₁=0, 1, 2, …, L由线性最小二乘估计, 当τ=1, 2, …, 30 d时, 长江下游地区低频降水(r_lcj)的1, 2, …, 30 d变化由上式预测。时滞线性回归式(1) 称为主成分低频多变量滞后线性回归模式, 反映了各个时间滞后τ主要分量t_j1和低频降水之间相互独立的回归关系; 式(1) 可变为

$ {{\mathit{r}}_{\mathit{lcj}}}\left(\mathit{i+}\tau \right)={{a}_{0}}\left(\tau \right)+\sum\limits_{{{j}_{1}}=1}^{\mathit{L}}{{{a}_{{{j}_{1}}}}\left(\tau \right)}{{\mathit{t}}_{{{j}_{1}}}}\left(i \right) $

(2)

对于初始时间的i=i₀的t_j1(i₀), 由时滞线性回归式(2) 可得到低频降水预测值r_lcy(i₀+τ), τ=1, 2, …, 30 d^*, 其中时滞相关结构[回归系数a_j(τ)随τ变化]体现了低频降水与各个环流低频分量之间的相互作用; 这种相互作用增强, 低频降水也增强, 其正位相或正负位相转换时产生强降水的概率明显增大。

* http://www.lcjrerf30.org/soft.asp

3 预测试验

本文采用3月15日至7月31日(样本长度N=139) 的长江下游地区逐日低频降水和850 hPa低频经向风场资料进行预测试验。由于PCA得到主要空间分布型与样本序列长度有关, 因此本文用限定记忆法, 保持子序列N₀不变, 用MLR预报模型滑动进行独立样本预测试验。具体方法如下:对于初始时间i₀, 使用i₀-N₀+1, i₀-N₀+2, …, i₀-1, i₀时次的子序列资料(N₀是子序列长度), 对850 hPa低频经向风场(经过带通滤波)资料(并标准化处理), 经过PCA得到主要空间分布型和逐日变化的主成分t_j(i), 用长江下游逐日低频降水r_lcj和低频经向风主成分t_j(i), 建立MLR; 由初始时间的i=i₀的t_j(i₀), 用30个MLR(分别对应于预测时间τ=1, 2, …, 30 d)预测未来30 d长江下游逐日低频降水${{\hat{r}}_{lcj}}\left( {{i}_{\rm{0}}}\rm{+1} \right)$, $ {{\hat{r}}_{lcj}}\left( {{i}_{\rm{0}}}\rm{+2} \right)$, … $ {{\hat{r}}_{lcj}}\left( {{i}_{\rm{0}}}\rm{+30} \right)$。本文计算中, 主分量个数L=4, 子序列长度N₀=78。从5月31日(i₀=78, 基于3月15日至5月31日的资料建立MLR，初始时间是5月31日，基于3月15日至5月31日期间资料的空间模态和变量之间的时滞联系，用建立的30个MLR预测6月1—30日长江下游降水20~30 d低频分量的变化)开始每隔5 d(固定子序列长度N₀=78) 滑动进行独立样本预测, 其中i₀=78, 83, …, 103, 108共作7次预测。当i₀=108时, 初始时间是6月30日，基于4月14日至6月30日期间资料的空间模态和变量之间的时滞联系(基本反映了初夏5—6月期间振荡模态特征)，预测7月1—30日长江下游降水20~30 d低频分量的变化。以上建模方案有利于适应方差贡献较大的主要低频振荡系统空间结构随时间的改变(从春季到夏季), 在一定程度上反映低频振荡时空结构变化的非线性特征。

图 1a给出2013年5—8月长江下游地区逐日降水的时间演变, 表现为5—7月上旬的3次集中强降水过程(5月17日、6月7—8日和7月5—6日)，降水过程中最大日平均雨量(区域)大于50 mm。7月中旬到8月下旬, 强降水过程减弱(7月20—21日、8月19—20和24—26日的强降水过程的最大日平均雨量在20~30 mm)。可见5—7月降水活跃期与中断期交替明显, 这种变化表现为极显著的季节内振荡, 进入8月以后, 这种振荡逐步减弱。图 1b是逐日降水时间序列的非整数波功率谱分析(Schickedanz et al, 1977)。通过α=0.10显著性水平检验的显著周期(非整数)是6~7 d, 10 d的高频振荡和较显著的13和25 d左右的低频振荡, 其中周期小于10 d的高频振荡和13 d左右的准双周振荡变化与强降水频数年际变化相关较小(主要与中等降水过程变化有关), 而25 d周期是与强降水频数密切相关的降水ISO, 这种20~30 d ISO的强度与长江下游夏季强降水频数的正相关最显著(杨秋明, 2009)。因此，2013年5—8月长江下游地区降水低频振荡中存在较显著的20~30 d振荡(尤其是5—7月), 有利于长江下游地区暴雨过程频繁发生, 夏季出现显著区域洪涝。下面对2013年初夏(6—7月)长江下游地区降水20~30 d低频分量进行延伸期预报试验。

过去的研究表明, ISO对长江下游地区强降水过程的作用具有时间尺度上的选择性，仅20~30 dISO强度与长江下游强降水(暴雨)存在最显著的稳定的正相关(杨秋明, 2009)，而且与全球环流20~30 d振荡型变化密切相关。这种降水20~30 d ISO活动和长江下游强降水频数均呈准2年周期变化(杨秋明, 2013), 与东亚梅雨强弱的年际变化关系不大(杨秋明等，2012)。全球大气20~30 d振荡是直接影响两半球热带外地区天气气候最主要的因子, 并且在20~30 d时间尺度上南北半球中纬度ISO之间存在显著的低频遥相关，其中南半球热带外绕球遥相关(the southern circumglobal teleconnection wave train, SCGT)波列的传播主要影响南半球中纬度地区降水变化, 这种20~30 d振荡也表现在半球尺度的涡动动能变化上, 即斜压环状模(baroclinic annular mode, BAM), 与斜压性和涡动热通量之间的双向反馈有关(Thompson et al, 2014a; 2014b)。通过大气遥相关, 可以间接影响长江下游地区20~30 d降水低频分量变化和强降水过程的形成, 是长江下游地区强降水延伸期预报的可预报性重要来源之一；而ENSO、QBO、MJO仅直接影响全球热带地区天气气候, 对热带外地区天气气候的影响具有一定的时空不稳定性和不确定性, 它们对长江下游地区强降水作用的直接影响也较小。上述20~30 d时间尺度的SCGT的传播是全球热带外地区主要气候振荡模态之一(主要由大气内部动力过程决定)，它与南非、澳大利亚、南美、印度地区和长江下游地区降水低频变化关系密切。通过南半球热带外类似SCGT的低频波列与热带印度洋地区环流的相互作用, 调控越赤道气流来影响东亚季风的强度和长江下游地区20~30 d低频降水变化, 是影响东亚副热带地区强降水过程的主要ISO型之一(杨秋明, 2009)。所以本文采用南半球中高纬度地区(35°~65°S、0°~360°)方差贡献较大的前4个低频主成分PC1~PC4(图 2给出了前4个主要空间模态(资料时间是2013年3月15日至5月31日), 解释方差分别是46.2%、39.1%、6.8%和5.8%, 均表现为沿纬向传播的波列结构, 其中第一、第二模态对应于SCGT的传播, 而第三、第四模态对应于另一类沿纬向传播的波列结构, 空间尺度较小, 仅在南大西洋和南印度洋地区传播, 反映了类似于SCGT的波列活动的多样性)作为因子构建预测模型。每次独立预测试验中, 将南半球中高纬度地区的低频经向风计算出的前4个主要空间模态对应的低频主成分PC1~PC4的观测值(子序列长度为78 d), 建立MLR模型作2013年6—7月长江下游低频降水变化的30 d预测。

图 3a~3g给出了7次预报结果(初始时间是5月31日, 6月5、10、15、20、25日和30日, 图中r是预测和实况之间的相关系数, 即预报技巧)。从图 3a~3g中可以发现这7次预报预报效果均很好, 预报技巧分别是0.91、0.93、0.86、0.90、0.80、0.86和0.87, 均大于0.60(考虑序列的持续性, 通过α=0.05的显著性水平检验), 最好的是初始时间是6月5日的30 d预报, 预报技巧达到0.93(图 3b, 表明这一时段长江下游地区20~30 d降水ISO与SCGT之间存在极显著的相关), 它们能较好地预测长江下游地区6月上旬的多雨过程(其中6月7—8日长江中下游的大暴雨过程(图 4a和4b))、中旬的干旱少雨和6月下旬多雨洪涝过程。另外, 初始时间是6月10、15、20和25日(图 3c~3f)的4次30 d预报, 均准确预报7月5日左右降水低频位相的正负转换(它反映了副热带低频系统迅速增强, 形成7月5—6日长江中下游地区的大暴雨过程，如图 4c和4d), 而初始时间是6月25日(图 3f)的30 d预报，则预测7月上旬后期到7月中旬进入高温少雨时段(低频分量的负位相)，7月下旬前期再次进入多雨时段(低频分量的正位相)。因此, MLR模型可以在6月10日提前25 d左右预报出7月5日的长江下游地区大暴雨过程。这表明南半球SCGT波列对于2013年5—7月长江下游地区强降水形成起主导作用(通过南半球中纬度地区环流和东亚环流在20~30 d时间尺度上的强相互作用), 是2013年初夏长江下游地区暴雨10~30 d延伸期预报的重要可预报性来源之一。另外，初始时间是6月20日的30 d预报(图 3e)，相关技巧有所降低(0.80)，主要是预测的16~30 d降水低频分量负位相超前2~4 d，表明这个时段建立的MLR模型15 d以后的预报误差增大，可能与低频降水和低频环流之间相互作用的非线性不稳定性增大有关，其机制有待深入研究。

上述预测试验表明, 基于2013年夏季长江下游地区降水20~30 d ISO与南半球地区ISO之间的相互关系, 选择合适关键区的低频主分量作为因子, 针对一种位相上的超前滞后关系, 构建MLR预报模型, 能有效预测10~30 d延伸期长江下游地区降水低频分量的变化, 较好地体现了低频降水与主要环流低频分量之间较强的相互作用, 为长江下游地区延伸期强降水过程发生的预报提供主要依据。

进一步利用多年(2001—2012年)初夏资料, 分别建立MLR对每年6—7月的长江下游地区降水20~30 d低频分量预测试验(35°~65°S、0°~360°的低频PC1~PC4作为4个预报因子)结果表明, 预报时效也达到30 d左右(表 1), 能很好地预测多强降水过程的正位相的变化。从表 1可以看出, 除了2005和2012年初夏的7次预报的平均相关技巧较小外[分别是0.45和0.22, 未通过0.05的显著性水平检验(考虑序列持续性的影响), 反映从春到夏, SCGT对长江下游低频降水的作用较弱或不稳定], 其余年份均能有效预测长江下游地区降水20~30 d低频分量的未来30 d变化趋势(尤其是2001—2004、2007—2009、2011和2013年, 平均相关预报技巧≥0.80, 达到0.02的显著性水平检验)。13年(2001—2013年)初夏的91次预测试验的平均预报技巧是0.78, 表明了这些20~30 d振荡偏强和正常年份中, 南半球中纬度地区SCGT与夏季长江下游地区20~30 d低频降水存在强相互作用, 是影响长江下游地区强降水形成和变化的主导因素。图 5a~5l是2001—2012年期间每年初始时间为6月30日的12次预报, 其中预报技巧最大的是2001和2009年, 均达到0.97(图 5a和5i), 最低的是2005和2012年(图 5e和5l)。13年中有11年的6月30日均能较准确预报7月1—30日长江下游地区低频降水分量变化, 特别是成功预测长江下游地区空梅的2002年, 2009年7月中旬后期到下旬20~30 d低频降水分量正位相的变化(图 5b和5i，对应于7月中下旬长江下游地区持续多雨, 强降水过程频繁发生)。另一方面, 非整数波功率谱分析也表明, 2012年仅存在40~50 d振荡(图略), 20~30 d振荡不存在(强降水过程明显偏少, 主要是中等降水过程), 所以20~30 d振荡分量很小, MLR模型容易出现显著的不稳定, 预报时效很短(图 5l); 而2005年5—6月20~30 d振荡较弱, 7—8月迅速增强(即20~30 d振荡强度存在季节内的突变, 图略), 因此20~30 d低频分量的MLR预测模型也出现一些不稳定, 预测误差明显增大(图 5e)。

通过以上多年预测试验表明，基于初夏长江下游地区低频降水和南半球热带外环流各个主要分量之间较高的超前(滞后)相关关系, 构建MLR预测模型(反映类似SCGT波列的变化的多样性), 可以很好地作出20~30 d振荡活动较明显和正常的大部分年份初夏长江下游低频降水分量30 d左右的延伸期预报。南半球热带外类似SCGT的波列强度和位置, 也有一定的年际变化(波列在35°~65°S之间东移), 本文采用滑动序列的建模方案中主要空间模态也随时间有所变化, 以适应不同年份这些波列强度和东移传播途径的一些变化。此外, 对于实时ISO预报, 由于用传统的滤波方法提取ISO特征的过程中会在序列的前后各损失10 d左右的资料, 对于初始时间, 如需采用ECMWF环流的10 d预报(Miller et al, 2010)延长资料序列(右端延长10 d), 滤波后再略去序列右端的10 d的资料, 得到右端精度较高的滤波序列, 建立MLR模型。用得到较精确的初始时间t₀的滤波值, 预报未来30 d低频分量变化, 使得预报时效稳定在25 d以上, 可以有效提高实时预报精度。值得注意的是，在延伸期预报中，本文中现时的带通滤波窗口取为20~30 d，当这种ISO分量显著增大或稳定时，振荡过程的正位相或正(负)位相转换时，伴随明显的长江下游地区强降水过程；而这种分量减弱时，强降水过程明显减少，只对应于一般降水过程(小到中雨或中到大雨过程)。若20~30 d ISO分量很弱, 这时10~20 d或30~50 d降水ISO可能起主导作用(如2012年)，可用选择对应的关键区的低频主分量, 建立MLR模型预测它们未来30 d的变化，估计多(少)雨时段; 但这两类ISO变化过程中，长江下游地区出现强降水的概率较小。因此，通过预测20~30 d降水ISO分量大小以及符号的变化, 可以更客观地估计未来30 d中长江下游地区强降水过程发生的敏感时段，提高极端天气事件预测的精度。

4 结论

本文分析了2013年长江下游地区5—8月降水的ISO变化特征, 并用多变量时滞回归模型(MLR), 对影响2013年初夏长江下游地区强降水过程的20~30 d降水低频分量进行10~30 d延伸期逐日变化预测试验, 同时对2001—2012年初夏的资料进行历史回报试验研究, 分析可预报性的年际变化, 得到如下结论:

(1) 2013年夏季长江下游地区降水存在较显著的25 d左右和13 d左右的振荡周期, 前者(20~30 d ISO)伴随着5—7月3次集中强降水过程(5月17日、6月7—8日和7月5—6日), 这种显著的20~30 d振荡位相变化是造成强降水(暴雨或大暴雨)过程的主要原因; 而后者(10~15 d ISO)的变化, 有利于中等强度的降水过程的形成。

(2) 基于与南半球热带外地区SCGT波列变化相关的850 hPa经向风低频主要模态, 建立MLR预报模型可以有效预测2013年6—7月长江下游降水20~30 d低频分量未来30 d的变化, 能为长江下游延伸期强降水过程发生的时段预测提供预报信号, 将长江下游强降水事件的预报时效延长到25~30 d。进一步利用多年(2001—2012年)6—7月资料的历史回报试验表明, 对于20~30 d振荡偏强和正常的大部分年份, 预报时效也稳定在25 d以上；仅20~30 d振荡不明显的2012年, 预报技巧明显偏低, 预报时效较短。另外, 2005年20~30 d振荡强度季节内突变(5—6月较弱, 7—8月显著增强), 预测模型出现一些不稳定, 预报技巧下降。在2001—2013年期间初夏的91次预测试验的平均预报技巧是0.78, 能较好预报6—7月长江下游地区多(少)雨时段的变化, 为强降水过程延伸期预测提供主要信息。

(3) 通过预测20~30 d降水ISO分量大小以及符号的变化, 可以估计长江下游强降水发生的时段。而这种分量很弱时(如2012年)，强降水过程明显偏少，正位相或正负位相转换时一般形成中等降水过程(或小到中雨过程)为主。这时10~20 d和30~50 d降水ISO也可能起主导作用，可用MLR模型预测它们未来30 d的变化，估计多(少)雨时段。

(4) 基于大量观测资料中隐藏的主要时滞相关关系的变化, 完全由动态数据驱动构建MLR预测模型, 可以有效地全面描述多样化的两半球中纬度大气20~30 d振荡的变化规律, 较精确地预测10 d以上的长江下游地区降水低频变化过程; 而且选择适当的子序列长度, 滑动进行动态建模的方法可将环流的复杂系统空间结构变化及其与低频降水之间时滞相关的变化信息反馈到预测模型中, 有利于提高预报精度。因此，以数据为基础的延伸期预报技术的发展有助于未来10 d以上天气预报精度的明显提高，大量观测资料中提取的关键数据的更新(动态数据)是构建时变简化的低频分量预测模型，并获得更长时效延伸期预报的最主要途径。