CN103577694A - 一种基于多尺度分析的水产养殖水质短期组合预测方法 - Google Patents

Abstract

本发明公开了一种基于多尺度分析的水产养殖水质短期组合预测方法，在线采集水质时间序列数据并进行修复；采用经验模态分解对所选的水质时间序列样本集数据分解成不同频率尺度的IMF分量和残差rn分量；将IMF和rn分量进行分类，根据分类特征分别选择人工蜂群优化最小二乘支持向量回归机、BP神经网络和自回归滑动平均模型进行预测，最后将各部分结果加权求和得到水质时间序列预测结果。本发明通过经验模态分解将原水质时间序列数据分解为不同时频尺度分量，更能准确的把握原水质序列内部的变化情况；并将人工蜂群优化最小二乘支持向量回归机、BP神经网络和自回归滑动平均模型模型优势互补相结合，有效的提高组合预测模型性能。

Description

一种基于多尺度分析的水产养殖水质短期组合预测方法

技术领域

本发明属于智能信息处理和水产养殖技术交叉领域，尤其涉及一种基于多尺度分析的水产养殖水质短期组合预测方法。

背景技术

集约化水产养殖水体是水产品栖息场所，养殖水质的好坏直接决定着水产品的生长状况及其产品质量。而水产养殖水质预测是谋划养殖区域水环境规划、现代渔业健康养殖精准化管理的重要基础性工作，准确的水质预测对养殖水质科学化调控、防范水质恶化、水产品疾病爆发具有十分重要的经济价值和现实意义。

目前，水质预测主要有基于机理预测模型和基于数值定量预测模型两大类。而基于机理预测模型要求测量水质参数较多，计算量大，误差累积率高，不满足水产养殖企业对水质短期预测的需求。基于数值定量预测方法常采用多项式回归、数理统计法、灰色系统理论法、神经网络模型法、水质模拟模型法等方法，它们各自有各自的研究特点和使用条件，虽然在水质预测取得了一定成效，但是单独预测的效果不是很突出，预测精度也比较低。此外，水产养殖水体受天气、人类活动、水生动植物、理化因子等多种因素作用，使养殖水质存在其非线性、大时滞、不确定性等特征。因此，如何根据集约化水产养殖水质特性，将不同的预测方法进行优势互补，构建一种快速、准确，且能为集约化水产养殖水质精准化调控管理提供决策支持的新型预测方法，是目前水产养殖业现代化管理亟需解决的重要问题之一。

发明内容

本发明实施例的目的在于提供一种基于多尺度分析的水产养殖水质短期组合预测方法，旨在解决现有的水质预测方法存在的预测精度低、不准确的问题。

本发明实施例是这样实现的，一种基于多尺度分析的水产养殖水质短期组合预测方法，该基于多尺度分析的水产养殖水质短期组合预测方法包括以下步骤：

步骤一，在线获取集约化水产养殖池塘的水质数据，按照时间顺序排列成水质时间序列数据X，水质时间序列数据X记为：X＝{x_t，t＝1，2，...，N}，其中，N为水质时间序列数据的个数且N为正整数；

步骤二，对在线采集的水质原始数据修复处理，根据数据特征对历史水质时间序列数据中的残缺值进行修复处理；

步骤三，经验模态本征分解：采用EMD算法对已经修复处理水质参数时间序列进行分解，得到n个IMF分量和一个残差r_n，其中n是大于1的自然数；

步骤四，IMF和r_n分量分类预测：根据不同尺度频率滤波的特点，将IMF分量分为高频滤波、中频滤波和低频滤波，分别使用人工蜂群优化最小二乘支持向量回归机、BP神经网络和自回归滑动平均模型进行建模预测；

将分解得到的几个较高频率的IMF分量作为高频滤波；将分解得到的中级几个一般频率的IMF分量作为中频滤波；将后几个频率较低和残差rn作为低频滤波，然后针对不同时频的滤波特性进行分类预测，对高频滤波部分，使用人工蜂群优化最小二乘支持向量回归机进行建模预测，对中频滤波部分，使用BP神经网络进行建模预测，对低频滤波部分，使用自回归滑动平均模型进行建模预测；

步骤五，组合模型权重计算：假定某一水质参数时间序列为(x₁，x₂，...，x_t)，则模型权重计算式为：

$w_t^m=\frac{\frac{1}{\sqrt{{2\mathrm{\π\σ}}_m}}{w}_{t-1}^m{e}^{-{\left[(x_{i,t}-{\hat{x}}_{i,t}^m/{\mathrm{\σ}}_m)\right]}^2}}{\underset{j=1}{\overset{J}{\mathrm{\Σ}}}\frac{1}{\sqrt{{2\mathrm{\π\σ}}_j}}{w}_{t-1}^m{e}^{{-\left[(x_{i,t}-{\hat{x}}_{i,t}^j/{\mathrm{\σ}}_j)\right]}^2}}$

其中，

为t+1时刻第j个基本预测模型在组合模型中的权重，j＝1，2，...，J，J为基本预测模型的个数，

为第j个模型第i个分量的预测值，σ_i为x_i，t

的方差；

步骤六，将各模型预测结果进行加权求和，得到最终的集约化水产养殖水质预测值

${\hat{x}}_{t+1}^{*}=\underset{j=1}{\overset{J}{\mathrm{\Σ}}}{w}_t^j{x}_{i,t+1}^j$

式中，

为第t+1时刻水质参数最终预测结果。

进一步，在步骤二中，残缺值进行修复处理包括；

残缺数据的修补处理：

若缺失数据的前后时间间隔不大，采用线性插值的方法补上；

$x_{k+i}=x_k+\frac{i\&CenterDot;(x_{k+j}-x_k)}{j},0<i<j$

式中，x_k和x_k+j分别为已知k时刻和k+j时刻的水质参数值，x_k+i为k+i时刻缺失的水质参数取值，若时间间隔较大，则采用天气类型相同或相近相邻几天中同时刻的数据来填补；

错误的不良数据处理：

因水质数据具有连续性和时序性，前后相邻时段监测的数据一般不会发生跳变，通过比较，若该时刻水质数据变化范围在其前后水质监测值的±10％以外，则认为该数据有误，可采用均值平滑法进行水平处理；

和

分别为相邻数据误差阈值。

进一步，在步骤三中，找出水质时间序列x(t)的所有极大值点和极小值点，将所有的极大值点用三次样条函数拟合一条数据x(t)的上包络线e_up(t)；同样的把所有极小值点采用三次样条函数拟合一条数据x(t)的下包络线e_down(t)，则上下包络线的平均值m₁(t)＝(e_up(t)+e_down(t))/2；将原始水质时间序列x(t)减去上下包络线的均值m₁(t)，得到一个新的数据序列h1(t)＝x(t)-m₁(t)，并进行k次筛选，直到h₁(t)满足IMF条件，记c₁(t)＝h₁(t)，则c₁(t)为水质时间序列的第1个本征模态函数IMF分量，它表示原始水质时间序列中的最高频成分，记为：h_1(k)(t)＝h_1(k-1)(t)-m_1(k)(t)，式中h_1(k)(t)为第k次筛选得到的水质时间序列数据，h_1(k-1)(t)为第k-1次筛选得到的水质时间序列数据；如此多次筛选，原水质时间序列被分解成若干个IMF分量和一个残差分量rn，随着分解次数的增加，分量的非平稳行为逐渐减少，并且不同分量代表水质序列在不同时频上的波动规律，更能准确的反映原水质序列内部的变化情况，并保留原水质序列的固有特征。

进一步，在步骤四中，建立改进人工蜂群优化最小二乘支持向量回归机数学模型，包括以下步骤：

第一步，设置人工蜂群算法参数，蜂群总数为NS，维数为D＝2，采蜜蜂个数为N_e，守望蜂个数为N_u，食物源数量FN，且N_e＝N_u＝FN＝NP/2，最大迭代次数为T_max，最大更新次数T_limit，设当前迭代次数为t＝1；每个食物源的当前更新次数t_Fi＝0；

第二步，生成待优化预测模型参数的初始解集θ：

$\mathrm{\&theta;}=[{\mathrm{\&theta;}}_1,{\mathrm{\&theta;}}_2,...,{\mathrm{\&theta;}}_{N_s}],{\mathrm{\&theta;}}_i=({\mathrm{\&theta;}}_{i,1},{\mathrm{\&theta;}}_{i,2})$

其中，i＝1，...，N_e，θ_i为第i组初始解；θ_i，j表示第i组解中的第j个值，对应为待确定预测模型中的一个参数，可通过下式产生θ_i，j的初始解：

θ_i，j＝θ_min+rand(0，1)*(θ_max-θ_min，j)

其中，i＝1，...，N_e，j＝1，2；rand(0，1)表示产生[0，1]之间的随机数，θ_{min， j}表示θ_i，j的最小值，θ_max，j表示θ_i，j的最大值；

利用下式计算所有解θ_i对应的适应度值，

${\mathrm{fitness}}_i\left({\mathrm{\&theta;}}_i\right)=\left\{\begin{array}{cc}\frac{1}{1+{\mathrm{fit}}_i\left({\mathrm{\&theta;}}_i\right)},& {\mathrm{fit}}_i\left({\mathrm{\&theta;}}_i\right)\&GreaterEqual;0\\ 1+\left|{\mathrm{fit}}_i\left({\mathrm{\&theta;}}_i\right)\right|,& {\mathrm{fit}}_i\left({\mathrm{\&theta;}}_i\right)<0\end{array}\right.$

其中，fit_i(θ_i)表示此处水质分解后分量数据拟合值与原始水质分量数据之间的误差，即：

${\mathrm{fit}}_i\left({\mathrm{\&theta;}}_i\right)={\mathrm{\&lambda;}}_1*\sqrt{\frac{j1}{1,\mathrm{tn}}\underset{i=1}{\overset{n}{\mathrm{\&Sigma;}}}{(x_{i,t}-{\hat{x}}_{i,t}^j)}^2}+{\mathrm{\&lambda;}}_2*\frac{\max \left\{\right|x_{1,t}-{\hat{x}}_{1,t}^j|,\&CenterDot;\&CenterDot;\&CenterDot;,|x_{n,t}-{\hat{x}}_{n,t}^j\left|\right\}}{\sqrt{\frac{1}{n}\underset{i=1}{\overset{n}{\mathrm{\&Sigma;}}}{(x_{i,t}-\stackrel{\&OverBar;}{x})}^2}}$

其中，x_i，t为水质分解后的第i个分量在t时刻的真实值，

为水质分解后的第j个模型第i个分量在t时刻的拟合值，

为分解后的第i个分量的平均值，λ₁和λ₂为调节因子，分别是0到1之间的随机数，n是水质分解后分量数据的个数；

第三步，按照下式计算各个解的适应度值finess_i(θ_i(t))，其中t是指第t次迭代，并将θ_i(t)的适应度值从小到大进行排序，适应度最大的θ_i(t)为侦查蜂，θpb_i(t)为每只采蜜蜂当前发现的最好解，表示为θpb_i(t)＝θ_i(t)，i＝1，2，...，N_e，并且θpb_i，j(t)＝θ_i，j(t)，其中θpb_i，j(t)为第i只采蜜蜂在第j维的最好解；θgb(t)设定为群体中的全局最好解，

并且 $f(\mathrm{\&theta;gb}\left(t\right)=\min \{f({\mathrm{\&theta;pb}}_1\left(t\right),...,f\left({\mathrm{\&theta;pb}}_{N_e}\left(e\right)\right))\};$

第四步，采蜜蜂按下式对食物源θ_i，j(t)进行邻域搜索产生新解θ_i，j(t+1)，并进行解的替代；

$\left\{\begin{array}{cc}{\mathrm{\&theta;}}_{i,j}(t+1)=w_{i,j}\left(t\right)*{\mathrm{\&theta;}}_{i,j}\left(t\right)+c_1*r_a*({\mathrm{\&theta;pb}}_{i,j}\left(t\right)-{\mathrm{\&theta;}}_{i,j}\left(t\right))+\mathrm{\&lambda;}\&CenterDot;N(0,1)+c_2*r_b*({\mathrm{\&theta;gb}}_j\left(t\right)-{\mathrm{\&theta;}}_{i,j}\left(t\right)),& f\left({\mathrm{\&theta;}}_{i,j}(t+1)\right)<f\left({\mathrm{\&theta;}}_{i,j}\left(t\right)\right)\\ {\mathrm{\&theta;}}_{i,j}(t+1)={\mathrm{\&theta;}}_{i,j}\left(t\right),& f\left({\mathrm{\&theta;}}_{i,j}(t+1)\right)>f\left({\mathrm{\&theta;}}_{i,j}\left(t\right)\right)\end{array}\right.$

$w_{i,j}=w_{\max }-\frac{(w_{\max }-w_{\min })t}{T_{\max }}$

其中，θ_i，j(t+1)为产生的新解，θ_i，j(t)为原来的旧解；θpb_i，j(t)为θ[b_i(t)在第j维最好解；θbg_j(t)为θbg(t)的第j维全局最好解，(i=1，2，...，N_e；j＝1，2，...，D)，c₁和c₂为学习因子，N(0，1)为产生高斯分布的随机数，W_max和W_min分为为惯性权重w的上下限，λ为增量系数，r_a和r_b是[-1，1]之间的随机数，t指第t次迭代；

第五步，各跟随蜂按照下式计算每个食物源被选择的概率，采用轮盘赌选择策略选择一个食物源并按照下式进行探索，产生新解并进行解的替代；

$P_i=\frac{{\mathrm{fitness}}_i\left({\mathrm{\&theta;}}_i\right)}{\underset{k=1}{\overset{\mathrm{Ne}}{\mathrm{\&Sigma;}}}{\mathrm{fitness}}_k\left({\mathrm{\&theta;}}_k\right)}$

$\left\{\begin{array}{cc}{\mathrm{\&theta;}}_{i,j}(t+1)={\mathrm{\&theta;}}_{i,j}\left(t\right)+{\mathrm{\&eta;}}_{i,j}({\mathrm{\&theta;}}_{i,j}\left(t\right)-{\mathrm{\&theta;}}_{k,j}\left(t\right))+{\mathrm{\&xi;}}_{i,j}({\mathrm{\&theta;gb}}_j\left(t\right)-{\mathrm{\&theta;}}_{k,j}\left(t\right)),& f\left({\mathrm{\&theta;}}_{i,j}(t+1)\right)<f\left({\mathrm{\&theta;}}_{i,j}\left(t\right)\right)\\ {\mathrm{\&theta;}}_{i,j}(t+1)={\mathrm{\&theta;}}_{i,j}\left(t\right),& f\left({\mathrm{\&theta;}}_{i,j}(t+1)\right)>f\left({\mathrm{\&theta;}}_{i,j}\left(t\right)\right)\end{array}\right.$

${\mathrm{\&xi;}}_{i,j}=\frac{2t}{T_{\max }+t}$

其中，θ_i，j(t+1)为产生的新解，θ_i，j(t)为原来的旧解；η_ij是[-1，1]之间的随机数，θ_k为跟随蜂选择的采蜜蜂，θ_k，j为采蜜蜂θ_k的第j为元素，并且i＝1，2，...，N_u；j＝1，2，...，D；t_Fi＝t_Fi+1，ξ_i，j为自适应系数；

第六步，若t_Fi超过预设最大更新次数T_limit时，探索得到的θ_i仍未能改进的食物源适应度，则派出侦察蜂采用下式进行对其重新初始化，获得新解并进行解的替代；

其中，i＝1；j＝1，2，...，D，θ_i，j(t+1)为产生的新解，θ_i，j(t)为原来的旧解；

是列维概率分布产生的随机搜索步长，θpb_i，j(t)和分别为解θ_i的第j维的最大值；

第七步，若迭代次数t小于最大迭代次数t_max，则返回第三步，否则结束并输出最优解θgb(t)，即得到模型正则化参数γ和核函数参数σ的最优解组合，以获得最佳的正则化参数γ和核函数参数σ，建立人工蜂群优化最小二乘支持向量回归机水质预测模型。

进一步，在步骤四中，对中频滤波部分，使用BP神经网络进行建模预测的具体步骤为：

第一步，接收中频滤波序列数据，网络初始化；

第二步，计算误差和梯度向量，进行模型的权值修正；

第三步，对精度和迭代次数是否满足算法终止条件作出判断，满足进行下一步，不满足返回第二步；

第四步，BP神经网络预测模型建立完成。

进一步，在步骤四中，对低频滤波部分，使用自回归滑动平均模型进行建模预测的具体方法为：

第一步，接收低频序列数据，进行模型的结构辨识；

第二步，进行模型参数估计、诊断和检验；

第三步，对模型是否合适作出判断，是则建立自回归滑动平均预测模型，否返回第一步。

本发明提供的基于多尺度分析的水产养殖水质短期组合预测方法，通过经验模态本征分解和人工蜂群优化最小二乘支持向量回归机、BP神经网络和自回归滑动平均模型相结合的基于多尺度分析的水产养殖水质短期组合预测方法；

本发明具有以下有益效果：

(1)通过经验模态本征分解，将由历史养殖水质时间序列数据依尺度分解成不同尺度的n个IMF分量和一个残差rn，消除了噪音的干扰；

(2)在人工蜂群优化算法中，采用高斯变异、Levy变异和自适应变异算子，有效提高了人工蜂群算法全局搜索的能力；

(3)以经验模态本征分解的IMF分量和一个残差rn划分为高频滤波、中频滤波和低频滤波，根据滤波特性分别选择不同的预测模型进行预测，各模型优势互补，大大提高了组合预测方法的性能；

(4)为水产养殖人员做出准确的水质优化调控决策提供可靠、有效的技术支撑；

本发明将多种预测方法组合起来，综合利用各种方法的预测优势，提供预测精度，有效提高了水产养殖水质预测精度，且实现容易、计算量小，操作便捷，能够为实现水产福利养殖和节能减排的水质最优智能控制提供决策依据。

附图说明

图1是本发明实施例提供的基于多尺度分析的水产养殖水质短期组合预测方法的流程图；

图2是本发明实施例提供的基于人工蜂群算法优化最小二乘支持向量回归机预测模型的流程图；

图3是本发明实施例提供的基于BP神经网络预测模型的流程图；

图4是本发明实施例提供的基于自回归滑动平均预测模型的流程图。

具体实施方式

为了使本发明的目的、技术方案及优点更加清楚明白，以下结合实施例，对本发明进行进一步详细说明。应当理解，此处所描述的具体实施例仅仅用以解释本发明，并不用于限定本发明。

下面结合附图及具体实施例对本发明的应用原理作进一步描述。

如图1所示，本发明实施例的基于多尺度分析的水产养殖水质短期组合预测方法包括以下步骤：

S101：在线采集水质时间序列数据并进行修复；

S102：采用经验模态分解对所选的水质时间序列样本集数据分解成不同频率尺度的IMF分量和残差rn分量；

S103：将IMF和rn分量进行分类，根据分类特征分别选择人工蜂群优化最小二乘支持向量回归机、BP神经网络和自回归滑动平均模型进行预测；

S104：最后将各部分结果加权求和得到水质时间序列预测结果。

本发明的具体步骤为：

步骤一，在线获取集约化水产养殖池塘的水质数据，按照时间顺序排列成水质时间序列数据X，水质时间序列数据x记为：X＝(x_t，t＝1，2，...，N}，其中，N为水质时间序列数据的个数且N为正整数；

步骤二，对在线采集的水质原始数据修复处理，根据数据特征对历史水质时间序列数据中的残缺值进行修复处理；

在实际工作中，为消除不规则的噪声数据和修补缺失数据，减少“不良数据”对预测结果的影响，对采集的数据进行修复处理；

1)残缺数据的修补处理

若缺失数据的前后时间间隔不大，采用线性插值的方法将其补上；

$x_{k+i}=x_k+\frac{i\&CenterDot;(x_{k+j}-x_k)}{j},0<i<j$

式中，x_k和X_k+j分别为已知k时刻和k+j时刻的水质参数值，X_k+i为k+i时刻缺失的水质参数取值，若时间间隔较大，则采用天气类型相同或相近相邻几天中同时刻的数据来填补；

2)错误的“不良数据”处理

因水质数据具有连续性和时序性，前后相邻时段监测的数据一般不会发生跳变，通过比较，若该时刻水质数据变化范围在其前后水质监测值的±10％以外，则认为该数据有误，可采用均值平滑法进行水平处理；

和

分别为相邻数据误差阈值；

步骤三，经验模态本征分解：采用EMD算法对已经修复处理水质参数时间序列进行分解，得到n个IMF分量和一个残差r_n，其中n是大于1的自然数；

在实际工作中，找出水质时间序列x(t)的所有极大值点和极小值点，将所有的极大值点用三次样条函数拟合一条数据x(t)的上包络线e_up(t)；同样的把所有极小值点采用三次样条函数拟合一条数据x(t)的下包络线e_down(t)，则上下包络线的平均值m₁(t)＝(e_up(t)+e_down(t))/2；将原始水质时间序列x(t)减去上下包络线的均值m₁(t)，得到一个新的数据序列h₁(t)＝x(t)-m₁(t)，并进行k次筛选，直到h₁(t)满足IMF条件，记c₁(t)＝h₁(t)，则c₁(t)为水质时间序列的第1个本征模态函数IMF分量，它表示原始水质时间序列中的最高频成分，记为：h_1(k)(t)＝h_1(k-1)(t)-m_1(k)(t)，式中h_1(k)(t)为第k次筛选得到的水质时间序列数据，h_1(k-1)(t)为第k-1次筛选得到的水质时间序列数据；如此多次筛选，原水质时间序列被分解成若干个IMF分量和一个残差分量rn，随着分解次数的增加，分量的非平稳行为逐渐减少，并且不同分量代表水质序列在不同时频上的波动规律，更能准确的反映原水质序列内部的变化情况，并保留原水质序列的固有特征；

步骤四，IMF和r_n分量分类预测：根据不同尺度频率滤波的特点，将IMF分量分为高频滤波、中频滤波和低频滤波，分别使用人工蜂群优化最小二乘支持向量回归机、BP神经网络和自回归滑动平均模型进行建模预测；

在实际工作中，将分解得到的几个较高频率的IMF分量作为高频滤波；将分解得到的中级几个一般频率的IMF分量作为中频滤波；将后几个频率较低和残差rn作为低频滤波，

然后针对不同时频的滤波特性进行分类预测，如图2所示，对高频滤波部分，使用人工蜂群优化最小二乘支持向量回归机进行建模预测，如图3所示，对中频滤波部分，使用BP神经网络进行建模预测，如图4所示，对低频滤波部分，使用自回归滑动平均模型进行建模预测；

步骤五，组合模型权重计算：假定某一水质参数时间序列为(x₁，X₂，...，X_t)，则模型权重计算式为：

$w_t^m=\frac{\frac{1}{\sqrt{{2\mathrm{\&pi;\&sigma;}}_m}}{w}_{t-1}^m{e}^{-{\left[(x_{i,t}-{\hat{x}}_{i,t}^m/{\mathrm{\&sigma;}}_m)\right]}^2}}{\underset{j=1}{\overset{J}{\mathrm{\&Sigma;}}}\frac{1}{\sqrt{{2\mathrm{\&pi;\&sigma;}}_j}}{w}_{t-1}^m{e}^{{-\left[(x_{i,t}-{\hat{x}}_{i,t}^j/{\mathrm{\&sigma;}}_j)\right]}^2}}$

其中，为t+1时刻第j个基本预测模型在组合模型中的权重，j=1，2，...，J，J为基本预测模型的个数，

为第j个模型第i个分量的预测值，σ_i为x_i，t 的方差；

步骤六，将各模型预测结果进行加权求和，得到最终的集约化水产养殖水质预测值

${\hat{x}}_{t+1}^{*}=\underset{j=1}{\overset{J}{\mathrm{\&Sigma;}}}{w}_t^j{x}_{i,t+1}^j$

式中，为第t+1时刻水质参数最终预测结果；

在步骤四中建立改进人工蜂群优化最小二乘支持向量回归机数学模型，包括以下步骤：

第一步，设置人工蜂群算法参数，蜂群总数为NS，维数为D＝2，采蜜蜂个数为N_e，守望蜂个数为N_u，食物源数量FN，且N_e＝N_u＝FN＝NP/2，最大迭代次数为T_max，最大更新次数T_limit，设当前迭代次数为t＝1；每个食物源的当前更新次数t_Fi＝0；

第二步，生成待优化预测模型参数的初始解集θ：

$\mathrm{\&theta;}=[{\mathrm{\&theta;}}_1,{\mathrm{\&theta;}}_2,...,{\mathrm{\&theta;}}_{N_s}],{\mathrm{\&theta;}}_i=({\mathrm{\&theta;}}_{i,1},{\mathrm{\&theta;}}_{i,2})$

其中，i＝1，...，N_e，θ_i为第i组初始解；θ_i，j表示第i组解中的第j个值，对应为待确定预测模型中的一个参数，可通过下式产生θ_i，j的初始解：

θ_i，j＝θ_min+rand(0，1)*(θ_max-θ_min，j)

其中，i＝1，...，N_e，j＝1，2；rand(0，1)表示产生[0，1]之间的随机数，θ_{min， j}表示θ_i，j的最小值，θ_max，j表示θ_i，j的最大值；

利用下式计算所有解θ_i对应的适应度值，

${\mathrm{fitness}}_i\left({\mathrm{\&theta;}}_i\right)=\left\{\begin{array}{cc}\frac{1}{1+{\mathrm{fit}}_i\left({\mathrm{\&theta;}}_i\right)},& {\mathrm{fit}}_i\left({\mathrm{\&theta;}}_i\right)\&GreaterEqual;0\\ 1+\left|{\mathrm{fit}}_i\left({\mathrm{\&theta;}}_i\right)\right|,& {\mathrm{fit}}_i\left({\mathrm{\&theta;}}_i\right)<0\end{array}\right.$

其中，fit_i(θ_i)表示此处水质分解后分量数据拟合值与原始水质分量数据之间的误差，即：

${\mathrm{fit}}_i\left({\mathrm{\&theta;}}_i\right)={\mathrm{\&lambda;}}_1*\sqrt{\frac{j1}{1,\mathrm{tn}}\underset{i=1}{\overset{n}{\mathrm{\&Sigma;}}}{(x_{i,t}-{\hat{x}}_{i,t}^j)}^2}+{\mathrm{\&lambda;}}_2*\frac{\max \left\{\right|x_{1,t}-{\hat{x}}_{1,t}^j|,\&CenterDot;\&CenterDot;\&CenterDot;,|x_{n,t}-{\hat{x}}_{n,t}^j\left|\right\}}{\sqrt{\frac{1}{n}\underset{i=1}{\overset{n}{\mathrm{\&Sigma;}}}{(x_{i,t}-\stackrel{\&OverBar;}{x})}^2}}$

其中，x_i，t为水质分解后的第i个分量在t时刻的真实值，

为水质分解后的第j个模型第i个分量在t时刻的拟合值，

为分解后的第i个分量的平均值，λ₁和λ₂为调节因子，分别是0到1之间的随机数，n是水质分解后分量数据的个数；

第三步，按照下式计算各个解的适应度值finess_i(θ_i(t))，其中t是指第t次迭代，并将θ_i(t)的适应度值从小到大进行排序，适应度最大的θ_i(t)为侦查蜂，θpb_i(t)为每只采蜜蜂当前发现的最好解，表示为θpb_i(t)＝θ_i(t)，i＝1，2，...，N_e，并且θpb_i，j(t)＝θ_i，j(t)，其中θpb_i，j(t)为第i只采蜜蜂在第j维的最好解；θgb(t)设定为群体中的全局最好解，

并且 $f(\mathrm{\&theta;gb}\left(t\right)=\min \{f({\mathrm{\&theta;pb}}_1\left(t\right),...,f\left({\mathrm{\&theta;pb}}_{N_e}\left(e\right)\right))\};$

第四步，采蜜蜂按下式对食物源θ_i，j(t)进行邻域搜索产生新解θ_i，j(t+1)，并进行解的替代；

$\left\{\begin{array}{cc}{\mathrm{\&theta;}}_{i,j}(t+1)=w_{i,j}\left(t\right)*{\mathrm{\&theta;}}_{i,j}\left(t\right)+c_1*r_a*({\mathrm{\&theta;pb}}_{i,j}\left(t\right)-{\mathrm{\&theta;}}_{i,j}\left(t\right))+\mathrm{\&lambda;}\&CenterDot;N(0,1)+c_2*r_b*({\mathrm{\&theta;gb}}_j\left(t\right)-{\mathrm{\&theta;}}_{i,j}\left(t\right)),& f\left({\mathrm{\&theta;}}_{i,j}(t+1)\right)<f\left({\mathrm{\&theta;}}_{i,j}\left(t\right)\right)\\ {\mathrm{\&theta;}}_{i,j}(t+1)={\mathrm{\&theta;}}_{i,j}\left(t\right),& f\left({\mathrm{\&theta;}}_{i,j}(t+1)\right)>f\left({\mathrm{\&theta;}}_{i,j}\left(t\right)\right)\end{array}\right.$

$w_{i,j}=w_{\max }-\frac{(w_{\max }-w_{\min })t}{T_{\max }}$

其中，θ_i，j(t+1)为产生的新解，θ_i，j(t)为原来的旧解；θpb_i，j(t)为θpb_i(t)在第j维最好解；θbg_j(t)为θbg(t)的第j维全局最好解，(i=1，2，...，N_e；j=1，2，...，D)，c₁和c₂为学习因子，N(0，1)为产生高斯分布的随机数，W_max和W_min分为为惯性权重w的上下限，λ为增量系数，r_a和r_b是[-1，1]之间的随机数，t指第t次迭代；

第五步，各跟随蜂按照下式计算每个食物源被选择的概率，采用轮盘赌选择策略选择一个食物源并按照下式进行探索，产生新解并进行解的替代；

$P_i=\frac{{\mathrm{fitness}}_i\left({\mathrm{\&theta;}}_i\right)}{\underset{k=1}{\overset{\mathrm{Ne}}{\mathrm{\&Sigma;}}}{\mathrm{fitness}}_k\left({\mathrm{\&theta;}}_k\right)}$

$\left\{\begin{array}{cc}{\mathrm{\&theta;}}_{i,j}(t+1)={\mathrm{\&theta;}}_{i,j}\left(t\right)+{\mathrm{\&eta;}}_{i,j}({\mathrm{\&theta;}}_{i,j}\left(t\right)-{\mathrm{\&theta;}}_{k,j}\left(t\right))+{\mathrm{\&xi;}}_{i,j}({\mathrm{\&theta;gb}}_j\left(t\right)-{\mathrm{\&theta;}}_{k,j}\left(t\right)),& f\left({\mathrm{\&theta;}}_{i,j}(t+1)\right)<f\left({\mathrm{\&theta;}}_{i,j}\left(t\right)\right)\\ {\mathrm{\&theta;}}_{i,j}(t+1)={\mathrm{\&theta;}}_{i,j}\left(t\right),& f\left({\mathrm{\&theta;}}_{i,j}(t+1)\right)>f\left({\mathrm{\&theta;}}_{i,j}\left(t\right)\right)\end{array}\right.$

${\mathrm{\&xi;}}_{i,j}=\frac{2t}{T_{\max }+t}$

其中，θ_i，j(t+1)为产生的新解，θ_i，j(t)为原来的旧解；η_ij是[-1，1]之间的随机数，θ_k为跟随蜂选择的采蜜蜂，θ_k，j为采蜜蜂θ_k的第j为元素，并且i＝1，2，...，N_u；j=1，2，...，D；t_Fi＝t_Fi+1，ξ_i，j为自适应系数；

第六步，若t_Fi超过预设最大更新次数T_limit时，探索得到的θ_i仍未能改进的食物源适应度，则派出侦察蜂采用下式进行对其重新初始化，获得新解并进行解的替代；

其中，i＝1；j=1，2，...，D，θ_i，j(t+1)为产生的新解，θ_i，j(t)为原来的旧解；是列维概率分布产生的随机搜索步长，θpb_i，j(t)和

分别为解θ_i的第j维的最大值；

第七步，若迭代次数t小于最大迭代次数t_max，则返回第三步，否则结束并输出最优解θgb(t)，即得到模型正则化参数γ和核函数参数σ的最优解组合，以获得最佳的正则化参数γ和核函数参数σ，建立人工蜂群优化最小二乘支持向量回归机水质预测模型。

如图3所示，在步骤四中，对中频滤波部分，使用BP神经网络进行建模预测的具体步骤为：

第一步，接收中频滤波序列数据，网络初始化；

第二步，计算误差和梯度向量，进行模型的权值修正；

第三步，对精度和迭代次数是否满足算法终止条件作出判断，满足进行下一步，不满足返回第二步；

第四步，BP神经网络预测模型建立完成。

如图4所示，在步骤四中，对低频滤波部分，使用自回归滑动平均模型进行建模预测的具体方法为：

第一步，接收低频序列数据，进行模型的结构辨识；

第二步，进行模型参数估计、诊断和检验；

第三步，对模型是否合适作出判断，是则建立自回归滑动平均预测模型，否返回第一步。

本发明将多种预测方法组合起来，综合利用各种方法的预测优势，提供预测精度，有效提高了水产养殖水质预测精度，且实现容易、计算量小，操作便捷，能够为实现水产福利养殖和节能减排的水质最优智能控制提供决策依据。

以上所述仅为本发明的较佳实施例而已，并不用以限制本发明，凡在本发明的精神和原则之内所作的任何修改、等同替换和改进等，均应包含在本发明的保护范围之内。

Claims (6)

1.一种基于多尺度分析的水产养殖水质短期组合预测方法，其特征在于，该基于多尺度分析的水产养殖水质短期组合预测方法包括以下步骤：

步骤一，在线获取集约化水产养殖池塘的水质数据，按照时间顺序排列成水质时间序列数据X，水质时间序列数据X记为：X＝{x_t，t＝1，2，...，N}，其中，N为水质时间序列数据的个数且N为正整数；

步骤二，对在线采集的水质原始数据修复处理，根据数据特征对历史水质时间序列数据中的残缺值进行修复处理；

步骤三，经验模态本征分解：采用EMD算法对已经修复处理水质参数时间序列进行分解，得到n个IMF分量和一个残差r_n，其中n是大于1的自然数；

步骤四，IMF和r_n分量分类预测：根据不同尺度频率滤波的特点，将IMF分量分为高频滤波、中频滤波和低频滤波，分别使用人工蜂群优化最小二乘支持向量回归机、BP神经网络和自回归滑动平均模型进行建模预测；

将分解得到的几个较高频率的IMF分量作为高频滤波；将分解得到的中级几个一般频率的IMF分量作为中频滤波；将后几个频率较低和残差rn作为低频滤波，然后针对不同时频的滤波特性进行分类预测，对高频滤波部分，使用人工蜂群优化最小二乘支持向量回归机进行建模预测，对中频滤波部分，使用BP神经网络进行建模预测，对低频滤波部分，使用自回归滑动平均模型进行建模预测；

步骤五，组合模型权重计算：假定某一水质参数时间序列为(x₁，x₂，...，X_t)，则模型权重计算式为：

$w_t^m=\frac{\frac{1}{\sqrt{{2\mathrm{\&pi;\&sigma;}}_m}}{w}_{t-1}^m{e}^{-{\left[(x_{i,t}-{\hat{x}}_{i,t}^m/{\mathrm{\&sigma;}}_m)\right]}^2}}{\underset{j=1}{\overset{J}{\mathrm{\&Sigma;}}}\frac{1}{\sqrt{{2\mathrm{\&pi;\&sigma;}}_j}}{w}_{t-1}^m{e}^{{-\left[(x_{i,t}-{\hat{x}}_{i,t}^j/{\mathrm{\&sigma;}}_j)\right]}^2}}$

其中，

为t+1时刻第j个基本预测模型在组合模型中的权重，j=1，2，...，J，J为基本预测模型的个数，

为第j个模型第i个分量的预测值，σ_i为

的方差；

步骤六，将各模型预测结果进行加权求和，得到最终的集约化水产养殖水质预测值

${\hat{x}}_{t+1}^{*}=\underset{j=1}{\overset{J}{\mathrm{\&Sigma;}}}{w}_t^j{x}_{i,t+1}^j$

式中，

为第t+1时刻水质参数最终预测结果。

2.如权利要求1所述的基于多尺度分析的水产养殖水质短期组合预测方法，其特征在于，在步骤二中，残缺值进行修复处理包括；

残缺数据的修补处理

若缺失数据的前后时间间隔不大，采用线性插值的方法补上；

$x_{k+i}=x_k+\frac{i\&CenterDot;(x_{k+j}-x_k)}{j},0<i<j$

式中，x_k和x_k+j分别为已知k时刻和k+j时刻的水质参数值，x_k+i为k+i时刻缺失的水质参数取值，若时间间隔较大，则采用天气类型相同或相近相邻几天中同时刻的数据来填补；

错误的不良数据处理

因水质数据具有连续性和时序性，前后相邻时段监测的数据一般不会发生跳变，通过比较，若该时刻水质数据变化范围在其前后水质监测值的±10％以外，则认为该数据有误，可采用均值平滑法进行水平处理；

和

分别为相邻数据误差阈值。

3.如权利要求1所述的基于多尺度分析的水产养殖水质短期组合预测方法，其特征在于，在步骤三中，找出水质时间序列x(t)的所有极大值点和极小值点，将所有的极大值点用三次样条函数拟合一条数据x(t)的上包络线e_up(t)；同样的把所有极小值点采用三次样条函数拟合一条数据x(t)的下包络线e_down(t)，则上下包络线的平均值m₁(t)＝(e_up(t)+e_down(t))/2；将原始水质时间序列x(t)减去上下包络线的均值m₁(t)，得到一个新的数据序列h₁(t)＝x(t)-m₁(t)，并进行k次筛选，直到h₁(t)满足IMF条件，记c₁(t)＝h₁(t)，则c₁(t)为水质时间序列的第1个本征模态函数IMF分量，它表示原始水质时间序列中的最高频成分，记为：h_1(k)(t)＝h_1(k-1)(t)-m_1(k)(t)，式中h_1(k)(t)为第k次筛选得到的水质时间序列数据，h_1(k-1)(t)为第k-1次筛选得到的水质时间序列数据；如此多次筛选，原水质时间序列被分解成若干个IMF分量和一个残差分量rn，随着分解次数的增加，分量的非平稳行为逐渐减少，并且不同分量代表水质序列在不同时频上的波动规律，更能准确的反映原水质序列内部的变化情况，并保留原水质序列的固有特征。

4.如权利要求1所述的基于多尺度分析的水产养殖水质短期组合预测方法，其特征在于，在步骤四中，建立改进人工蜂群优化最小二乘支持向量回归机数学模型，包括以下步骤：

第一步，设置人工蜂群算法参数，蜂群总数为NS，维数为D＝2，采蜜蜂个数为N_e，守望蜂个数为N_u，食物源数量FN，且N_e＝N_u＝FN＝NP/2，最大迭代次数为T_max，最大更新次数T_limit，设当前迭代次数为t＝1；每个食物源的当前更新次数t_Fi＝0；

第二步，生成待优化预测模型参数的初始解集θ：

$\mathrm{\&theta;}=[{\mathrm{\&theta;}}_1,{\mathrm{\&theta;}}_2,...,{\mathrm{\&theta;}}_{N_s}],{\mathrm{\&theta;}}_i=({\mathrm{\&theta;}}_{i,1},{\mathrm{\&theta;}}_{i,2})$

其中，i＝1，...，N_e，θ_i为第i组初始解；θ_i，j表示第i组解中的第j个值，对应为待确定预测模型中的一个参数，可通过下式产生θ_i，j的初始解：

θ_i，j＝θ_min+rand(0，1)*(θ_max-θ_min，j)

其中，i＝1，...，N_e，j＝1，2；rand(0，1)表示产生[0，1]之间的随机数，θ_{min， j}表示θ_i，j的最小值，θ_max，j表示θ_i，j的最大值；

利用下式计算所有解θ_i对应的适应度值，

${\mathrm{fitness}}_i\left({\mathrm{\&theta;}}_i\right)=\left\{\begin{array}{cc}\frac{1}{1+{\mathrm{fit}}_i\left({\mathrm{\&theta;}}_i\right)},& {\mathrm{fit}}_i\left({\mathrm{\&theta;}}_i\right)\&GreaterEqual;0\\ 1+\left|{\mathrm{fit}}_i\left({\mathrm{\&theta;}}_i\right)\right|,& {\mathrm{fit}}_i\left({\mathrm{\&theta;}}_i\right)<0\end{array}\right.$

其中，fit_i(θ_i)表示此处水质分解后分量数据拟合值与原始水质分量数据之间的误差，即：

${\mathrm{fit}}_i\left({\mathrm{\&theta;}}_i\right)={\mathrm{\&lambda;}}_1*\sqrt{\frac{j1}{1,\mathrm{tn}}\underset{i=1}{\overset{n}{\mathrm{\&Sigma;}}}{(x_{i,t}-{\hat{x}}_{i,t}^j)}^2}+{\mathrm{\&lambda;}}_2*\frac{\max \left\{\right|x_{1,t}-{\hat{x}}_{1,t}^j|,\&CenterDot;\&CenterDot;\&CenterDot;,|x_{n,t}-{\hat{x}}_{n,t}^j\left|\right\}}{\sqrt{\frac{1}{n}\underset{i=1}{\overset{n}{\mathrm{\&Sigma;}}}{(x_{i,t}-\stackrel{\&OverBar;}{x})}^2}}$

其中，x_i，t为水质分解后的第i个分量在t时刻的真实值，

为水质分解后的第j个模型第i个分量在t时刻的拟合值，

为分解后的第i个分量的平均值，λ₁和λ₂为调节因子，分别是0到1之间的随机数，n是水质分解后分量数据的个数；

第三步，按照下式计算各个解的适应度值finess_i(θ_i(t))，其中t是指第t次迭代，并将θ_i(t)的适应度值从小到大进行排序，适应度最大的θ_i(t)为侦查蜂，θpb_i(t)为每只采蜜蜂当前发现的最好解，表示为θpb_i(t)＝θ_i(t)，i＝1，2，...，N_e，并且θpb_i，j(t)＝θ_i，j(t)，其中θpb_i，j(t)为第i只采蜜蜂在第j维的最好解；θgb(t)设定为群体中的全局最好解，

并且 $f(\mathrm{\&theta;gb}\left(t\right)=\min \{f({\mathrm{\&theta;pb}}_1\left(t\right),...,f\left({\mathrm{\&theta;pb}}_{N_e}\left(e\right)\right))\};$

第四步，采蜜蜂按下式对食物源θ_i,j(t)进行邻域搜索产生新解θ_i,j(t+1)，并进行解的替代；

$\left\{\begin{array}{cc}{\mathrm{\&theta;}}_{i,j}(t+1)=w_{i,j}\left(t\right)*{\mathrm{\&theta;}}_{i,j}\left(t\right)+c_1*r_a*({\mathrm{\&theta;pb}}_{i,j}\left(t\right)-{\mathrm{\&theta;}}_{i,j}\left(t\right))+\mathrm{\&lambda;}\&CenterDot;N(0,1)+c_2*r_b*({\mathrm{\&theta;gb}}_j\left(t\right)-{\mathrm{\&theta;}}_{i,j}\left(t\right)),& f\left({\mathrm{\&theta;}}_{i,j}(t+1)\right)<f\left({\mathrm{\&theta;}}_{i,j}\left(t\right)\right)\\ {\mathrm{\&theta;}}_{i,j}(t+1)={\mathrm{\&theta;}}_{i,j}\left(t\right),& f\left({\mathrm{\&theta;}}_{i,j}(t+1)\right)>f\left({\mathrm{\&theta;}}_{i,j}\left(t\right)\right)\end{array}\right.$

$w_{i,j}=w_{\max }-\frac{(w_{\max }-w_{\min })t}{T_{\max }}$

其中，θ_i，j(t+1)为产生的新解，θ_i，j(t)为原来的旧解；θpb_i，j(t)为θpb_i(t)在第j维最好解；θbg_j(t)为θbg(t)的第j维全局最好解，(i＝1，2，...，N_e；j=1，2，...，D)，c₁和c₂为学习因子，N(0，1)为产生高斯分布的随机数，W_max和W_min分为为惯性权重w的上下限，λ为增量系数，r_a和r_b是[-1，1]之间的随机数，t指第t次迭代；

第五步，各跟随蜂按照下式计算每个食物源被选择的概率，采用轮盘赌选择策略选择一个食物源并按照下式进行探索，产生新解并进行解的替代；

$P_i=\frac{{\mathrm{fitness}}_i\left({\mathrm{\&theta;}}_i\right)}{\underset{k=1}{\overset{\mathrm{Ne}}{\mathrm{\&Sigma;}}}{\mathrm{fitness}}_k\left({\mathrm{\&theta;}}_k\right)}$

$\left\{\begin{array}{cc}{\mathrm{\&theta;}}_{i,j}(t+1)={\mathrm{\&theta;}}_{i,j}\left(t\right)+{\mathrm{\&eta;}}_{i,j}({\mathrm{\&theta;}}_{i,j}\left(t\right)-{\mathrm{\&theta;}}_{k,j}\left(t\right))+{\mathrm{\&xi;}}_{i,j}({\mathrm{\&theta;gb}}_j\left(t\right)-{\mathrm{\&theta;}}_{k,j}\left(t\right)),& f\left({\mathrm{\&theta;}}_{i,j}(t+1)\right)<f\left({\mathrm{\&theta;}}_{i,j}\left(t\right)\right)\\ {\mathrm{\&theta;}}_{i,j}(t+1)={\mathrm{\&theta;}}_{i,j}\left(t\right),& f\left({\mathrm{\&theta;}}_{i,j}(t+1)\right)>f\left({\mathrm{\&theta;}}_{i,j}\left(t\right)\right)\end{array}\right.$

${\mathrm{\&xi;}}_{i,j}=\frac{2t}{T_{\max }+t}$

其中，θ_i，j(t+1)为产生的新解，θ_i，j(t)为原来的旧解；η_ij是[-1，1]之间的随机数，θ_k为跟遗蜂选择的采蜜蜂，θ_k，j为采蜜蜂θ_k的第j为元素，并且i＝1，2，...，N_u；j＝1，2，...，D；t_Fi=t_Fi+1，ξ_i，j为自适应系数；

第六步，若t_Fi超过预设最大更新次数T_limit时，探索得到的θ_i仍未能改进的食物源适应度，则派出侦察蜂采用下式进行对其重新初始化，获得新解并进行解的替代；

其中，i＝1；j＝1，2，...，D，θ_i，j(t+1)为产生的新解，θ_i，j(t)为原来的旧解；

是列维概率分布产生的随机搜索步长，θpb_i，j(t)和

分别为解θ_i的第j维的最大值；

第七步，若迭代次数t小于最大迭代次数t_max，则返回第三步，否则结束并输出最优解θgb(t)，即得到模型正则化参数γ和核函数参数σ的最优解组合，以获得最佳的正则化参数γ和核函数参数σ，建立人工蜂群优化最小二乘支持向量回归机水质预测模型。

5.如权利要求1所述的基于多尺度分析的水产养殖水质短期组合预测方法，其特征在于，在步骤四中，对中频滤波部分，使用BP神经网络进行建模预测的具体步骤为：

第一步，接收中频滤波序列数据，网络初始化；

第二步，计算误差和梯度向量，进行模型的权值修正；

第三步，对精度和迭代次数是否满足算法终止条件作出判断，满足进行下一步，不满足返回第二步；

第四步，BP神经网络预测模型建立完成。

6.如权利要求1所述的基于多尺度分析的水产养殖水质短期组合预测方法，其特征在于，在步骤四中，对低频滤波部分，使用自回归滑动平均模型进行建模预测的具体方法为：

第一步，接收低频序列数据，进行模型的结构辨识；

第二步，进行模型参数估计、诊断和检验；

第三步，对模型是否合适作出判断，是则建立自回归滑动平均预测模型，否返回第一步。

Images