CN116933028A - 基于b样条包络分析的遥测数据异常检测方法 - Google Patents

Abstract

本发明公开了基于B样条包络分析的遥测数据异常检测方法，具体为：将遥测数据{x₁，x₂，…，x_M}分为k+1个片段；求取每个数据片段上最大值P_max,i与最小值P_min,i，构成数据点集P_max，P_min；将每个片段的最大值P_max,i与最小值P_min,i进行移动，得到新的数据点集P_max'，P_min'；将数据点集P_max'，P_min'中的点进行B样条曲线拟合，将P_max'的曲线作为上包络线，将P_min'的曲线作为下包络线，将上包络线和下包络线作为判断数据是否异常的动态阈值，完成数据的异常检测。使用B样条拟合遥测数据的自适应包络作为阈值，有助于提高信号的质量和准确性，减少噪声的干扰。使用自适应包络能够更好地适应信号的动态变化，可以提高故障检测的灵敏度和准确性，对降低航天器故障损失和保障航天任务成功具有重大意义。

Description

基于B样条包络分析的遥测数据异常检测方法

技术领域

本发明属于航空航天数据处理技术领域，具体涉及基于B样条包络分析的遥测数据异常检测方法。

背景技术

在航天领域中，遥测是一项关键技术，能够向地面基站传回很多关键设备的遥测数据。近些年来，随着航空探测的需要，外太空航天器的系统功能快速增加，设备整体结构日趋复杂，使得设备发生故障的可能越来越大，这些故障的产生往往会导致巨大的经济损失和航天任务失败，因此通过对关键设备遥测数据的分析，达到设备故障排查的目的，对航天任务、科学研究等方面有重大意义。

在实际中，使用固定阈值的故障诊断方法通常不需要复杂的计算或算法，更容易实施。阈值的选择可以基于经验或专家知识，使得诊断结果能贴合实际情况并且该方法适用于许多故障诊断场景，适用性广泛。然而，固定阈值方法缺乏适应性，其是根据经验或理论设定的固定数值，无法自动适应信号的动态变化，并且容易受到主观因素的影响，可能导致漏检或误报。当系统出现未知的故障或变化时，固定阈值可能无法准确地检测和识别这些故障。

理论与实践均已表明，包络线分析是工程信号处理中行之有效的方法，在机械故障诊断和振动机械信号分析中发挥了重要作用。自适应包络能够根据实时信号的变化自动调整阈值，更好地适应信号的变化特性，并且可以根据故障信号的强度来调整阈值，提高故障的检测率。

发明内容

本发明的目的是提供基于B样条包络分析的遥测数据异常检测方法，能够准确识别出航天器遥测数据中的异常数据，从而实现航天器的故障诊断。

本发明所采用的技术方案是，基于B样条包络分析的遥测数据异常检测方法，具体按照以下步骤实施：

步骤1、将遥测数据{x₁，x₂，…，x_M}分为k+1个片段；

步骤2、求取每个数据片段上最大值P_max,i与最小值P_min,i，构成数据点集P_max，P_min；

步骤3、将每个片段的最大值P_max,i与最小值P_min,i进行移动，得到新的数据点集P_max'，P_min'；

步骤4、将数据点集P_max'，P_min'中的点进行B样条曲线拟合，将P_max'的曲线作为上包络线，将P_min′的曲线作为下包络线，将上包络线和下包络线作为判断数据是否异常的动态阈值，完成数据的异常检测。

本发明的特点还在于，

步骤1中，遥测数据为{x₁，x₂，…，x_M}，M为遥测数据长度，将遥测数据分为k个长度为d的数据片段，剩余长度不足d的视为第k+1个数据片段；这样，整个数据被分为k+1个片段，即{x₁₁，x₁₂，…，x_1d}，{x₂₁，x₂₂，…，x_2d}，…，{x_k1，x_k2，…，x_kd}，{x_(k+1)1，x_(k+1)2，…，x_(k+1)m}。

步骤2中，将遥测数据{x₁，x₂，…，x_M}划分为k+1个片段，找出每个数据片段的最大值P_max，i与最小值P_min，i(i＝1，2，…，k+1)，并将这些点作为数据点集P_max＝{P_max，1，P_max，2，…，P_max，k+1}，P_min＝{P_min，1，P_min，2，…，P_min，k+1}。

步骤3中，具体为：步骤1中的前k个长度为d的数据片段为{x₁₁，x₁₂，…，x_1d}，{x₂₁，x₂₂，…，x_2d}，…，{x_k1，x_k2，…，x_kd}，第i(i＝1，2，…，k)个片段的均值u_i、标准差σ_i计算公式如下：

最后一个包含m个剩余数据的数据片段为{x_(k+1)1，x_(k+1)2，…，x_(k+1)m}，其均值μ_k+1、标准差σ_k+1，计算如下：

将每个区间的最大值P_max，i，最小值P_min，i分别向上、向下移动标准差σ_i，得到P_max′＝{P_max，1′，P_max，2'，…，P_max，k+1′}，P_max′＝{P_max，1′，P_max，2′，…，P_max，k+1′}，作为新的数据点集。

步骤4中，具体为：

步骤4.1，数据点集P_max′中的点可表示为：(t_max，j，P_max，j′)，j＝1，2，…，k+1，t_max,j为P_max,j′所对应的时间，采用均匀B样条曲线进行拟合；

步骤4.2，数据点集P_min'中的点可表示为：(t_min，j，P_min，j′)，j＝1，2，…，k+1，t_min,j为P_min,j′所对应的时间，采用均匀B样条曲线进行拟合；

步骤4.3，由P_max'数据点集得到的B样条拟合曲线作为数据的上包络线，P_min'数据点集得到的B样条拟合曲线作为数据的下包络线，将数据的上包络线和下包络线作为判断数据是否异常的“动态阈值”，对于任意一组遥测数据，如果数据在阈值范围内，则判定为正常数据，若不在阈值范围内，则判定为异常数据。

步骤4.1中，拟合过程为：

步骤4.11，将T_max,0、T_max,,pmax作为数据点集P_max'的两个时间端点，T_max,0＝t_max,1，T_max,pmax＝t_max,k+1，并计算T_max,i，T_max，i＝T_max，0+ih_max，i＝0，1，…，p_max，其中p_max为拟合分段数，h_max为每个分段的时间长度，

步骤4.12，将(t_max,j,P_max,j′),j=1,2,…,k+1和T_max,i代入三次B样条方程中，得到矩阵T_max，如下式所示；

步骤4.13，将T_max代入B样条最优拟合逼近方程，如下式所示；

T_max ^TT_maxα_max＝T_max ^TX_max

其中，X_max为(t_max,j,P_max,j′),j=1,2,…,k+1组成的矩阵，α_max为需要求解的系数矩阵：

其中，

得到方程的解为α_max＝(T_max ^TT_max)^-1T_max ^TX_max；

步骤4.14，计算得到数据点集P_max'的B样条拟合曲线，其中α_i为系数矩阵α_max的第i个元素，t为数据点集P_max'所在的时间区间的任意值。

步骤4.2中，拟合过程为：

步骤4.21，将T_min,0、作为数据点集P_min'的两个时间端点，即T_min,0＝t_min,1，并计算T_min,i，T_min,i=T_min,0+ih_min,i=0,1,…,p_min，其中p_min为拟合分段数，h_min为每个分段的时间长度，

步骤4.22，将(t_min,j,P_min,j),j=1,2,…,k+1和T_min,i代入三次B样条方程中，得到矩阵T_min，如下式所示；

步骤4.23，将T_min代入B样条最优拟合逼近方程，如下式所示；

T_min ^TT_minα_min＝T_min ^TX_min

其中，X_min为(t_min,j,P_min,j′),j=1,2…,k+1组成的矩阵，α_min为需要求解的系数矩阵；

其中，

得到方程的解为α_min＝(T_min ^TT_min)^-1T_min ^TX_min；

步骤4.24，计算得到数据点集P_min'的B样条拟合曲线，其中α_i为系数矩阵α_min的第i个元素，t为数据点集P_min′所在的时间区间任意值。

本发明的有益效果是，使用B样条拟合遥测数据的自适应包络作为阈值，有助于提高信号的质量和准确性，减少噪声的干扰。使用自适应包络能够更好地适应信号的动态变化，可以提高故障检测的灵敏度和准确性，能够更好地捕捉到遥测数据的特征变化和趋势，提前发现潜在的故障迹象，并进行相应的维护和修复，对降低航天器故障损失和保障航天任务成功具有重大意义。

附图说明

图1是本发明基于B样条数据包络的遥测数据异常检测方法的流程图；

图2是卫星遥测数据曲线图；

图3是卫星遥测数据等间隔数据片段的最大最小值结果图；

图4是卫星遥测数据等间隔数据片段的最大最小值处理后的结果图；

图5是使用固定阈值进行异常检测的结果图；

图6是使用多项式拟合自适应包络的曲线图；

图7是使用B样条拟合自适应包络的曲线图；

图8是使用B样条拟合所得自适应包络作为阈值进行异常检测的曲线图。

具体实施方式

下面结合附图和具体实施方式对本发明进行详细说明。

实施例1

本发明基于B样条包络分析的遥测数据异常检测方法，如图1所示，具体按照以下步骤实施：

步骤1、设遥测数据序列为{x₁,x₂,…,x_M}，M为遥测数据长度，将遥测数据序列分为k个长度为d的数据片段，剩余长度不足d的视为第k+1个数据片段；这样，整个数据被分为k+1个片段，即{x₁₁,x₁₂,…,x_1d},{x₂₁,x₂₂,…,x_2d},…,{x_k1,x_k2,…,x_kd},{x_(k+1)1,x_(k+1)2，…，x_(k+1)m}；

对于遥测数据序列{x₁，x₂，…，x_M}，取长度为d的数据为一个片段，用d除以M的商为k，余数为m，将整个数据分为间隔为d的k个数据片段以及最后不足d的m个数据作为第k+1个片段，这样，整个数据被分为k+1个片段；

步骤2、在每个数据片段(即第i个数据片段)上求取最大值P_max,i与最小值P_min,i，将这些点分别构成数据点集P_max，P_min；

步骤1中将给定的遥测数据序列{x₁，x₂，…，x_M}划分为k+1个片段，找出每个数据片段的最大值P_max,i与最小值P_min，i(i＝1，2，…，k+1)，并将这些点作为数据点集P_max＝{P_max，1，P_max，2，…，P_max，k+1}，P_min＝{P_min，1，P_min，2，…，P_min，k+1}。

步骤3、求取每个数据片段上的标准差σ_i，其大小反映了数据的波动程度，将该片段对应最大值P_max,i与最小值P_min,i分别向上、向下移动，移动的距离标准差σ_i，得到新的数据点集P_max'，P_min'；

步骤1中给出的前k个长度为d的数据片段为{x₁₁，x₁₂，…，x_1d}，{x₂₁，x₂₂，|…，x_2d}，…，{x_k1，x_k2，…，x_kd}，第i＝(i＝1，2，…，k)个片段的均值u_i、标准差σ_i计算公式如下：

最后一个包含m个剩余数据的数据片段为{x_(k+1)1，x_(k+1)2，…，x_(k+)m}，其均值μ_k+1、标准差σ_k+1，计算如下：

将每个区间的最大值P_max,i，最小值P_min,i分别向上、向下移动标准差σ_i，得到P_max′＝{P_max，1′，P_max，2′，…，P_max，k+1′}，P_max′＝{P_max，1′，P_max，2′，…，P_max，k+1′}，作为新的数据点集。

步骤4、将数据点集P_max'，P_min'中的点进行B样条曲线拟合,对应于P_max'的曲线作为数据的上包络线，对应于P_min'的曲线作为数据的下包络线,将数据的上包络线和下包络线作为判断数据是否异常的“动态阈值”，从而完成对数据的异常检测；

步骤4.1，数据点集P_max'中的点可表示为：(t_max，j，P_max，j′)，j＝1，2，…，k+1，t_max,j为P_max,j′所对应的时间，采用均匀B样条曲线进行拟合；具体为：

步骤4.11，将T_max,0、作为数据点集P_max'的两个时间端点，T_max,0＝t_max,1，T_max,pmax＝t_max,k+1，并计算T_max,i，T_max，i＝T_max，0+ih_max，i＝0，1，…，p_max，其中p_max为拟合分段数(p_max<k)，h_max为每个分段的时间长度，

步骤4.12，将(t_max，j，P_max，j′)，j＝1，2，…，k+1和T_max,i代入三次B样条方程中，得到矩阵T_max，如下式所示；

三次B样条定义如下：

步骤4.13，将T_max代入B样条最优拟合逼近方程，如下式所示；

T_max ^TT_maxα_max＝T_max ^TX_max

其中，X_max为(t_max，j，P′_maxj)，j＝1，2，…，k+1组成的矩阵，α_max为需要求解的系数矩阵：

其中，

得到方程的解为α_max＝(T_max ^TT_max)^-1T_max ^TX_max；

步骤4.14，计算得到数据点集P_max'的B样条拟合曲线，其中α_i为系数矩阵α_max的第i个元素，t为数据点集P_max'所在的时间区间的任意值，即t_max,1≤t≤t_max,k+1；

步骤4.2，数据点集P_min'中的点可表示为：(t_min，j，P′_minj)，j＝1，2，…，k+1，t_min,j为P_min,j′所对应的时间，采用均匀B样条曲线进行拟合，具体如下：

步骤4.21，将T_min,0、T_min,pmin作为数据点集P_min'的两个时间端点，即T_min,0＝t_min,1，T_min,pmin＝t_min,k+1，并计算T_min,i，T_min，i＝T_min，0+ih_min，i＝0，1，…，p_min，其中p_min为拟合分段数(p_min<k)，h_min为每个分段的时间长度，

步骤4.22，将(t_min，j，P_min，j′)，j＝1，2，…，k+1和T_min,i代入三次B样条方程中，得到矩阵T_min，如下式所示；

步骤4.23，将T_min代入B样条最优拟合逼近方程，如下式所示；

T_min ^TT_minα_min＝T_min ^TX_min

其中，X_min为(t_min，j，P_min，j′)，j＝1，2，…，k+1组成的矩阵，α_min为需要求解的系数矩阵；

其中，

得到方程的解为α_min＝(T_min ^TT_min)^-1T_min ^TX_min；

步骤4.24，计算得到数据点集P_min'的B样条拟合曲线，其中α_i为系数矩阵α_min的第i个元素，t为数据点集P_min'所在的时间区间任意值，即t_min,1≤t≤t_min,k+1。

步骤4.3，由P_max'数据点集得到的B样条拟合曲线作为数据的上包络线，P_min'数据点集得到的B样条拟合曲线作为数据的下包络线，将数据的上包络线和下包络线作为判断数据是否异常的“动态阈值”，对于任意一组遥测数据，如果数据在阈值范围内，则判定为正常数据，若不在阈值范围内，则判定为异常数据；

传统的固定阈值异常检测，只适用于正常工况或者某一特殊工况，阈值设置过大，不能精确地检测异常值，阈值设置过小，则会造成误判，本发明提出的动态阈值能够提高对异常值的敏感程度和检测能力，并且具有一定的自适应性，对异常值检测的准确率更高，更接近于工程实际问题。

实施例2

本发明基于B样条包络分析的遥测数据异常检测方法，对遥测数据先进行数据片段划分，并求取每个数据片段的最大值和最小值作为上、下包络线的数据点集，将该数据点集分别向上、向下移动，移动的距离为对应数据片段的标准差，得到新的数据点集，使用B样条曲线对新的数据点集进行拟合得到上、下包络线，最终得到带有一定阈值的数据包络线，该包络线可以跟随遥测数据的趋势变化，从而可以准确地捕捉到遥测数据是否发生异常。

实施例3

本实施例采用某卫星的遥测数据进行分析实验，数据如图2所示；遥测数据(t_i，x_i)，i＝1，2，...，M在正常情况下的变化曲线具有一定的规律性，如果发生故障，则会明显不符合曲线的变化趋势，因此可以设定一定的阈值来对故障是否发生进行判断。如果设置一个固定阈值T来判断是否出现异常，阈值设置过大，不能精确地检测异常值，阈值设置过小，则会造成误判，因此固定阈值不能够准确地捕捉数据异常；

首先将该数据划分为k+1个数据片段，求取k+1个区间数据的最大值和最小值作为上、下包络线的数据点集，结果如图3所示，将该数据点集分别向上，向下移动，移动距离为对应数据片段数据的标准差，得到新的数据点集，结果如图4所示；

使用依据3σ原则取得的固定阈值进行异常检测出现漏检，结果如图5所示，使用多项式拟合对新的数据点集进行拟合得到的上、下包络图拟合曲线，结果如图6所示。使用B样条算法对新的数据点集进行拟合，得到了上、下包络线的拟合曲线，结果如图7所示。使用B样条拟合出遥测数据自适应包络作为阈值进行异常检测，结果如图8所示。

现在定义如下指标，对传统多项式拟合求取的数据包络线与本发明所提出方法进行对比，指标定义如下：

1.区间覆盖率

区间覆盖率ICP指真实值落入区间的概率，ICP的值越大，模型构造的区间越能最大限度地覆盖所有真实值，ICP可用于评价区间的准确度和可靠性。

其中，N代表样本数，CP_i是布尔值，用于表示第i个样本是否落入区间内。

其中，x_i是样本真实值，U_i和L_i是区间的上下界。

2.区间平均宽度百分比

宽度百分比IMWP是区间宽度和样本真实值的百分比的平均值；

其中，N表示样本数，x_i是样本真实值，U_i和L_i是区间的上下界。

该指标衡量了每一个样本值对应的区间宽度占实际观测值的百分比的平均数，IMWP的数值越小，模型越能满足实际需求，区间宽度一般与样本值呈相关性，因此，使用IMWP这个指标能对区间宽度给予更合理的评价。

3.平均带宽误差

平均带宽误差是用于衡量动态阈值预测准确度的重要指标，MWD的数值越小，模型的准确度越高。如果MWD的数值足够的小意味着样本真实值落入动态阈值区间的可能性较大，异常检测效果也相应的有所提高。

平均带宽误差用如下公式来衡量：

其中N表示样本数，x_i是样本真实值，U_i和L_i是区间的上下界，WD_i是样本的带宽误差。

表1为用多项式拟合方法和用本发明的方法在区间覆盖率，区间平均宽度百分比，平均带宽误差三个指标方面的结果对比：

表1两种方法的指标对比

	区间覆盖率	区间平均宽度百分比	平均带宽误差
				多项式拟合法	0.7559	0.1744	0.0482
B样条拟合法	0.9862	0.1648	0.0041

从表1可以看出本发明在三种指标方面均优于传统多项式拟合方法，说明本发明对于遥测数据的异常有更好的检测能力。通过以上分析，本发明所提出一种基于B样条包络分析据异常检测方法，可以较为精确实现遥测数据的异常检测。

Claims (7)

1.基于B样条包络分析的遥测数据异常检测方法，其特征在于，具体按照以下步骤实施：

步骤1、将遥测数据{x₁，x₂，…，x_M}分为k+1个片段；

步骤2、求取每个数据片段上最大值P_max，i与最小值P_min，i，构成数据点集P_max，P_min；

步骤3、将每个片段的最大值P_max，i与最小值P_min，i进行移动，得到新的数据点集P_max′，P_min′；

步骤4、将数据点集P_max′，P_min′中的点进行B样条曲线拟合，将P_max′的曲线作为上包络线，将P_min′的曲线作为下包络线，将上包络线和下包络线作为判断数据是否异常的动态阈值，完成数据的异常检测。

2.根据权利要求1所述的基于B样条包络分析的遥测数据异常检测方法，其特征在于，所述步骤1中，遥测数据为{x₁，x₂，…，x_M}，M为遥测数据长度，将遥测数据分为k个长度为d的数据片段，剩余长度不足d的视为第k+1个数据片段；这样，整个数据被分为k+1个片段，即{x₁₁，x₁₂，…，x_1d}，{x₂₁，x₂₂，…，x_2d}，…，{x_k1，x_k2，…，x_kd}，{x_(k+1)1，x_(k+1)2，…，x_(k+1)m}。

3.根据权利要求1所述的基于B样条包络分析的遥测数据异常检测方法，其特征在于，所述步骤2中，将遥测数据{x₁，x₂，…，x_M}划分为k+1个片段，找出每个数据片段的最大值P_max，i与最小值P_min，i(i＝1，2，…，k+1)，并将这些点作为数据点集P_max＝{P_max，1，P_max，2，…，P_max，k+1}，P_min＝{P_min，1，P_min，2，…，P_min，k+1}。

4.根据权利要求3所述的基于B样条包络分析的遥测数据异常检测方法，其特征在于，所述步骤3中，具体为：步骤1中的前k个长度为d的数据片段为{x₁₁，x₁₂，…，x_1d}，{x₂₁，x₂₂，…，x_2d}，…，{x_k1，x_k2，…，x_kd}，第i(i＝1，2，…，k)个片段的均值u_i、标准差σ_i计算公式如下：

最后一个包含m个剩余数据的数据片段为{x_(k+1)1，x_(k+1)2，…，x_(k+1)m}，其均值μ_k+1、标准差σ_k+1，计算如下：

将每个区间的最大值P_max，i，最小值P_min，i分别向上、向下移动标准差σ_i，得到P_max′＝{P_max，1′，P_max，2′，…，P_max，k+1′}，P_max′＝{P_max，1，P_max，2′，…，P_max，k+1′}，作为新的数据点集。

5.根据权利要求4所述的基于B样条包络分析的遥测数据异常检测方法，其特征在于，所述步骤4中，具体为：

步骤4.1，数据点集P_max′中的点可表示为：(t_max，j，P_max，j′)，j＝1，2，…，k+1，t_max，j为P_max，j′所对应的时间，采用均匀B样条曲线进行拟合；

步骤4.2，数据点集P_min′中的点可表示为：(t_min，j，P_min，j′)，j＝1，2，…，k+1，t_min，j′为P_min，j′所对应的时间，采用均匀B样条曲线进行拟合；

步骤4.3，由P_max′数据点集得到的B样条拟合曲线作为数据的上包络线，P_min′数据点集得到的B样条拟合曲线作为数据的下包络线，将数据的上包络线和下包络线作为判断数据是否异常的“动态阈值”，对于任意一组遥测数据，如果数据在阈值范围内，则判定为正常数据，若不在阈值范围内，则判定为异常数据。

6.根据权利要求5所述的基于B样条包络分析的遥测数据异常检测方法，其特征在于，所述步骤4.1中，拟合过程为：

步骤4.11，将T_max，0、作为数据点集P_max′的两个时间端点，T_max，0＝t_max，1，T并计算T_max，i，T_max，i＝T_max，0+ih_max，i＝0，1，…，p_max，其中p_max为拟合分段数，h_max为每个分段的时间长度，

步骤4.12，将(t_max，j，P_max，j′)，j＝1，2，…，k+1和T_max，i代入三次B样条方程中，得到矩阵T_max，如下式所示；

步骤4.13，将T_max代入B样条最优拟合逼近方程，如下式所示；

T_max ^TT_maxα_max＝T_max ^TX_max

其中，X_max为(t_max，j，P_max，j′)，j＝1，2，…，k+1组成的矩阵，α_max为需要求解的系数矩阵：

其中，

得到方程的解为α_max＝(T_max ^TT_max)^-1T_max ^TX_max；

步骤4.14，计算得到数据点集P_max′的B样条拟合曲线，其中α_i为系数矩阵α_max的第i个元素，t为数据点集P_max′所在的时间区间的任意值。

7.根据权利要求5所述的基于B样条包络分析的遥测数据异常检测方法，其特征在于，所述步骤4.2中，拟合过程为：

步骤4.21，将T_min，0、作为数据点集P_min′的两个时间端点，即T_min，0＝t_min，1，并计算T_min，i，T_min，i＝T_min，0+ih_min，i＝0，1，…，p_min，其中p_min为拟合分段数，h_min为每个分段的时间长度，

步骤4.22，将(t_min，j，P_min，j′)，j＝1，2，…，k+1和T_min，i代入三次B样条方程中，得到矩阵T_min，如下式所示；

步骤4.23，将T_min代入B样条最优拟合逼近方程，如下式所示；

T_min ^TT_minα_min＝T_min ^TX_min

其中，X_min为(t_min，j，P_min，j′)，j＝1，2，…，k+1组成的矩阵，α_min为需要求解的系数矩阵；

其中，

得到方程的解为α_min＝(T_min ^TT_min)^-1T_min ^TX_min；

步骤4.24，计算得到数据点集P_min′的B样条拟合曲线，其中α_i为系数矩阵α_min的第i个元素，t为数据点集P_min′所在的时间区间任意值。