CN106872658A - 一种基于向量时间序列模型的污水cod负荷预测的方法 - Google Patents

Abstract

本发明公开了一种基于向量时间序列模型的污水COD负荷预测的方法，包括下列步骤：利用数据挖掘建模过程进行挖掘目标的定义；对获取的污水进水相关数据进行数据的探索和预处理以达到向量自回归模型建模的要求；选用最大似然估计模型的参数，同时选用信息准则确定模型的阶次，然后利用多元混成统计进行模型检验；进一步通过目标参数来简化模型，从而建立简洁有效的预测模型。所述采集的进水数据已确定进水量和进水COD，同时构造了污水COD负荷这一属性；利用得到的预测模型进行测试数据的预测，模型的输出则是污水COD负荷相关变量的预测结果。本发明的方法，模型简洁而且可以同时预测多个变量，预测精度比较高同时所需要的时间短。

Description

一种基于向量时间序列模型的污水COD负荷预测的方法

技术领域

本发明涉及污水处理厂进水的污水负荷预测的技术领域，具体涉及一种基于向量时间序列模型的污水COD负荷预测的方法。

背景技术

随着工业化程度的不断提高和人口的增长，城市污水排放量迅速增长，对环境的影响非常巨大,世界各国都投入大量的资金用于城市污水处理技术的研究,开发出了许多新技术和新工艺,对改善水环境状况起了巨大的作用。目前,建造污水处理厂是公认的解决水污染现状的一条有效途径,各国都在大力兴建。

但污水处理是一项能量密集型的综合技术,随着污水厂的大规模兴建,污染物的排放得到一定程度的控制,但是随之问题也浮现出来——污水处理设施的高投资、高运行成本在一定程度上阻碍了污水处理事业的发展,许多中小型城市由于资金问题没有修建污水处理厂,已建成的一些污水厂也因运行成本过高处于停产和半停产状态。纵观影响城市污水处理厂建设发展的因素,能耗问题是最主要的影响因素,如果不能及时有效的解决,污水处理的发展会受到很大的阻碍。我国相比于国外同种工艺类型的城市污水处理厂,能耗普遍偏高而且出水水质和污泥处理均达不到发达国家的处理水平。清华大学的钱易教授曾指出传统活性污泥法流程复杂、基建费用高、运行中需要很大的能源消耗、运行费用很可观,我国经济实力和发达国家相比较弱,不应该直接照搬发达国家的经验和技术来解决自己城市污水处理问题。针对这种现象,通过对处理工艺的深入研究,找出制约因素并提出相应的改造方案,节能的空间是非常大的,这是我国污水处理行业急需攻克的技术难题。有效降低污水处理系统的能耗、合理分配能源成为决定污水厂正常运行的关键因素。因此,研究和探索污水处理工艺的节能技术和途径具有非常重要的现实意义。

同时，污水的实际处理量与能耗的关系非常密切,随着污水处理量的增加,能耗也随之增加,但是能耗增幅的大小和该污水处理厂的设备的实际运行状况有着非常大的关系。随着污水处理量的不断增长以及出水排放标准的日趋严格，如何充分应对污水处理厂进水负荷和运行环境的变化，实现污染物的稳定达标并降低运行成本，成为污水处理运行管理所要面对的关键问题。自上世纪70年代以来，依赖于监测设备(硬件)及控制算法(软件)的发展，污水处理仪表、控制和自动化技术(Instrumentation Control andAutomation,ICA)在污水处理厂日常运行中发挥着越来越重要的作用。研究及实践表明，根据工艺的运行原理及特点，提高污水处理系统ICA技术的应用程度是加强污水处理厂脱氮除磷效果，实现节能降耗目标的有效手段。

对污水处理厂而言，进水流量及污染物浓度的变化是控制系统所要应对的扰动的主要来源。不同于其他工业生产过程，污水处理厂进水负荷一般不可调控，其突变性强，波动范围大，尤其在降雨的影响下，对生化处理系统具有较强的冲击性。

国内现有的污水处理中的相关模型技术，主要集中在污水处理后的水质COD、BOD等的预测和软测量中，而忽略了进水的重要性。实际上，对于污水处理的控制系统而言，出水的污染物浓度达到国家许可的排放标准是基本目标，进水的污染物负荷对于控制系统而言则是主要考虑的输入影响因素。

发明内容

本发明的目的是为了解决现有技术中的上述缺陷，提供一种基于向量时间序列模型的污水COD负荷预测的方法，该方法考虑污水进水的水量、水质COD等关联变量。

本发明的目的可以通过采取如下技术方案达到：

一种基于向量时间序列模型的污水COD负荷预测的方法，包括下列步骤：

S1、基于数据建模目标的变量选取：利用控制系统采集得到的污水处理厂进水的变量进行分析，所述变量包括进水量、污水COD、NH₄N、PH和进水温度T，通过定性分析采集变量与污水COD负荷的关联和影响程度，选取出对污染物负荷有影响的变量；

S2、对选取变量进行数据探索：对于步骤S1中选取的对污水COD负荷相关的变量进行数据结构的分析，包括：数据质量分析和数据特征分析；

S3、对探索后的数据进行数据预处理，达到向量自回归建模的要求：通过对步骤S2分析的数据进行数据序列的清洗、滤波数据的集成、数据的变换和新属性的构造、数据属性的规约、数据序列的平稳化和零均值化的处理，通过对建模变量的筛选进行数据的规约，进行规约筛选确定选取进水量、进水COD和污水COD负荷为模型的变量；

S4、针对满足建模要求的数据建立污水COD负荷的预测模型并进行模型的检验：利用经过数据预处理达到建模要求的数据采用向量时间序列建立污水COD负荷的VAR预测模型，并通过采用AIC信息准则对模型定阶，进而建立污水COD负荷的VAR预测模型，然后利用多元混成统计进行模型的有效性检验；

S5、对建立的污水COD负荷的预测模型简化并优化模型的参数：通过目标参数简化得到简化模型，同样采用最大似然估计得到模型的参数得到优化模型；

S6、利用准备的数据测试优化模型，通过分析预测结果进行模型的评价：将规约好的进水数据作为污水COD负荷的预测模型的输入，得到进水量、进水COD和污水COD负荷的预测结果，进行结果的分析并给出模型的相关评价。

进一步地，所述数据质量分析用于检查原始数据中是否存在脏数据以及无法进行直接分析的数据，所述脏数据包括缺失值、异常值以及含有特殊字符的数据；

所述数据特征分析通过绘制图表、计算指定特征量进行变量的分布分析来发现数据的分布特征和分布类型，通过统计量分析进行数据的统计描述，以及周期性分析。

进一步地，利用单变量的散点图和箱形图初步进行所述异常值的判断，同时也利用统计学的3σ法原则进行所述异常值的排查。

进一步地，所述数据序列的清洗具体为：针对所述脏数据的处理中少量数据的缺失采用单序列的局部回归拟合插值，针对所述脏数据的处理中时间段数据的缺失采用所采集的数据当前时间点的均值填补，针对所述脏数据的处理中含有特殊字符的数据变换为常数值，针对所述脏数据的处理中异常值采用全局均值进行替换；

所述滤波数据的集成具体为：选取卡尔曼滤波对数据序列进行滤波处理，接着将所述步骤S1中选取出对污染物负荷有影响的变量集成为一张数据变量表；

所述数据的变换和新属性的构造具体为：采用滤波后的进水量和污水COD的乘积得到污水COD负荷这个新属性，新属性的构造计算如下：

污水COD负荷＝单位时间进水量x进水COD (1)；

所述数据属性的规约具体为：通过对集成的数据进行变量的相关性分析，选取对污水COD负荷影响明显的变量进行属性的规约，剔除对污水COD负荷影响不显著的属性，最后确定模型采用的数据是进水量、进水COD和污水COD负荷，

其中，相关性分析中采用皮尔逊相关系数：

式中，cov(X,Y)表示两个边路间的互协方差，σ_X，σ_Y则分别代表自协方差，E代表的是各个变量的期望；

所述数据序列的平稳化和零均值化处理具体为：首先，采用单位根检验法判断序列的平稳性：考虑k维VAR(p)模型，假定表达式为：

式中，是一个p次矩阵多项式，其中，φ₀是一个k维常数向量，并且对于i＞0，φ_i是kxk维矩阵且φ_i≠0，Ι_k是kxk维的单位矩阵，Βⁱ是i阶次的后移算子，a_t是独立同分布随机向量序列，其均值为0，协方差矩阵∑_a为正定矩阵，

假定序列是平稳的，求期望有：

(Ι_k-φ₁-…-φ_p)μ＝[φ(1)]μ＝φ₀ (4)

其中，μ＝Ε(z_t)，因此有μ＝[φ(1)]^-1φ₀，则有|Ι_k-φΒ|＝|φ(Β)|，因此序列VAR(p)弱平稳的充分和必要条件是：行列式方程|Ι_k-φΒ|＝0的所有解的绝对值都大于1，也可以说解的模大于1或者所有的解在单位圆外，

如果序列经判断不符合平稳性，则需要进行差分、取对数的数据变换方法来达到平稳性的要求，然后求取序列的均值，得到零均值的平稳序列。

进一步地，所述步骤S4具体包括：

S41、模型参数估计方法的选用：采用ML估计方法估计VAR(p)模型的参数，对于给定的数据集{z₁,…,z_T}，VAR(p)模型的极大似然函数为

式中，是关于φ_i的估计值，T是时间序列的样本数，p是阶次，k是模型的维数；

S42、模型阶次的选择：利用不同的信息准则来判定适合的模型阶次，所述信息准包括AIC准则、BIC准则和HQ准则，在正态分布假设下，VAR(p)模型的AIC准则、BIC准则和HQ准则分别如下：

其中，则是上式讨论的Σ_a的ML估计；

S43、模型的多元混成统计检验：进行模型检验，对模型的准确性进行判断，其中，模型的准确性判断准则包括：

<1>、所有拟合参数都是统计显著的，即在给定的显著水平α下；

<2>、残差没有显著序列或横截面积相关；

<3>、不存在结构变化或反常值；

<4>、残差没有破坏数据序列的原始分布假设。

本发明相对于现有技术具有如下的优点及效果：

1、本发明采取的数据处理和建模的一系列步骤是一种常用且重要的数据挖掘建模的步骤，而向量自回归模型(VAR)又是多元时间序列模型中最常用的模型，结合数据挖掘建立时间序列的模型对拥有大量数据的污水处理厂来说，是非常适合采用的方法。

2、采用的向量自回归模型(VAR)，模型相对容易估计，可以用最小二乘(LS)方法、极大似然(ML)方法或者贝叶斯方法；同时VAR模型与广泛应用在多元统计分析中的多因变量多元线性回归模型相似，易于理解和采用。

3、不同于单序列的时间序列模型(ARMA等)，向量自回归模型不仅考虑了自身的影响，同时考虑了对其具有影响的变量，在当前具有大量数据的影响下对于建立更加精确的模型是相当有助益的。

4、本发明利用向量自回归模型(VAR)建立的污水COD负荷模型，在对进水COD负荷预测的同时，也同时预测了进水量和进水COD，不仅即时反映了污水进水的相关情况，也对于污水处理厂的控制系统的实施提供了便利，同时为污水处理厂进一步节能降耗提供了前提。

附图说明

图1是本发明公开的基于向量自回归模型VAR(p)的污水COD负荷预测的模型建立的流程图；

图2是经过初步数据预处理达到建模要求的进水COD负荷的一阶差分的时间序列趋势图；

图3是经过初步数据预处理达到建模要求的进水量的一阶差分的时间序列趋势图；

图4是经过初步数据预处理达到建模要求的进水COD的一阶差分的时间序列趋势图；

图5是建模序列VAR模型不同阶次的三种信息准则随阶次增大的变化趋势图；

图6是建模序列初始VAR(5)的Q_k(m)统计量的p值和m的散点图；

图7是建模序列简化VAR(5)的Q_k(m)统计量的p值和m的散点图；

图8是图1所述方法的进水COD负荷的一阶差分序列的预测残差的结果图；

图9是图1所述方法的进水量一阶差分序列的预测残差的结果图；

图10是图1所述方法的进水COD的一阶差分序列的预测残差的结果图。

具体实施方式

为使本发明实施例的目的、技术方案和优点更加清楚，下面将结合本发明实施例中的附图，对本发明实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例是本发明一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有作出创造性劳动前提下所获得的所有其他实施例，都属于本发明保护的范围。

实施例一

本实施例公开了一种基于向量自回归模型(VAR)的多元时间序列的预测方法，该方法是一种采用相互关联变量过去的趋势对当前和未来影响的判断方法，包含以下步骤：

S1、基于数据建模目标的变量选取：利用控制系统采集得到的污水处理厂进水的变量进行分析，所述变量包括进水量、污水COD、NH₄N、PH和进水温度T，通过定性分析采集变量与污染物负荷的关联和影响程度，选取出对污染物负荷有影响的变量；

该步骤基于“污水COD负荷预测”这一数据建模目标进行变量选取，利用控制系统采集得到的污水处理厂进水的变量进行初步分析，结合污水处理厂实际情况采集进水量、污水COD、NH₄N、PH和进水温度T等变量，通过定性分析采集变量与污水COD负荷的关联和影响程度，选取出上述对污水COD负荷有影响的变量。

通过从工艺的角度和案例选取的污水厂实际情况选取了关于污水的进水量、进水COD、NH₄N、PH和进水温度T等变量作为初步选取的对污水COD负荷影响的输入变量。

S2、对选取与污水COD负荷相关变量进行数据探索：对于步骤S1中选取的对污水COD(化学需氧量，Chemical Oxygen Demand)负荷相关的变量进一步进行数据结构的分析，主要包括：数据质量分析和数据特征分析；

其中，数据质量分析用于检查原始数据中是否存在脏数据(脏数据主要指的是缺失值、异常值以及含有特殊字符的数据等)，以及无法进行直接的相应分析的数据；数据特征分析通过绘制图表、计算某些特征量等手段进一步进行变量的分布分析来发现数据的分布特征和分布类型，统计量分析进行数据的统计描述，以及周期性分析等。

具体应用中，步骤S2对进水数据的探索包含以下的几个详细的步骤：

步骤S21、数据质量的分析，主要是查找序列中的脏数据：由于控制系统传感器的敏感性等各种干扰原因导致的少量数据的缺失，传感器的损坏或者信号中断等原因导致的时间段数据的缺失，部分含有单位等特殊字符的变量数据；

步骤S22、主要根据数据的图形进行初步的数据特征分析：利用单变量的散点图和箱形图初步进行异常值的判断，同时也利用统计学的3σ法原则进行异常值的排查。

S3、对探索后的数据进行数据预处理，达到向量自回归建模的要求：数据预处理一方面是要提高数据的质量，另一方面是可以更加贴切地适应所采用的挖掘模型。通过对步骤S2分析的数据进行数据序列的清洗、滤波数据的集成以及数据的变换和新属性的构造、数据属性的规约来得到需要的数据，主要针对脏数据进行缺失值的填补、异常值的替换，同时为了排除干扰进行滤波处理，为了满足采用模型的要求需要进行的数据序列的平稳化和零均值化的处理，以及通过对建模变量的筛选进行数据的规约，进行规约筛选确定选取进水量、进水COD和污水COD负荷为模型的变量；

该步骤具体包括：

步骤S31、数据序列的清洗：通过步骤S2对数据质量和数据特征的分析，针对性地进行数据的清洗。步骤S21中脏数据的处理——少量数据的缺失，这里采用单序列的局部回归拟合插值；时间段数据的缺失则采用所采集的数据当前时间点的均值填补；含有特殊字符的数据变换为常数值。步骤S22的异常值采用全局均值进行替换；

步骤S32、滤波数据的集成：为了排除测量的干扰需要对数据序列进行滤波处理，这里提供采用均值滤波、中值滤波或卡尔曼滤波等不同的滤波方法的效果对比，结合实际案例的工艺情况选取卡尔曼滤波为案例滤波方法，接着将步骤S1选取的变量集成为一张数据变量表；

S33、数据的变换和新属性的构造：考虑污水COD负荷这一变量无法进行数据的采集，需要进行属性的构造——采用滤波后的进水量和污水COD这2个属性的乘积得到污水COD负荷这个新的属性，新属性的构造计算：

污水COD负荷＝单位时间进水量x进水COD (1)

S34、数据属性的规约：通过对集成的数据进行变量的关联性分析——相关性分析，选取对污水COD负荷影响比较明显的变量进行属性的规约，剔除对污水COD负荷影响不够显著的属性，最后确定模型采用的数据是进水量、进水COD和污水COD负荷这三个属性。采用的是熟知的皮尔逊相关系数：

式中，cov(X,Y)表示两个边路间的互协方差，σ_X，σ_Y则分别代表自协方差，E代表的是各个变量的期望。

S35、数据序列的平稳化和零均值化处理：平稳序列一般满足下列两个条件——(a)、序列的均值是与时间无关的常数；(b)、序列的协方差矩阵是常数正定矩阵，即一个弱平稳序列的均值和协方差不依赖于时间T，从而构成弱平稳时间序列。这里采用常用的单位根检验法判断序列的平稳性：考虑k维VAR(p)模型，假定表达式为：

式中，是一个p次矩阵多项式。其中，φ₀是一个k维常数向量，并且对于i＞0，φ_i是kxk维矩阵且φ_i≠0，Ι_k是kxk维的单位矩阵，Βⁱ是i阶次的后移算子，a_t是独立同分布随机向量序列，其均值为0，协方差矩阵∑_a为正定矩阵。

假定序列是平稳的，求期望有：

(Ι_k-φ₁-…-φ_p)μ＝[φ(1)]μ＝φ₀ (4)

其中，与式(1)中同样的表达式代表一样的含义，而μ＝Ε(z_t)，因此有μ＝[φ(1)]^-1φ₀

则有|Ι_k-φΒ|＝|φ(Β)|，因此序列VAR(p)弱平稳的充分和必要条件是行列式方程|Ι_k-φΒ|＝0的所有解的绝对值都大于1，也可以说解的模大于1或者它们所有的解在单位圆外。

如果序列经判断不符合平稳性，则需要进行差分(通常是一阶差分或者二阶差分就满足平稳性)、取对数等数据变换来达到平稳性的要求，然后求取序列的均值，得到零均值的平稳序列。

S4、针对满足建模要求的数据建立污水COD负荷的预测模型并进行模型的检验：利用经过数据预处理达到建模要求的数据采用向量时间序列建立污水COD负荷的VAR预测模型，并通过采用AIC等信息准则对模型定阶，进而建立污水COD负荷的VAR预测模型，然后利用多元混成统计进行模型的检验；

该步骤具体包括：

S41、模型参数估计选用的方法：一个VAR(p)模型可通过最小二乘(LS)、极大似然(ML)或者贝叶斯估计方法来估计模型的参数。在多元正态假设下，z_t服从k维正态分布，VAR(p)模型的ML估计和LS估计近似相等，这里采用ML估计方法。

对于给定的数据集{z₁,…,z_T}，VAR(p)模型的极大似然函数为

式中，是关于φ_i的估计值，T是时间序列的样本数，p是阶次，k是模型的维数。

S42、模型阶次的选择：设定VAR模型并选择了模型参数的估计方法后，需要参考不同的方法来确定模型的阶次p。通常模型阶次是不断试探出来的，在不断试探的过程中取确定比较适合的模型阶次。通常采用的方法就是利用不同的信息准则来判定适合的模型阶次，信息准则在选择一个统计模型中是非常有效的，所有的信息准则都包含两部分(第一部分是关于数据模型的拟合优度，而第二部分是惩罚更复杂模型)。最著名的信息准则有3个(AIC准则、BIC准则和HQ准则)，在正态分布假设下，VAR(p)模型的3个准则：

其中的参数同以上各式的含义一样，则是上式讨论的Σ_a的ML估计。AIC是Akaike(1973)提出来的赤池信息准则，BIC代表Schwarz(1978)的贝叶斯信息准则，HQ是由Han和Quinn(1979)以及Quinn(1980)提出的。AIC用因子2来惩罚每个参数，而BIC和HQ采用的惩罚因子依赖于样本数——T越大，BIC和HQ更加严重地惩罚复杂模型。

S43、模型的多元混成统计检验：初步确定模型的阶次之后，需要进行模型检验，也就是众所周知的诊断检验或者说残差分析，这在模型建立中发挥着重要的作用。模型检验主要包括：(a)、确保拟合模型的准确性；(b)、如果有实际需要，通过分析给出模型进一步改进的方向——根据选择的信息准则来判断模型的准确性，会依赖分析结果，模型的准确性主要从以下几点来说明：

<1>、所有拟合参数都是统计显著的(在给定的显著水平α下)；

<2>、残差没有显著序列或横截面积相关；

<3>、不存在结构变化或反常值；

<4>、残差没有破坏数据序列的原始分布假设(诸如多元正态分布)等。

准确模型的残差应该是白噪声序列，因而检验残差连续性以及交叉相关性就成为模型检验的重要组成部分。

令R_p是a_t的理论之后p的交叉相关矩阵，模型检验的相关假设为：

H₀:R₁＝…＝R_m与H_a:R_j≠0，对于某些1≤j≤m (9)

其中，m是预先设定的正整数，进而得到残差序列的混成统计量：

它是服从自由度为(m-p)k²的卡方分布。

S5、对建立的模型简化并优化模型的参数：通过步骤S4的多元混成统计检验模型的有效性，然后通过目标参数简化得到的模型，同样采用最大似然估计得到模型的参数，从而进一步得到优化模型；

对模型的简化：事实上，经常观测到的有些参数在给定的显著水平下统计不显著，而剔除不显著参数有利于简化模型，尤其是当没有先验知识用于支持这些参数时，就需要剔除这些参数。给定特定的显著水平，例如α＝0.05，可以识别要剔除的目标参数。所谓目标参数，是指这些参数各自的T比率小于可允许的α的正态分布的临界值。

S6、利用准备的数据测试优化模型，通过分析预测结果进行模型的评价：将规约好的进水数据作为污水COD负荷的预测模型的输入，得到进水量、进水COD和污水COD负荷的预测结果，进行结果的分析并给出模型的相关评价。

该步骤中通过利用测试样本进行数据预测，给出预测结果的展示和分析，来判断模型的适用性。

得到较为准确的模型并进行模型检验之后，接着进行步骤S6用测试样本数据来进行模型测预测效果检测，并给出预测结果的展示和分析，最终得到确定模型适用性的结论。

上述方法中，主要考虑的是步骤S3中的数据预处理和步骤S4中的模型定阶两个方面：(a)、数据预处理不够好，影响模型的探索和模型的准确性；(b)、模型阶次判断的不够准确，则直接影响模型的稳定性和预测的效果。

实施例二

如图1，一种基于向量时间序列的污水COD负荷的预测方法，包含以下的建模和模型评价步骤：

1、通过从某污水处理厂的控制系统获取的污水进水量的数据库，里面包含2016年第二季度5月和6月的所有数据，结合污水处理的A2O工艺和数据库采集的数据变量，选取了与污水污染物负荷相关联的诸如进水量、进水COD、进水NH₄N、进水PH和进水水温T等变量；

2、接着进行数据探索的过程，首先查看数据的质量：

A、对于选取的5个变量，最先开始也是易于查看的就是数据的缺失值，通过对5月和6月数据的初步查看，发现5月的数据缺失严重，从而选择6月比较全的部分(2016年6月4日到2016年6月23日)数据作进一步的分析。

同时，短期COD测量的值比较稳定，从而考虑采集频率1h的数据，获得的原始数据经过初步的整理，特殊字符的加工后的结果由附表1可以查看；

B、附表1中的进水量代表的是安装测量表以来总的进水处理量，为了得到单位1h内的进水量，我们需要进行一阶差分获取单位时间的进水量这一属性；同时，对于各变量序列中局部少量缺失值进行局部线性拟合插补，局部缺失值连续量过多的部分，则采用单位时间点的数据的均值插补；由于无法直接采集污水COD负荷这一属性，因而我们需要采用式(1)来构造和获取单位时间污水COD负荷这一属性，从而获取了如附表2的集成的初步样本数据。

C、对于集成的样本数据z_t＝{l_t,f_t,c_t,n_t,p_t,t_t}(其中z_t表示向量时间序列，括号内则分别代表是污水COD负荷、进水量、COD、NH4N、PH和温度T的单时间序列)，利用式(2)的皮尔逊相关系数进行变量间的相关性分析，变量间的相关性如表1所示。

表1.污水COD负荷变量间的相关系数

相关系数	COD负荷	进水量	进水COD	进水NH4N	进水PH	进水水温T
							COD负荷	1	0.663	0.639	0.126	0.215	0.167
进水量	0.663	1	-0.151	0.275	0.234	0.039
							进水COD	0.639	-0.151	1	-0.113	0.052	0.185
进水NH4N	0.126	0.275	-0.113	1	0.450	-0.015
							进水PH	0.215	0.234	0.052	0.450	1	-0.193
进水水温T	0.167	0.039	0.185	-0.015	-0.193	1

D、由相关系数表进行属性规约，从而确定模型采用污水COD负荷、进水量、进水COD三个属性为模型的输入z_t＝{l_t,f_t,c_t}(t＝1,…,n)。

E、接着进行模型的平稳性检验，根据式(3)和式(4)利用单位根ADF检验法检验，向量时间序列不满足平稳性，而一阶差分z_t0'＝z_t-z_t-1(t＝2,…,n)的向量时间序列则满足向量时间序列的平稳性，进一步对平稳化的向量时间序列进行零均值化处理z_t'＝z_t0'-φ₀，从而得到符合建模要求的向量时间序列(按照3:1的比例，数据分为建模数据和测试数据)，满足建模要求的污水COD负荷向量时间序列的趋势图见附图2至图4。

3、设定向量时间序列的模型为VAR(p)模型，接着需要确定的是模型的阶次，对于向量时间序列，考虑公式(3)展开的VAR(p)模型：

z_t＝φ₀+φ₁z_t-1+…+φ_pz_t-p+a_t,t＝p+1,…,T (11)

设定不同的阶次p，然后利用公式(5)的ML估计计算模型的参数，求取不同阶次下的信息准则的大小，得到不同阶次下的信息准则如表2所示。

表2.不同阶次p对应的模型的信息准则表

p	AIC	BIC	HQ	P值
					0	-0.739	-0.739	-0.739	0
1	-1.770	-1.691	1.739	0
					2	-1.809	-1.653	-1.748	0.0001
3	-1.870	-1.635	-1.778	0
					4	-2.180	-1.867	-2.057	0
5	-2.300	-1.908	-2.146	0
					6	-2.305	-1.835	-2.121	0.0233
7	-2.304	-1.755	-2.088	0.0651
					8	-2.734	-1.747	-2.127	0
9	-2.385	-1.680	-2.108	0.0104
					10	-2.429	-1.645	-2.121	0.0001
11	-2.432	-1.570	-2.09	0.0410
					12	-2.489	-1.548	-2.119	0
13	-2.538	-1.519	-2.138	0

应用三种不同的信息准则，允许选取的最大阶为13：AIC准则最优的阶次p＝13，而BIC和HQ准则选取的最优阶次p＝5，表2总结了这些统计量。这表明，对于多元时间序列，不同的准则将选取不同的阶次，然而这些统计量是估计的，所以某种程度上的取值就不是太严格。附图5给出了3个信息准则的时序图，对于p∈{1,2,3,4,5}，AIC一直在减小，而BIC和HQ的趋势相同，在p＝5显示最小，之后的阶次BIC和HQ的p值不断增大。而观察后续的阶次，发现AIC随着阶次的增大不断较小，但是并没有很明显的优势，同时考虑到阶次p越高，待估计的参数越多，参数估计的难度也随之增大，不确定增加，从而VAR(5)将作为三维时间序列的开始模型。

利用ML对参数进行估计，得到VAR(5)模型的参数矩阵：

从而可以得到初始VAR(5)模型的表达式：

z_t＝φ₀+φ₁z_t-1+φ₂z_t-2+φ₃z_t-3+φ₄z_t-4+φ₅z_t-5+a_t (12)

其中，残差协方差矩阵为

初步模型的AIC＝-2.146,BIC＝-1.754,HQ＝-1.991。同时，得到初步的模型之后，需要进行模型检验：利用多元混成统计量来检验残差交叉相关性，附图6详细给出了统计量的Q₅(m)的p值。因为有48个参数，所以Q_k(m)卡方分布的自由度为9m-48。因此，如果m≥5，近似卡方分布成立。但是根据附图4，Q₅(m)统计量并非如此，反而在m≥6显示残差没有强序列和交叉相关性。说明残差检验并不合格，需要进一步进行模型的参数优化，也表明求取的模型确实存在某些参数在给定的α＝0.05的显著水平上是统计不显著的的。

根据卡方统计量，用α＝0.05来识别目标参数，得到简化模型的参数矩阵：

此时的残差协方差矩阵为

查看简化后的模型，发现简化后的模型参数由原来的48个减少到37，大大降低了模型估计的难度。同时简化后的模型中AIC＝-1.187,BIC＝-0.865,HQ＝-1.060。对于本案例，简化模型的所有3个信息准则反而有所增加的值。而此时污水COD负荷的简化模型VAR(5)为

再次利用多元混成统计来检验残差交叉相关性。如附图7详细给出了统计量Q₅(m)的p值，观察此时的图形发现，Q_k(m)统计量表明VAR(5)模型已经不存在强序列或交叉相关性了。

从而，上述污水COD负荷简化模型VAR(5)对于负荷的预测是合适的。4.利用得到的模型对测试的数据进行预测，同时检验模型的预测能力，从而进一步对模型进行评价。

利用测试数据得到预测结果。利用简化后的VAR(5)模型得到的预测残差很小，残差效果如附图8至图10所示，模型预测的精度比较高，而且同时也达到了预测多个变量的便利和要求。

从整体上来看，建立的模型能够使用污水COD负荷的要去，而且依赖项的进水量和进水COD通过在线传感器的测量也比较方便，可以说模型是很适用的。

上述实施例为本发明较佳的实施方式，但本发明的实施方式并不受上述实施例的限制，其他的任何未背离本发明的精神实质与原理下所作的改变、修饰、替代、组合、简化，均应为等效的置换方式，都包含在本发明的保护范围之内。

Claims (5)

1.一种基于向量时间序列模型的污水COD负荷预测的方法，其特征在于，包括下列步骤：

S1、基于数据建模目标的变量选取：利用控制系统采集得到的污水处理厂进水的变量进行分析，所述变量包括进水量、污水COD、NH₄N、PH和进水温度T，通过定性分析采集变量与污水COD负荷的关联和影响程度，选取出对污染物负荷有影响的变量；

S2、对选取变量进行数据探索：对于步骤S1中选取的对污水COD负荷相关的变量进行数据结构的分析，包括：数据质量分析和数据特征分析；

S3、对探索后的数据进行数据预处理，达到向量自回归建模的要求：通过对步骤S2分析的数据进行数据序列的清洗、滤波数据的集成、数据的变换和新属性的构造、数据属性的规约、数据序列的平稳化和零均值化的处理，通过对建模变量的筛选进行数据的规约，进行规约筛选确定选取进水量、进水COD和污水COD负荷为模型的变量；

S4、针对满足建模要求的数据建立污水COD负荷的预测模型并进行模型的检验：利用经过数据预处理达到建模要求的数据采用向量时间序列建立污水COD负荷的VAR预测模型，并通过采用AIC信息准则对模型定阶，进而建立污水COD负荷的VAR预测模型，然后利用多元混成统计进行模型的有效性检验；

S5、对建立的污水COD负荷的预测模型简化并优化模型的参数：通过目标参数简化得到简化模型，同样采用最大似然估计得到模型的参数得到优化模型；

S6、利用准备的数据测试优化模型，通过分析预测结果进行模型的评价：将规约好的进水数据作为污水COD负荷的预测模型的输入，得到进水量、进水COD和污水COD负荷的预测结果，进行结果的分析并给出模型的相关评价。

2.根据权利要求1所述的一种基于向量时间序列模型的污水COD负荷预测的方法，其特征在于，

所述数据质量分析用于检查原始数据中是否存在脏数据以及无法进行直接分析的数据，所述脏数据包括缺失值、异常值以及含有特殊字符的数据；

所述数据特征分析通过绘制图表、计算指定特征量进行变量的分布分析来发现数据的分布特征和分布类型，通过统计量分析进行数据的统计描述，以及周期性分析。

3.根据权利要求2所述的一种基于向量时间序列模型的污水COD负荷预测的方法，其特征在于，

利用单变量的散点图和箱形图初步进行所述异常值的判断，同时也利用统计学的3σ法原则进行所述异常值的排查。

4.根据权利要求2所述的一种基于向量时间序列模型的污水COD负荷预测的方法，其特征在于，

所述数据序列的清洗具体为：针对所述脏数据的处理中少量数据的缺失采用单序列的局部回归拟合插值，针对所述脏数据的处理中时间段数据的缺失采用所采集的数据当前时间点的均值填补，针对所述脏数据的处理中含有特殊字符的数据变换为常数值，针对所述脏数据的处理中异常值采用全局均值进行替换；

所述滤波数据的集成具体为：选取卡尔曼滤波对数据序列进行滤波处理，接着将所述步骤S1中选取出对污染物负荷有影响的变量集成为一张数据变量表；

所述数据的变换和新属性的构造具体为：采用滤波后的进水量和污水COD的乘积得到污水COD负荷这个新属性，新属性的构造计算如下：

污水COD负荷＝单位时间进水量x进水COD (1)；

所述数据属性的规约具体为：通过对集成的数据进行变量的相关性分析，选取对污水COD负荷影响明显的变量进行属性的规约，剔除对污水COD负荷影响不显著的属性，最后确定模型采用的数据是进水量、进水COD和污水COD负荷，

其中，相关性分析中采用皮尔逊相关系数：

${\mathrm{\&rho;}}_{X,Y}=\frac{\mathrm{cov}(X,Y)}{{\mathrm{\&sigma;}}_X{\mathrm{\&sigma;}}_Y}=\frac{E\left(\right(X-{\mathrm{\&mu;}}_X\left)\right(Y-{\mathrm{\&mu;}}_Y\left)\right)}{{\mathrm{\&sigma;}}_X{\mathrm{\&sigma;}}_Y}=\frac{E\left(XY\right)-E\left(X\right)E\left(Y\right)}{\sqrt{E\left(X^2\right)-E^2\left(X\right)}\sqrt{E\left(Y^2\right)-E^2\left(Y\right)}}---\left(2\right)$

式中，cov(X,Y)表示两个边路间的互协方差，σ_X，σ_Y则分别代表自协方差，E代表的是各个变量的期望；

所述数据序列的平稳化和零均值化处理具体为：首先，采用单位根检验法判断序列的平稳性：考虑k维VAR(p)模型，假定表达式为：

式中，是一个p次矩阵多项式，其中，φ₀是一个k维常数向量，并且对于i＞0，φ_i是kxk维矩阵且φ_i≠0，Ι_k是kxk维的单位矩阵，Βⁱ是i阶次的后移算子，a_t是独立同分布随机向量序列，其均值为0，协方差矩阵∑_a为正定矩阵，

假定序列是平稳的，求期望有：

(Ι_k-φ₁-…-φ_p)μ＝[φ(1)]μ＝φ₀ (4)

其中，μ＝Ε(z_t)，因此有μ＝[φ(1)]^-1φ₀，则有|Ι_k-φΒ|＝|φ(Β)|，因此序列VAR(p)弱平稳的充分和必要条件是：行列式方程|Ι_k-φΒ|＝0的所有解的绝对值都大于1，也可以说解的模大于1或者所有的解在单位圆外，

如果序列经判断不符合平稳性，则需要进行差分、取对数的数据变换方法来达到平稳性的要求，然后求取序列的均值，得到零均值的平稳序列。

5.根据权利要求1所述的一种基于向量时间序列模型的污水COD负荷预测的方法，其特征在于，所述步骤S4具体包括：

S41、模型参数估计方法的选用：采用ML估计方法估计VAR(p)模型的参数，对于给定的数据集{z₁,…,z_T}，VAR(p)模型的极大似然函数为

$L(\widehat{\mathrm{\&beta;}},\widehat{{\&Sigma;}_a}|z_{1:p})={\left(2\mathrm{\&pi;}\right)}^{-k(T-p)/2}{\left|\widehat{{\&Sigma;}_a}\right|}^{-(T-p)/2}\exp \&lsqb;\frac{-k(T-p)}{2}\&rsqb;---\left(5\right),$

式中，是关于φ_i的估计值，T是时间序列的样本数，p是阶次，k是模型的维数；

S42、模型阶次的选择：利用不同的信息准则来判定适合的模型阶次，所述信息准包括AIC准则、BIC准则和HQ准则，在正态分布假设下，VAR(p)模型的AIC准则、BIC准则和HQ准则分别如下：

$AIC\left(p\right)=\ln \left|\widehat{{\&Sigma;}_{a,p}}\right|+\frac{2}{T}{\mathrm{pk}}^2---\left(6\right)$

$BIC\left(p\right)=\ln \left|\widehat{{\&Sigma;}_{a,p}}\right|+\frac{\ln \left(T\right)}{T}{\mathrm{pk}}^2---\left(7\right)$

$HQ\left(p\right)=\ln \left|\widehat{{\&Sigma;}_{a,p}}\right|+\frac{2\ln \&lsqb;\ln \left(T\right)\&rsqb;}{T}{\mathrm{pk}}^2---\left(8\right)$

其中，则是上式讨论的Σ_a的ML估计；

S43、模型的多元混成统计检验：进行模型检验，对模型的准确性进行判断，其中，模型的准确性判断准则包括：

<1>、所有拟合参数都是统计显著的，即在给定的显著水平α下；

<2>、残差没有显著序列或横截面积相关；

<3>、不存在结构变化或反常值；

<4>、残差没有破坏数据序列的原始分布假设。