CN101852871A - 基于经验模态分解和数值集合预报的短期气候预测方法 - Google Patents

Abstract

本发明公布了一种基于经验模态分解和数值集合预报的短期气候预测方法，本发明使用数值集合预报技术与均生函数逐步回归模型相整合的方式，并结合了经验模态分解(EMD)这种处理数据序列的新方法，首先将非平稳的气候数据序列分解成平稳的、多尺度特征的本征模态函数分量(IMF)，再对每一个IMF使用集合预报与逐步回归分析相结合的方式构建不同的预报模型，最后线性拟合成预报结果。使用本发明的***进行短期气候预测时，用户可以根据实际数据需要截取指定的序列长度和预报长度，并且针对性的选取集合预报过程中的预报模型参数。本发明相对于直接预测和单一预测方法，对于气候的变化趋势以及突发性气候具有更好的预报能力。

Description

基于经验模态分解和数值集合预报的短期气候预测方法

技术领域

本发明涉及一种基于经验模态分解和集合预报技术对气候序列数据(20-50年)的短期预测(5-10年)，适用于对气候突发性变化需要敏感把握的情况，如霜冻(降温)、暴雨等灾害性天气。

背景技术

气候预测是大气科学的一个重要分支，而气候变化的预测则又是各级气象台站开展业务短期气候预测的最主要内容之一。

气候***是一种耗散的、具有多个不稳定源的高阶非线性***，其复杂的内部相互作用和自由变化导致了气候的可变性和复杂性，也导致了气候变化的多尺度特性^[1-3]。由于目前我们对岩石圈、水圈、冰雪圈、生物圈(包括人类活动)的了解甚少，所积累的资料精度和完备性十分不足，因此就目前人类对气候***的了解而言，要建立起精确的反映气候***变化的非线性动力学方程组是不可能的^[4]，这就给气候变化预测的动力学理论研究带来了极大的困难。

就目前来说，我们赖以进行短期气候预测工作的根据有两个方面：一是对气候变化进行了严密的监视，利用各种仪器在全球范围进行了定时观测，从而积累了关于它过去和现在情况的大量数据；二是气候***为一物理***，通过物理学和数学的研究，获知了一些气候变化应遵循的物理规律，并可以用数学语言来表示^[5，6]。第二种方法实际上就是一种时间序列预测方法。

时间序列是指同一种现象在不同时间上的相继观察值排列而成的一组数字序列。时间序列预测方法的基本思想是：在预测一个现象的未来变化时，用该现象的过去行为来预测未来。即通过时间序列的历史数据揭示现象随时间变化的规律，将这种规律延伸到未来，从而对该现象的未来做出预测。随着计算机技术的发展，特别是自从著名统计学家Box和Jenkins于1968年发展了一整套随机时间序列的模型识别、参数估计和诊断检验的经典建模^[7]方法以来，时间序列预测的发展非常迅速，在气象、天文、电力、医学、生物、经济、金融和计算机等各个领域已有广泛的应用，并日益显示出强大的生命力和重要性，已成为一门独立的重要数学分支。迄今已有不少专家、学者从事这方面的研究，并建立了一套比较完整的理论和应用体系^[8，9]。

本发明是一种针对气候的预测方法，属于时间序列预测的一种，而考虑到气候***的特殊性，在预测过程中加入了两种对数据序列的处理，即经验模态分解(Empirical ModeDecomposition，EMD)和集合预报技术。EMD方法从本质上讲是对一个信号(或其导数，视所需的分解精度而定)进行平稳化处理，其结果是将信号中不同尺度的波动或趋势逐级分解开来，产生一系列具有不同特征尺度的数据序列，每一个序列称为一个本征模函数(Intrinsic Mod Function，IMF)分量，不同的分量成为了平稳信号，使得后面的数据预测等处理更加合理、有效。集合预报是针对单一模式预报不确定性而引入到数值预报中的，它过滤了预报中的随机成分，从而能获得一个相对稳定的预报结论，本发明中将序列值根据预先指定的成员数进行处理，重新组合再求集合平均，以期得到不稳定解的概率降低。

参考文献

[1]林振山，杨修群.理论气候学[M].南京：南京大学出版社，1996.

[2]林振山.气候建模、诊断和预测的研究[M].北京：气象出版社，1996.

[3]尤卫红.气候变化的多尺度诊断分析和预测的多种技术方法研究[M].北京：气象出版社，1998.

[4]林振山，邓自旺，尤卫红.理论气候学中的几个基本问题[J].热带气象学报，1995，11(2)：187-192.

[5]胡增臻，黄荣辉.长期天气预报业务和方法研究的最新进展[J].气象科技，1993，(1)：1-10.

[6]尤卫红.对云南省气象台近年来长期天气预报的检验[J].云南气象，1992，(2)：30-31.

[7]Pawlak Z.Rough Sets[J].International Journal of Computer and information Science，1982，11(5)：341-356.

[8]Hu Xiao Hua，Cercone N.Learning in relational databases：a rough set approaeh[J].Computational Intelligence，1995，11(2)：323-337.

[9]王国胤，于洪，杨大春.基于条件信息熵的决策表约简[J].计算机学报，2002，25(7)：759

发明内容

本发明目的是针对气候数据序列的不稳定和非线性特征提供一种复合方式的气候预测方法。将经验模态分解方法与均生函数逐步回归模型以及集合预报技术相结合，进行气候预测。该发明对于气候的变化趋势以及突发性气候具有更好的预报能力。

本发明为实现上述目的，采用如下技术方案：

本发明基于经验模态分解和数值集合预报的短期气候预测方法，包括以下几个步骤：

(1)对数据进行z-score标准化处理：将原始输入序列A的原始值v使用z-score标准化到标值v′，即

$v^{\′}=\frac{v-\stackrel{\‾}{A}}{{\mathrm{\σ}}_A},$

其中

和σ_A分别为原始输入序列的均值和标准差；

(2)将步骤(1)所述的标值v′用EMD算法分解，得出n个IMF分量和一个趋势分量，其中n为大于1的自然数；

(3)将步骤(2)所述分解出的n个IMF分量分别用数值集合预报方法进行处理，其中行数即为成员数s；随机指定每一行去掉的元素，预报元素逐行减少，最终构成一个拥有s行的变长序列组；

(4)最后对步骤(3)所述的由IMF分量处理所得的序列组构建均生函数模型，预测计算并记录结果，然后将每组的结果拟合成一个IMF分量的预报结果，最后再将n个包含不同特征的预报结果线性拟合成最终预报解。

优选地，步骤(2)所述的经验模态分解算法EMD，在分解的筛选过程中，采用设置两个连续的迭代结果之间的标准差SD的大小在0.2到0.3之间，且最大迭代次数为200作为追加的停止条件，采用极值延拓法来拟合序列的端点极值。

优选地，步骤(2)所述的n个IMF分量和一个趋势分量中，第一个IMF分量对应随机成分，第一至第n个IMF分量对应周期成分，趋势分量即残差Rn对应趋势成分。

优选地，步骤(3)中所述每一行中的元素，零值不参与预报。

本发明基于经验模态分解和集合预报技术的短期气候预测方法，根据气候时间序列的非线性、非平稳特性，首先利用经验模态分解技术对时间序列进行平稳化处理，在保留时间序列本身特征的基础上减少序列间的干扰或耦合信息，再利用集合预报技术与均生函数逐步回归模型相结合的方式进行时间序列预测。相对于单一使用均生函数逐步回归的方法，可以有效提高预测准确度，特别适合处理逐年降水或温度变化的非平稳气候时间序列。

附图说明

图1是基于经验模态分解和数值集合预报技术的气候时序预测方法的方法流程图；

图2是对气候时间序列进行经验模态分解的具体处理流程图；

图3是使用集合预报和均生函数逐步回归模型相结合的预测方法流程图。

具体实施方式

下面结合附图对发明的技术方案进行详细说明：

本发明基于经验模态分解这种全新的处理序列数据的方法，首先将非平稳非线性的数据序列分解成若干个代表一组特征尺度的新数据序列，先将原始数据序列分解为各种不同特征波形的叠加。其中最关键的一步是通过信号的极值点拟合信号包络线，本发明采用的是应用最广泛的三次样条插值函数法。三次样条函数需要信号两端数据的一阶或二阶导数作为其边界已知条件，而由EMD算法的原理可知，无法直接获得两端点对应的极值，针对本发明在短期气候预测方面的应用，采用极值延拓法作为端点极值拟合方法。本发明具体包括以下步骤：

(1)对数据进行z-score标准化处理，这种方法基于原始数据的均值(mean)和标准差进行数据的标准化。将A的原始值v使用z-score标准化到v′，即

$v^{\′}=\frac{v-\stackrel{\‾}{A}}{{\mathrm{\σ}}_A}$

其中v为原值，

和σ_A分别为序列的均值和标准差。

(2)将已预处理的数据用EMD算法分解，得出n个IMF分量。时间序列一般包含随机成分、周期成分和趋势成分。一般来说时间序列中随机成分对应一些突发事件或噪声，周期成分对应气温变化周期，而趋势成分则对应长期大的变化规律。经过EMD分解，原序列被分解成若干个IMF分量和一个趋势分量，IMF1(第一个分量)对应随机成分，IMF2、IMF3、……IMFn等对应各周期成分，残差Rn对应趋势成分。并且IMF1的振幅最大，波长最短，平稳性最差，随着分解次数的增加，分量的非平稳行为逐渐减小，并且不同分量代表序列在不同时间周期上的波动规律，更能准确反映原序列的变化，并保留原序列的本身特征；

(3)将步骤(2)所述分解出的若干IMF分量分别用数值集合预报方法进行处理，其中行数即为成员数s。随机指定每一行去掉的元素(零值不参与预报)，预报元素逐行减少，最终构成一个拥有s行的变长序列组。

(4)最后对步骤(3)所述的由IMF分量处理所得的序列组构建均生函数模型，预测计算并记录结果，然后将每组的结果拟合成一个IMF分量的预报结果，最后再将n个包含不同特征的预报结果线性拟合成最终预报解。

步骤2所述的经验模态分解算法，在分解的筛选过程中，采用设置两个连续的迭代结果之间的标准差SD的大小在0.2到0.3之间，且最大迭代次数为200作为追加的停止条件，采用极值延拓法来拟合序列的端点极值。

步骤(3)中所述的集合预报属于数值模式集合预报的一种，相对单一确定性预报，考虑了数值模式集合预报中不同组合方案对预报造成的误差。

实施例：

如图1所示，本发明实例包括如下步骤：

步骤10将输入的气候时间序列，通过经验模态分解算法进行多时间尺度分解，得到若干个本征模态函数分量(Intrinsic Mode Function，IMF)和一个趋势分量。不同分量代表了序列在不同时间周期上的波动规律，更能准确反映原序列的变化，并保留时间序列本身特征。

如图2所示，经验模态分解的处理流程，具体包括以下步骤：

步骤101用序列X(t)代表气候时间序列，t为年份，并对其进行Z-SCORE标准化处理，得到序列S(t)。

步骤101a判断S(t)的极值数是否大于2。如果是则说明序列有平稳化的必要，执行步骤102，否则说明序列是平稳的，无须分解或者分解已经结束，执行步骤108。

步骤102令H(t)＝S(t)，开始在H(t)中进行提取IMF分量的迭代操作。

步骤103找出序列H(t)中所有的局部极大值点和极小值点。

步骤104对所有的极大值和极小值点，通过三次样条插值法形成上包络线H_max(t)和下包络线H_min(t)，使H(t)满足H_min(t)≤H(t)≤H_max(t)。其中对于边界点极值的拟合问题采用端点的极值延拓法解决。

步骤105计算上下包络线的算术平均值

并从H(t)中减去，即H(t)＝H(t)-m(t)。

步骤105a判断H(t)是否已满足成为IMF分量的判定条件，包括以下三点：

①极值点数和过零点数一致或者至多相差一个；

②上下包络线的算术平均值接近于0；

③限制两次连续的处理结果之间的标准差SD的取值在0.2～0.3之间，其中

$\mathrm{SD}=\underset{t=0}{\overset{T}{\mathrm{\Σ}}}\left[\frac{{|H_{(k-1)}\left(t\right)-H_k\left(t\right)|}^2}{{H^2}_{(k-1)}\left(t\right)}\right]$

其中，H_(k-1)(t)和H_k(t)是在步骤103至步骤105中连续两次的处理结果，k代表这一过程的处理次数。SD的取值范围通常限制在0.2～0.3之间。

若同时满足上述三个条件，执行步骤106，否则返回继续执行步骤103。

步骤106成功提取一个IMF分量，即IMF_i(t)＝H_k(t)。

步骤107从S(t)中减去提取的IMF分量，即S(t)＝S(t)-IMF₁(t)，并执行步骤101a判断S(t)的极值数。

步骤108此时所有的IMF分量都被提取出，剩余的S(t)则表现为一个单调或近似单调的趋势项，称为趋势分量R_n(t)。

步骤109时间序列X(t)实现经验模态分解，即

$X\left(t\right)=\underset{i=1}{\overset{n}{\mathrm{\Σ}}}{\mathrm{IMF}}_i\left(t\right)+R_n\left(t\right)$

步骤20和步骤30如图3所示，具体说明如下：

步骤20，对步骤10所述的各分量用数值集合预报方法进行处理，设原序列长度为N，要求预测步长为M，成员变量为S，具体步骤如下：

(1)首先通过一定规则构建出集合预报的若干预报成员，成员数小于等于N，继而组合成一个预报矩阵。本发明采用对X(t)预报量取S次不同的变长序列，构成一组S行N列的集合预报成员，由用户指定要修改的成员变量S大小。

(2)每一组待预报序列依照上述数学模型进行动态建模，然后通过步骤30的均生函数逐步回归模型得出一定长度的预测序列，包括原序列拟合数值和预报数值。

(3)将各组预报序列组成预报矩阵，将相应年份的预测或拟合数值通过一定算法调整到相同列，进行线性拟合，得出一个长度为(N+M)的预测结果，用于后面的IMF分量拟合。

步骤30采用步骤20处理所得的若干组预报矩阵，建立均生函数逐步回归模型，设时间序列为X(t)其数学模型如下：

X(t)＝{X(1)，X(2)，......X(N)}(301)

式中n为样本量。对于周期T，定义其均生函数为

$X_T\left(i\right)=T/N\underset{j=1}{\overset{T}{\mathrm{\Σ}}}X(i+\mathrm{jT})---\left(302\right)$

式中的i＝1，2，......，T；1≤T≤M；

对X_T(t)作周期性延拓，则可得到外延序列，即得到M-1个预报因子序列

公式(302)与(303)构成多元回归问题：f_T(t)为预报因子，X(t)为预报量。按逐步回归步骤得到如下预测模型

$X\left(t\right)=a_0+\underset{T=2}{\overset{M}{\mathrm{\Σ}}}{a}_i{f}_T\left(t\right)+e\left(t\right),(t=\mathrm{1,2},...,N)---\left(304\right)$

其中a₀，a_i为待定系数，e(t)为白噪声。用逐步回归分析方法求得a₀和a_i系数后，得到如下拟合预报方程。对于拟合有：

$X^{\′}\left(t\right)=a_0+\underset{i=2}{\overset{M}{\mathrm{\Σ}}}{a}_i{f}_i\left(t\right)---\left(305\right)$

对于预报有：

$X^{\′}(N+p)=a_0+\underset{i=2}{\overset{M}{\mathrm{\Σ}}}{a}_i{f}_i(N+p)---\left(306\right)$

p＝1，2，......为预报点数。

步骤40将步骤20和步骤30所述的分量预测结果通过线性组合的方式得到原气候时间序列的最终预测结果。预测结果与直接用均生函数逐步回归模型预测的结果相比，平稳性增强，不可预报性降低，根据对实际气候数值的验证，其对气候变化趋势的反映比直接单一预测更加准确，对气候的突然变化把握得更好。

Claims (4)

1.一种基于经验模态分解和数值集合预报的短期气候预测方法，其特征在于包括以下几个步骤：

(1)对数据进行z-score标准化处理：将原始输入序列A的原始值v使用z-score标准化到标值v′，即

$v^{\′}=\frac{v-\stackrel{\‾}{A}}{{\mathrm{\σ}}_A},$

其中

和σ_A分别为原始输入序列的均值和标准差；

(2)将步骤(1)所述的标值v′用EMD算法分解，得出n个IMF分量和一个趋势分量，其中n为大于1的自然数；

(3)将步骤(2)所述分解出的n个IMF分量分别用数值集合预报方法进行处理，其中行数即为成员数s；随机指定每一行去掉的元素，预报元素逐行减少，最终构成一个拥有s行的变长序列组；

(4)最后对步骤(3)所述的由IMF分量处理所得的序列组构建均生函数模型，预测计算并记录结果，然后将每组的结果拟合成一个IMF分量的预报结果，最后再将n个包含不同特征的预报结果线性拟合成最终预报解。

2.根据权利要求1所述的基于经验模态分解和数值集合预报的短期气候预测方法，其特征在于，步骤(2)所述的经验模态分解算法EMD，在分解的筛选过程中，采用设置两个连续的迭代结果之间的标准差SD的大小在0.2到0.3之间，且最大迭代次数为200作为追加的停止条件，采用极值延拓法来拟合序列的端点极值。

3.根据权利要求1所述的基于经验模态分解和数值集合预报的短期气候预测方法，其特征在于，步骤(2)所述的n个IMF分量和一个趋势分量中，第一个IMF分量对应随机成分，第一至第n个IMF分量对应周期成分，趋势分量即残差Rn对应趋势成分。

4.根据权利要求1所述的基于经验模态分解和集合预报技术的短期气候预测方法，其特征在于，步骤(3)中所述每一行中的元素，零值不参与预报。

Images