一种基于函数型主成分分析和小波变换的变压器运行数据缺失点修补方法与流程

本发明涉及电力设备数据清洗领域，尤其是一种基于函数型主成分分析和小波变换的变压器运行数据缺失点修补方法。

背景技术：

变压器是电网中最重要的设备之一，其运行数据对后续的大数据分析至关重要。然而在数据采集和运输过程中，可能由于一些故障和人为因素导致数据缺失，这不利于后续的数据分析和数据挖掘，因此有必要对缺失数据进行填补。

目前常用的缺失数据填补方法有：人工填充、插值、回归等。但是人工填充费时费力，要对每日产生的巨量数据进行人工填充也是不切实际的；插值方法的可信度不高，且插值方法难以适应大量连续数据缺失的情况；回归方法对函数模型的准确性要求很高，但是要确定一个合理的函数模型难度很大，而且模型的适用性很小。因此有必要建立一种新的方法，通过数据内在的特征，对缺失数据进行填补。

技术实现要素：

本发明的目的是提供一种基于函数型主成分分析和小波变换的变压器运行数据缺失点修补方法，只需对已知的数据点进行学习，就可以找到数据内在的特征，拟合出整个时间序列上的函数，通过这个函数在对应时间点上的值就可以对缺失点的值做出修补。

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

一种基于函数型主成分分析和小波变换的变压器运行数据缺失点修补方法，包括以下步骤：

利用fpca方法对已采集的运行数据点进行分析，拟合出整个时间序列上的运行数据函数xi(t)；

对原数据点和通过fpca得到的数据点做差，得到残差函数ε(t)；

利用小波对残差函数ε(t)变换去噪，得到ε′(t)；

估计函数将缺失点处的时间t0代入得到的作为缺失点处的修补值。

进一步地，所述利用fpca方法对已采集的运行数据点进行分析，拟合出整个时间序列上的运行数据函数xi(t)步骤之前，还包括：

取该变压器n天内所有有效采样时间和对应的有效采样值，得到一个时间记录矩阵t和采样值记录矩阵y，t的第i行第j列元素tij表示第i天第j个采样点的时间，y的第i行第j列元素yij表示第i天第j个采样点的数值。

进一步地，所述利用fpca方法对已采集的运行数据点进行分析，拟合出整个时间序列上的运行数据函数xi(t)，具体包括：

将采集的变压器运行数据按时间顺序排列，十天数据作为一组进行分析；

根据kl定理，变压器每天的运行数据可以表示为，xi(t)表示第i天变压器的运行数据关于时间的函数，μ(t)是变压器十天里样本下的均值函数，αik是一组依赖于天数的系数，是十天样本下的一组基函数，估计均值函数μ(t)，估计基函数估计系数αik，选取前k个特征函数来表示xi(t)，得到xi(t)的初步估计，即

进一步地，所述估计均值函数μ(t)，具体包括：利用加权最小二乘法估计均值函数μ(t)：

优化函数分别对β0和β1求偏导，ni表示第i天的有效采样点数目，k(·)采用高斯核函数，hμ是带宽，通过gcv的方法选择带宽；

令两个偏导数为零，得到其中对μ(t)的估计为：

进一步地，所述估计基函数具体包括：

估计样本协方差函数g(s,t)；

根据通过离散化估计或是数值估计的方法得到特征函数和特征值λk，其中，γ表示xi(t)的定义域，是样本协方差函数g(s,t)的估计函数，g(s,t)的表达式为：g(s,t)＝cov(x(s),x(t))。

进一步地，所述估计样本协方差函数g(s,t)具体包括：

a.在s≠t处，通过加权最小二乘法得到g(s,t)估计值，即优化函数：分别对β0、β11和β12求偏导数，并令三个偏导数等于零，得到，

其中

这里k2(·,·)是一个二元核函数，其表达式为hg为带宽，通过gcv的方法选择带宽，

b.在s＝t处，先对坐标轴顺时针旋转45°，即有

再最小化如下的优化函数：

得到

进一步地，所述确定系数αik具体包括：

通过计算条件期望来估计系数αik：

其中

这里δjl在j＝l时为1，j≠l时为0，由下式给出

其中是由步骤a中得到的估计函数在对角线上的值，是由步骤b中得到的估计函数在对角线上的值。

进一步地，所述选取前k个特征函数来表示xi(t)，具体包括：

选择交叉验证的方式来选择k，即通过使下式的值最小：

其中表示将第i天的数据全部移除后，得到的第i天的估计函数，即

其中，表示将第1天的数据全部移除后，得到的第1天的估计函数，

表示将第1天的数据全部移除后，得到的第1天的特征值，

表示将第1天的数据全部移除后，得到的第1天的特征函数。

进一步地，所述对原数据点和通过fpca得到的数据点做差，得到残差函数ε(t)，具体包括：

将变压器第i天tj时刻的测量数据yij与做差得到残差函数εi(t)在tj时刻的值

根据业务逻辑，可以将缺失点分为两类，一类是小块的缺失，即连续缺失数据的时间间隔小于等于十分钟；一类是大块的缺失，即连续缺失数据的时间间隔大于十分钟，对残差函数缺失点的处理方法如下：

a.对小块的缺失，用移动平均的方式预测缺失点残差函数的值:如果第k分钟的点缺失，则前十分钟的缺失无法用这个公式进行计算，因此在边界处，对的估计为：

b.对大块的缺失，将大块的缺失处的残差函数的值设为零。

进一步地，所述利用小波对残差函数ε(t)变换去噪，得到ε′(t)，具体包括：

用matlab对残差函数ε(t)做小波变换，去除噪声：小波类型取db4小波，分解层数为3，阈值方法选择固定阈值估计法，噪声结构选择unscaledwhitenoise，得到去噪后的残差函数εi′(t)。

发明内容中提供的效果仅仅是实施例的效果，而不是发明所有的全部效果，上述技术方案中的一个技术方案具有如下优点或有益效果：

1、本发明无需预先根据数据的业务逻辑建立函数模型，只要提供足够多的测量点，就能根据数据的内在特征拟合出函数曲线。

2、对缺失点的分布没有要求，既可以对单个的缺失值做出预测，也可以对连续的缺失值做出预测。

3、本发明采用fpca和小波变换相结合的方式，既能在整体上把握数据集的特征，又能提高局部的拟合度，较传统的拟合手段，本发明的预测值可信度更高。

具体实施方式

为能清楚说明本方案的技术特点，下面通过具体实施方式，对本发明进行详细阐述。此外，本发明可以在不同例子中重复参考数字和/或字母。这种重复是为了简化和清楚的目的，其本身不指示所讨论各种实施例和/或设置之间的关系。本发明省略了对公知组件和处理技术及工艺的描述以避免不必要地限制本发明。

基于函数型主成分分析和小波变换的变压器运行数据缺失点修补方法，具体步骤如下：

1)采集到的数据可能同时包含多个变压器的信息，排列也可能是杂乱无序的。因此首先需要对这些数据进行整理。首先选取其中某个变压器的所有数据，然后按照时间顺序排列。此时该变压器排列好的数据一般情况下将包含12类不同的数据，它们分别是高压侧-u，中压侧-u，低压侧-u，高压侧-i，中压侧-i，低压侧-i，高压侧-p，中压侧-p，低压侧-p，高压侧-q，中压侧-q，低压侧-q。选取这12类数据中的任意一种数据进行分析，以高压侧-u为例，取该变压器n天内高压侧-u的所有有效采样时间和对应的有效采样值，得到一个时间记录矩阵t和采样值记录矩阵y。t的第i行第j列元素tij表示第i天第j个采样点的时间，y的第i行第j列元素yij表示第i天第j个采样点的数值。

2)用fpca方法(functionalprincipalcomponentanalysis)和kl定理(karhunen–loèvetheorem)对变压器每天的高压侧-u数据建模。

假设变压器每天的高压侧-u数据可以表示成函数xi(t)，那么根据kl表示，函数可以表示为：

这里xi(t)表示第i天变压器的高压侧-u数据关于时间的函数，μ(t)是变压器n天里所有高压侧-u数据的均值函数，αik是一组依赖于天数的系数，是n天样本下的一组基函数。考虑到测量时可能存在一定误差，可以认为真实的变压器数据情况为:

其中ε是均值为0，方差为σ²的随机误差项。每天的高压侧-u数据对应着yi(t)上不同时间点的采样：

3)用加权最小二乘的方法估计均值函数μ(t)。利用优化函数：

其中k(·)是一个核函数，常见的核函数有线性核函数，多项式核函数等。本发明采用高斯核函数，其表达式为：

hμ是带宽，带宽的选择会影响到均值函数估计的好坏，我们这里采用广义交叉检验的方法(generalizedcross-validation,gcv)选择带宽，一般也可根据经验决定带宽大小。

ni表示第i天的有效采样点数目，tij和yij的含义同第一点中所述，β0和β1是两个关于t的函数，则我们对μ(t)的估计为：

具体的计算方法为，将优化函数分别对β0和β1求偏导，令两个偏导数为零，得到

其中

4)然后我们需要估计函数基根据karhunen–loève定理，是样本协方差函数的特征函数，即满足

其中γ表示xi(t)的定义域，在我们这里的叙述中γ为[1,1440]。是样本协方差函数g(s,t)的估计函数，g(s,t)的表达式为：

g(s,t)＝cov(x(s),x(t))

因此在估计函数基之前，我们首先需要估计样本协方差函数g(s,t)。

a.在s≠t处，g(s,t)的估计是通过加权最小二乘得到，即最小化下面的优化函数：

这里k2(·,·)是一个二元核函数，其表达式为hg为带宽，仍可以通过gcv的方法选择带宽，据此，我们得到估计具体的计算方法同(3)，仍可以通过分别对β0、β11和β12求偏导数，并令三个偏导数等于零得到。即有，

其中

b.在s＝t处，由于在垂直对角线方向上局部二次拟合比局部一次拟合更接近协方差曲面的形状，因此我们需要修改一下优化函数。具体的，我们先对坐标轴顺时针旋转45°，即有

我们再最小化如下的优化函数：

得到

这样我们就得到了样本协方差函数g(s,t)的估计

5)通过得到的估计特征函数基我们根据以下等式得到特征函数和特征值λk：

其中γ的含义同(4)中所述。具体的计算我们可以通过离散化估计或是数值估计的方法求得。

6)确定系数αik。

我们通过计算条件期望来估计系数αik：

其中

这里δjl在j＝l时为1，j≠l时为0.我们还需要对做出估计，由下式给出

其中是由(4).a中得到的估计函数在对角线上的值，是由(4).b中得到的估计函数在对角线上的值。

7)选取前k个特征函数来表示xi(t)，我们需要确定k的数量。可以选择交叉验证的方式来选择k，具体的，使下式的值最小：

其中表示将第i天的数据全部移除后，再重复(1)-(6)步得到的第i天的估计函数，即

其中的含义和的含义相同。

结合(1)-(7)，我们就通过fpca得到了xi(t)的初步估计，即

接下来我们将用小波变换对xi(t)的估计做进一步的完善。

8)估计残差函数。

将变压器第i天tj时刻的测量数据yij与做差得到残差函数εi(t)在tj时刻的值但由于变压器的测量数据有缺失，因此不能得到残差函数εi(t)在所有tj时刻的值。为了方便后面进行小波变换去噪，需要对残差函数进行缺失点的补足。对残差函数缺失点的处理方法如下：

a.根据业务逻辑，可以将缺失点分为两类，一类是小块的缺失，即连续缺失数据的时间间隔小于等于十分钟；一类是大块的缺失，即连续缺失数据的时间间隔大于十分钟。

b.对小块的缺失，用移动平均的方式预测缺失点残差函数的值。具体的，如果第k分钟的点缺失，则要注意的是，前十分钟的缺失无法用这个公式进行计算，因此在边界处，我们对的估计为：

c.对大块的缺失，由于残差函数的绝对值本身就比较小，且函数值受到噪声的干扰，因此对残差函数大段的缺失进行预测意义不大，反而会使结果更加不可信。因此将大块的缺失处的残差函数的值设为零。

9)用matlab对残差函数做小波变换，去除噪声。

根据正交试验结果，选择最优的参数组合为：小波类型取db4小波，分解层数为3，阈值方法选择固定阈值估计法，噪声结构选择unscaledwhitenoise，得到去噪后的残差函数ε′i(t)。

10)将去噪后的残差函数ε′i(t)加到上，得到最终的估计函数，即

综上，我们得到最终的预测，即第i天第k分钟的变压器高压侧-u数据预测值为用该预测值进行缺失点的数值填补。

上述虽然对本发明的具体实施方式进行了描述，但并非对本发明保护范围的限制，所属领域技术人员应该明白，在本发明的技术方案的基础上，本领域技术人员不需要付出创造性劳动即可做出的各种修改或变形仍在本发明的保护范围以内。