本发明涉及天气预测领域，尤其涉及一种基于显式遗传算法和奇异谱分析的气象水文要素预报方法。适用于利用历史站点气象水文要素时间序列数据（包括观测数据或预报数据），开展单站气象水文要素的预报。气象水文要素包括风（包括风速或纬向风和径向风）、浪、海表温度、相对湿度、降水、海表气压异常、潜流等。
背景技术：
海面风为一个重要的海洋学物理变量。海面风预报包括两方面，一是海面风速预报，另一方面为海面风分量（又可称为风矢量）预报（包括纬向风和径向风）。海面风速的预报对许多海洋活动的计划、实施等都非常重要；同时，海面风分量为海洋大气模式进行海洋状态预报的重要变量和输入信息。
目前，对于海面风分量预报，主要依赖于数值模式，并取得了很大成功，但数值模式仍然存在包括物理过程难以准确描绘，初始状态变量场不准确等缺陷，影响着海面风场预报的精度。尤其是在我们仅关于某一特定站点海面风的预报时，数值模式预报计算量很大，并且需要包括初始时刻的湿度场、温度场等在内的大量输入信息。
因此，发展不依赖于复杂数值模式且计算量较小的海面风预报方法，一直以来受到气象人员的广泛关注。国内外学者先后开展了多种方法研究，包括线性回归方法、人工神经网络方法，遗传算法（GA）等等。其中，GA算法为目前一种先进有效的预报算法，其主要基于达尔文进化理论（Szpiro，1997；Alvarea，2003），为一种有效的非线性数据拟合方法。GA算法的一个重要优点在于其仅需要利用过去若干时刻的观测资料信息，同时，其可以提供显式的预报方程（Sheng，2012）。
目前，GA算法在气象领域开展了包括阿拉伯海海表温度（Alvarez等，2000），印度夏季降雨（Kishtawal等，2003），Ligurian海海表温度和海平面异常（Alvarez等，2004），北印度洋波高（Basu等，2005；Sinha等，2013）和涌浪（Remaya等，2012）等变量在内的预报应用研究。同时，GA算法亦应用于北印度洋的站点in situ资料和散射计资料海面风预报，结果表明GA算法相对持续模型预报在1到3天预报中具有明显的优越性（Basu等，2005）。由于GA算法的特性，不同站点和不同海域对应的海面风预报模型并不相同，均需要重新构建（Basu等，2005；Sinha等，2013）。
本发明中，针对GA算法在中国南海（SCS）站点海面风预报中开展应用研究。Basu等（2005）中，构建了北印度洋站点海面风的GA算法预报模型。本发明研究相对Basu等（2005）研究工作，不仅在于研究海域不同，同时，在实施GA算法前，本发明采用奇异谱分析（SSA）方法对海面风时间序列实施噪音滤除，这在Basu等（2005）研究中并未采用。其基本原因是由于GA算法理论上仅能用于预报确定性不包含噪音的时间序列。然而，散射计海面风资料由于其观测误差和反演误差的存在，不可避免的存在一定程度的噪音，这些噪音的存在会在一定程度上影响GA算法预报方程的精度，进而影响预报效果。因此，在对海面风时间序列实施GA算法构建预报模型前实施噪音滤除非常有必要，而SSA方法为目前削弱噪音的最有效方法（Sharma等，2007）。
技术实现要素：
本发明所要解决的技术问题是针对上述现有技术提供一种较高优越性的基于显示遗传算法和奇异谱分析的站点海面风预报方法，该方法仅利用站点气象水文要素过去若干时次的资料时间序列数据，奇异谱分析分离要素时间序列的噪音，遗传算法对噪音分离后的要素滤波和残余时间序列构建相应的显式预报模型，最终实现利用站点过去若干时次资料开展气象水文要素预报。该方法计算量小，精度可靠，可为气象预报保障提供一种新的方法和资料。
本发明解决上述问题所采用的技术方案为：一种基于显式遗传算法和奇异谱分析的气象水文要素预报方法，首先，对不同站点的海面风时间序列（包括风速时间序列和风分量时间序列）实施SSA方法滤除噪音；之后，对SSA滤波后的滤波非线性时间序列和扰动时间序列分别实施GA算法进行预报；最后，利用预报的滤波非线性时间序列和扰动时间序列构建预报风时间序列，并与真实风进行比较以便评估预报质量。
GA算法的基本理论在前人工作中有详细说明（Alvarez，2003；Basu等，2005），在此仅作简单介绍。
简单来说，对于时间序列，存在投影满足，
其中，表示嵌入维数，表示单位时间间隔步长。方程（1）为进行时步预报的预报方程，其中，投影利用过去的观测时间序列得到。
GA算法主要通过生物进化理论对投影进行有效估计。其首先通过随机方式产生个种群，之后进行如下几个主要步骤操作，包括：（1）初始化；（2）计算拟合度；（3）种群基因排序；（4）种群选择；（5）复制和交叉；（6）变异。关于该算法的具体实施过程和细节见Sharma等（2007）。
第个方程种群拟合度计算可表示为
其中，，表示训练样本总的个数。
第个种群强度指数为
其中，
其中，表示训练样本的平均值。为预报值和观测值之间相关系数的平方，其可表示为第个种群训练样本的方差百分比。可知，越趋近于1，第个种群就越好。上述过程不断进行循环，直到满足停止条件，如强度指数不再增加等。当得到最佳种群后，相应的强度指数表示为。
由于散射计海面风资料存在一定程度观测误差和反演误差，因此其站点海面风时间序列不可避免存在一定程度噪音。由于GA算法主要基于Takens理论，理论上仅可应用于时间序列无噪音情形，时间序列存在一定程度噪音时将影响算法的有效性，因此在进行GA算法预报前削弱噪音影响非常有必要。SSA方法为目前广泛应用于削弱时间序列噪音的最有效方法，其详细理论说明具体见Alvarez（2003）。
其基本原理简单表述为，对于时间序列，可构建数组，具体形式为，其中，表示时间序列总的个数，表示滤波窗口长度。然后，计算的协方差矩阵，并计算相应的特征值和特征矩阵。最后，对特征值按大小进行排序，选择前若干个特征值及其对应的特征矩阵可重构得到滤波后新的时间序列，认为其余特征值表示噪音。本发明中滤波窗口，用以包括年循环。同时，选择前80个特征值重构滤波时间序列。
本发明的原理是利用数学领域先进的非线性时间序列预报方法——显式遗传算法，应用于站点气象水文预报领域，在采用奇异谱分析有效分离要素时间序列数据噪音的基础上，利用显式遗传算法构建气象水文要素时间序列的预报模型。
基于显式遗传算法和奇异谱分析的站点气象水文要素预报方法包括以下几个步骤：
步骤一：对站点非线性气象水文要素时间序列实施奇异谱分析分离噪音（分解为要素滤波时间序列和残余时间序列）
对于站点非线性气象水文要素时间序列，简化为。
利用时间序列数据，构建滞后自协方差矩阵，具体形式为：
其中，表示时间序列总的个数，表示滤波窗口长度。
计算矩阵的协方差矩阵，具体形式为：
计算协方差矩阵对应的特征值（按由大到小进行排序）和特征向量。
计算的状态向量在特征向量上的投影，具体形式为：
其中，特征向量反映了序列的时间演变型，称为时间（简称）；是所反映的时间演变型在原序列的时段的权重，称为时间主成分（简称）。
重构滤除噪音后的要素滤波时间序列：选择前若干个特征值和对应的特征矩阵，重构要素滤波时间序列，其中，表示特征值截断数。由第个和计算得到，具体形式为：
计算要素残余时间序列：要素非线性时间序列与要素滤波时间序列之差，即。
步骤二：对要素滤波和残余时间序列实施遗传算法构建显示预报模型
对于站点要素滤波或残余时间序列，
时间序列预报模型一般形式为：
其中，表示预报模型，表示时间序列时间间隔，表示提前预报的天数，为嵌入维数，表示利用过去资料的天数。
遗传算法用以最优估计预报方程，主要包括以下几个过程：
初始化：通过随机方式产生个种群，每个种群包含的系数及其运算符号（包括基本算子及其它运算算子）。
方程初始种群可表示为
其中，为实数或变量，表示四个基本算数算子之一，即加、减、除和乘。算子和变量随机指定：数值为一系列实数，其服从在之间的均一分布。
拟合度计算：其具体计算表达式如下
其中，表示数据训练样本个数。
种群排序：利用计算的拟合度，对种群按照拟合度降序的方法对其进行排序。
种群选择：对种群进行组合。按照计算的拟合度大小，拟合度最大的种群优先选择另一种群组成组合，选择过程随机选择，其选择的概率与其拟合度成比例。之后，拟合度排第二的种群进行选择组合，并按照拟合度大小一次进行选择组合，直到前个种群选择配对完毕。剩余的个种群消失。共产生个组合。
复制和交叉：父代信息通过以下方式传递给下一代：种群选择中每个组合（共个组合）中均具有4个父代，下一代中前两个子代跟父代相同，即复制；后两个子代通过后两个父代算子或变量的交换组合而成，算子或变量的交换是随机的，即交叉。
变异：一定数量的种群（百分比很小）会产生变异。种群的选择是随机的，同时变异的可能只是单个算子或变量。按拟合度排序排在前面的种群不会发生变异，以便保证他们的信息不会产生丢失。
遗传算法通过对上述拟合度计算、种群排序、种群选择、复制和交叉以及变异五个过程不断进行迭代循环，直到达到迭代终止条件。迭代终止条件为，达到最大迭代次数，或者强度指数不再增加。
强度指数具体计算表达式为：
其中，
遗传算法迭代终止后，得到最优预报模型。
针对要素滤波和残余时间序列，设定合理的遗传算法参数设置，包括方程数即种群数、变量和算子可允许总数、利用过去历史资料长度范围、变异概率、最大迭代次数等，分别实施遗传算法，得到滤波和残余时间序列的预报模型，分别记为和。
步骤三：利用构建的预报模型开展站点气象水文要素预报
利用显式的气象水文要素滤波和残余时间序列预报模型，构建要素预报模型。
基于气象水文要素预报模型，利用过去若干时次要素数据，得到要素预报时间序列，计算预报模型的均方根预报误差
与现有技术相比，本发明的优点在于：
本发明涉及的遗传算法能提供站点气象水文要素的显式预报方程，突破传统遗传算法隐式预报方程的局限；同时，针对站点气象水文要素数据不可避免存在噪音特点，结合奇异谱分析分离时间序列的噪音，进一步确保了构建的预报模型的有效性和精度。
本发明仅利用站点气象水文要素过去若干时次的资料，只需构建一次预报模型，计算量小，与目前成功的数值预报技术相比，输入信息和计算量均有极大程度的减小，同时试验结果亦表明本发明涉及算法预报精度可靠，具有很强的业务化应用前景。
具体实施方式
以下结合实施例对本发明作进一步详细描述。
本实施例中的一种基于显式遗传算法和奇异谱分析的气象水文要素预报方法，采用数据为日平均的QuikSCAT散射计海面风场数据，包括风分量（径向风和纬向风）及风速（风速由径向和纬向风计算得到）。该散射计海面风场数据由法国海洋数据开发中心（IFREMER）发布。与欧洲中期天气预报中心（ECMWF）预报海面风场进行日平均、周平均和月平均对比分析表明，该散射计海面风场数据可很好呈现全球海域的平均风场特征。同时，与太平洋和大西洋的热带大气/海洋（TAO）和美国国家浮标数据中心（NDBC）日平均浮标数据对比分析表明数据间均方根差异较小。这充分表明该散射计海面风场数据的高精度特性和有效性（Basu等, 2005）。风速预报对于海洋活动，如船舶航行路线规划等非常重要，同时，风分量预报对于数值模式预报海洋状态等亦非常重要。
本发明以站点海面风（纬向风和径向风）预报展示具体实施过程。选择中国南海海域某一站点的散射计海面风时间序列观测数据，其对应经纬度为，时间序列范围为1999年7月20日-2009年11月21日的日平均数据，共3753个数据。散射计数据包含十年数据，可更好考虑海面风的年际变化等特征。具体包括以下步骤：
步骤一、对站点非线性海面风时间序列实施奇异谱分析分离噪音（分解为海面风滤波时间序列和残余时间序列）
设置滤波窗口数为和特征值截断数为，对于海面风时间序列实施奇异谱分析。
奇异谱分析过程如下：
计算海面风时间序列的滞后自协方差矩阵；
计算矩阵的协方差矩阵和其对应的特征值（按由大到小进行排序）和特征向量；
计算的状态向量在特征向量上的投影；
重构滤除噪音后的海面风滤波时间序列，，；
计算海面风残余时间序列。
奇异谱分析将非线性海面风（包括纬向风和径向风）时间序列分解为海面风滤波和残余时间序列，分别记为、和、。
输出结果为：
海面纬向风滤波时间序列
海面纬向风残余时间序列
海面径向风滤波时间序列
海面径向风残余时间序列
步骤二、对海面风滤波和残余时间序列实施遗传算法构建显示预报模型
利用上步骤中的海面风滤波和残余时间序列，设定合理的遗传算法参数设置，包括方程数即种群数设置为60、变量和算子可允许总数设置为20、利用过去历史资料长度范围设置为2-25、变异概率设置为0.01、最大迭代次数设置为5000，分别实施遗传算法，可构建1到5天海面风滤波和残余时间序列预报模型，在此仅展示1天预报模型。
遗传算法过程包括：
初始化：通过随机方式产生个种群，每个种群包含要素过去变量值的系数及其运算符号（包括基本算子及其它运算算子）。方程初始种群可表示为，其中，为实数或变量，表示四个基本算数算子之一，即加、减、除和乘。算子和变量随机指定：数值为一系列实数，其服从在之间的均一分布。
拟合度计算：计算表达式为
种群排序：利用计算的拟合度，对种群按照拟合度降序的方法对其进行排序。
种群选择：按照计算的拟合度大小，对种群进行组合。共产生个组合。
复制和交叉：种群选择中每个组合（共个组合）中均具有4个父代，下一代中前两个子代跟父代相同，即复制；后两个子代通过后两个父代算子或变量的交换组合而成，算子或变量的交换是随机的，即交叉。
变异：一定数量的种群（百分比很小）会产生变异。种群的选择是随机的，同时变异的可能只是单个算子或变量。按拟合度排序排在前面的种群不会发生变异，以便保证他们的信息不会产生丢失。
遗传算法通过对上述拟合度计算、种群排序、种群选择、复制和交叉以及变异五个过程不断进行迭代循环，直到达到迭代终止条件。迭代终止条件为，达到最大迭代次数，或者强度指数不再增加。强度指数计算表达式为：，其中，。
遗传算法构建站点海面风（包括纬向风和径向风）滤波时间序列和残余时间序列的预报模型，分别记为、和、。
输出结果为（仅给出1天预报模型）：
海面纬向风滤波时间序列预报模型
海面纬向风残余时间序列预报模型
海面径向风滤波时间序列预报模型
海面径向风残余时间序列预报模型
其中，表示利用预报时刻过去第个时间间隔的资料变量值，、和含义类同。
步骤三、利用显式的海面风非线性和残余时间序列预报模型开展站点海面风预报
利用上步骤中的显式海面风非线性和残余时间序列预报模型，构建海面风预报模型，计算公式为
构建站点海面风（包括纬向风和径向风）预报模型，分别记为和。
输出结果为：
海面纬向风预报模型
海面径向风预报模型
利用上述预报模型，代入散射计海面风时间序列数据，得到海面风预报时间序列。计算海面风预报时间序列的整个时间区间内的均方根预报误差，计算公式为
计算站点海面风（包括纬向风和径向风）预报均方根误差，分别记为和。
海面纬向风均方根预报误差
海面径向风均方根预报误差。
除上述实施例外，本发明还包括有其他实施方式，凡采用等同变换或者等效替换方式形成的技术方案，均应落入本发明权利要求的保护范围之内。