卫星在轨故障传播和波及效应建模与预测方法及系统

摘要

本发明涉及卫星状态监测技术领域，公开一种卫星在轨故障传播和波及效应建模与预测方法及系统，以分析出卫星的故障传播路径并确保预测的精度。方法包括：收集卫星遥测参数历史数据；使用格兰杰因果关系模型判断卫星遥测参数之间是否具有格兰杰因果关系；基于互相关函数判断参数之间的格兰杰因果关系的相关强度、相关方向、滞后时间，剔除相关强度低于阈值的相关关系；根据剔除之后剩余的遥测参数之间的关系构造邻接矩阵，画出参数之间的因果关系图；根据遥测参数之间的邻接矩阵，应用解释结构模型，建立卫星遥测参数故障传播层级图；根据格兰杰因果关系模型，定量构造故障传播模型，根据故障传播模型来判断卫星故障发生的波及效应。

说明书

技术领域

本发明涉及卫星状态监测技术领域，尤其涉及一种卫星在轨故障传播和波及效应建模与预测方法及系统。

背景技术

目前关于故障传播和波及效应的建模与预测方法，主要有两种。一种是直接从系统架构出发，直接构建系统部件传播关系，系统分析方法用树和图来表示系统的结构，需要熟悉系统原理，较为复杂，需要很多的先验知识，难以建立。另一种是从数据的角度出发，从数据来分析出系统的故障传播关系，比如说贝叶斯网络，这是目前用的比较多的方法。但是这些方法在故障传播路径分析中并不总是成功的，而且难以建立起正确的故障传播路径。

由于卫星在轨普遍存在传播特性和波及效应。无论是部件还是分系统故障，若不加控制任其发展，都可能会降低卫星性能、功能或寿命。比如，某卫星发生陀螺的缓变故障后，导致卫星姿态逐渐偏离，控制系统失去自主故障判断机会；随着陀螺持续饱和，陀螺测量输出持续超差，卫星进入全姿态捕获模型，导致其无法执行部分功能。然而，关于故障传播和波及效应的先验数据是缺乏的，需要充分利用遥测参数之间的关系，分析在轨故障的传播特性和波及效应，预测故障及其处置措施对系统或整星性能或功能的影响。实践证明，卫星在轨故障是可以预测的，特别是一些缓变故障，可以通过遥测数据进行早期检测和预防。对于突发故障，虽然难以预测其发生的时机，但是可以通过遥测数据的分析预测其传播过程和波及效应，以便采取有效的故障初值措施降低故障影响。

发明内容

本发明目的在于公开一种卫星在轨故障传播和波及效应建模与预测方法及系统，以分析出卫星的故障传播路径，在卫星发生故障时采取措施降低故障影响，并确保预测的精度。

为达上述目的，本发明公开一种卫星在轨故障传播和波及效应建模与预测方法，包括以下步骤：

步骤S1、收集卫星遥测参数历史数据，使用前值填充的方法对缺失值进行填充；

步骤S2、使用格兰杰因果关系模型判断卫星遥测参数之间是否具有格兰杰因果关系；

步骤S3、基于互相关函数判断参数之间的格兰杰因果关系的相关强度、相关方向、滞后时间，剔除相关强度低于阈值的相关关系；

步骤S4、根据剔除之后剩余的遥测参数之间的关系构造邻接矩阵，画出参数之间的因果关系图；

步骤S5、根据遥测参数之间的邻接矩阵，应用解释结构模型，建立卫星遥测参数故障传播层级图；

步骤S6、根据格兰杰因果关系模型，定量的构造步骤S3剔除之后剩余的卫星遥测参数之间的故障传播模型，根据故障传播模型来判断卫星故障发生的波及效应。

优选地，所述步骤S2中，遥测参数之间的格兰杰因果关系的确定包括以下步骤：

步骤S21、对遥测数据进行去均值处理；

步骤S22、选用ADF对遥测数据进行单位根检验，进行ADF检验的步骤如下：

假设时间序列{y

y

式中ε

Ψ(B)y

令：

ρ＝φ

ζ

将滞后多项式Ψ(B)分解为：

Ψ(B)＝(1-φ

转化为：

ψ(B)y

整理可得：

y

若服从AR(P)过程的时间序列有且只有一个单位根，则其特征方程为：

1-φ

其中，z代表特征根，有且只有一个根则：

Ψ(1)＝1-φ

等价于：ρ＝1；

故对于服从AR(P)过程的时间序列单位根检验等价于检验公式(7)中是否有ρ＝1；使得ADF检验根据以下四种情况进行单位根检验：

情况一：随机游走过程判定，时间序列不含常数项不含趋势项，由公式(7)检验存在单位根ρ＝1；

情况二：单位根过程判定，时间序列含常数项不含趋势项，在公式(10)估计模型中检验ρ＝1，其中，常数项为α；

y

情况三：时间序列含趋势项不含常数项，在公式(10)估计模型中检验ρ＝1；

情况四：趋势非平稳过程判定，时间序列含常数项含趋势项，在公式(11)估计模型中检验ρ＝1，Lt表示趋势项；

y

对于情况二，在ρ＝1成立时，对公式(10)进行最小二乘估计，可得ρ的超极限估计

其中，N()代表正态分布函数；W()代表标准维纳过程；r代表0到1之间的参数；以最小二乘估计

对于检验ρ＝1的T统计量，有如下极限分布：

其中，

对于情况一和情况四，Dickey和Fuller同样证明检验ρ＝1的统计量Z

对于情况三，T统计量的极限分布为常规的T分布，因此可用常规T检验，查T分布表；在实际应用中，用如下三种回归模型进行ADF检验：

故检验单位根的假设ρ＝1在上述三种回归模型中变为回归系数

如果数据没有通过单位根检验，则需要进行差分处理，直到数据平稳为止，所述差分处理为取时间序列内连续相邻两项之差；

步骤S23、根据BIC准则，确定格兰杰因果关系检验中的最大延迟，BIC可由下式求得：

BIC＝kln(n)-2ln(L) 公式(18)

式中，k为模型参数个数，n为样本数量，L为似然函数，取BIC最大时的个数作为格兰杰因果关系检验的最大延迟；

确定最大迟延后，采用多项式拟合的方法构建回归方程，多项式拟合的步骤如下；

多项式定义为：

其中，M是多项式的阶数，w0,…,w

其中，E()代表误差，x

步骤S24、进行杜滨—沃特森检验，检测回归后的残差是否服从正态分布，如果没有通过检验，则参数之间不存在格兰杰因果关系；

所述杜滨—沃特森检验的步骤如下：

通过构建统计量d，依据d值所落的区间范围判定残差间是否存在自相关；d统计量可由下式求得：

式中e

为了改进无法判定的情况，提出新统计量d

d

式中a、b满足：

式中，E(d

P＝2(n-1)-tr[(X′AX)(X′X)

Q＝2(3n-4)-2tr[X′A

其中，n为样本容量，k为参数个数，tr为求矩阵的迹，X为样本的值，X′为X的转置；判定：

若：d

若：4-d

步骤S25、在进行格兰杰因果关系检验时，给定两个时间序列{X

首先验证X是否是Y变化的格兰杰原因，该检验要求估计以下的无约束的回归模型：

以及有约束的回归模型：

其中，α

分别计算所述无约束的回归模型和所述有约束的回归模型的残差平方和为：RSS

构造F统计量，n为样本容量：

检验零假设，如果F≥F

交换X与Y的位置，以同样的方式检验Y是否是引起X变化的格兰杰原因；

对一对遥测参数进行格兰杰因果关系检验后，得到它们的相关关系，从而得到节点之间边的关系。

优选地，所述遥测参数之间格兰杰因果关系强弱的确定包括以下步骤：

对遥测时间序列数据{x

其中，μ为均值，δ为标准差；对于两条归一化的曲线G＝{g

上式中，0的个数是|i|个，其中，-s＜i＜s；定义G

于是：

则最小值或者最大值的指标分别是：

sc

sc

令：

FCC(G,H)为最终得到的衡量遥测参数之间相关性的元组，FCC(G,H)[0]∈[-1,1],FCC(G,H)[0]的绝对值越接近于1表示参数之间的相关强度越强，正值的FCC(G,H)[0]表示参数之间是正相关，负值的FCC(G,H)[0]表示参数之间是负相关，,FCC(G,H)[1]代表参数之间影响的滞后时间，当|FCC(G,H)[0]|＜ε时，ε为阈值，代表参数之间的相关关系不成立，因此它们之间关系被剔除，这样就得到筛选之后的参数之间的相关关系；其中，FCC(G,H)[0]表示FCC(G,H)数组的第一个元素，FCC(G,H)[1]表示FCC(G,H)数组的第二个元素。

优选地，所述卫星遥测参数故障传播层级图的确定如下：

根根筛选之后的卫星参数之间的格兰杰因果关系，然后将其转化成矩阵的形式，构建要素之间的直接关系矩阵，最后利用Warshall算法计算可达矩阵，并进一步计算求解，得到系统的多级递阶层次结构模型。

优选地，所述故障传播模型的构造如下：

根据FCC(G,H)[1]得到的参数间的滞后时间，利用格兰杰因果关系模型，估计如下的回归：

公式中的e和f为FCC(G,H)[1]得到的卫星遥测参数之间的滞后时间，α

为达上述目的，本发明还公开一种基于改进的格兰杰因果关系模型的卫星在轨故障传播和波及效应建模与预测系统，包括存储器、处理器以及存储在存储器上并可在处理器上运行的计算机程序，所述处理器执行所述计算机程序时实现上述方法的步骤。

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

本发明基于互相关函数判断参数之间的格兰杰因果关系的相关强度、相关方向、滞后时间，剔除相关强度低于阈值的相关关系，解决了传统的因果关系建模无法判断参数之间相关关系的强弱问题，并基于剔除相关强度低于阈值的相关关系之后剩余的遥测参数得到故障传播层级图，分析遥测参数之间的层级关系，解决了有关故障传播图无法得出参数的层次结构的问题，并确保了根据故障传播模型来判断卫星故障发生的波及效应的精度。

下面将参照附图，对本发明作进一步详细的说明。

附图说明

构成本申请的一部分的附图用来提供对本发明的进一步理解，本发明的示意性实施例及其说明用于解释本发明，并不构成对本发明的不当限定。在附图中：

图1是本发明实施例的卫星在轨故障传播和波及效应建模与预测方法流程示意图。

图2是本发明实施例DW检验统计量d取值区间检验结果示意图。

图3是本发明实施例改进DW检验统计量d取值区间检验结果示意图。

图4是本发明实施例初始卫星参数之间的格兰杰因果关系图。

图5是本发明实施例筛选之后的卫星参数之间的故障传播图。

具体实施方式

以下结合附图对本发明的实施例进行详细说明，但是本发明可以由权利要求限定和覆盖的多种不同方式实施。

实施例1

本实施例公开一种卫星在轨故障传播和波及效应建模与预测方法，如图1所示，包括以下步骤：

步骤S1、收集卫星遥测参数历史数据，使用前值填充的方法对缺失值进行填充。

步骤S2、使用格兰杰因果关系模型判断卫星遥测参数之间是否具有格兰杰因果关系。

该步骤中，遥测参数之间的格兰杰因果关系的确定包括以下步骤：

步骤S21、对遥测数据进行去均值处理。以此保证卫星遥测数据满足格兰杰因果关系检验中对数据围绕y轴波动的要求。

步骤S22、选用ADF对遥测数据进行单位根检验。为确定卫星遥测数据是否满足格兰杰因果关系检验中待检测数据为平稳序列的要求，需对遥测数据进行单位根检验。进行ADF检验的步骤如下：

假设时间序列{y

y

式中ε

Ψ(B)y

令：

ρ＝φ

ζ

将滞后多项式Ψ(B)分解为：

Ψ(B)＝(1-φ

转化为：

ψ(B)y

整理可得：

y

若服从AR(P)过程的时间序列有且只有一个单位根，则其特征方程为：

1-φ

其中，z代表特征根，有且只有一个根则：

Ψ(1)＝1-φ

等价于：ρ＝1；

故对于服从AR(P)过程的时间序列单位根检验等价于检验公式(7)中是否有ρ＝1；使得ADF检验根据以下四种情况进行单位根检验：

情况一：随机游走过程判定，时间序列不含常数项不含趋势项，由公式(7)检验存在单位根ρ＝1；

情况二：单位根过程判定，时间序列含常数项不含趋势项，在公式(10)估计模型中检验ρ＝1，其中，常数项为α；

y

情况三：时间序列含趋势项不含常数项，在公式(10)估计模型中检验ρ＝1；

情况四：趋势非平稳过程判定，时间序列含常数项含趋势项，在公式(11)估计模型中检验ρ＝1，Lt表示趋势项；

y

对于情况二，在ρ＝1成立时，对公式(10)进行最小二乘估计，可得ρ的超极限估计

其中，N()代表正态分布函数；W()代表标准维纳过程；r代表0到1之间的参数；以最小二乘估计

对于检验ρ＝1的T统计量，有如下极限分布：

其中，

对于情况一和情况四，Dickey和Fuller同样证明检验ρ＝1的统计量Z

对于情况三，T统计量的极限分布为常规的T分布，因此可用常规T检验，查T分布表；在实际应用中，用如下三种回归模型进行ADF检验：

故检验单位根的假设ρ＝1在上述三种回归模型中变为回归系数

如果数据没有通过单位根检验，则需要进行差分处理，直到数据平稳为止，所述差分处理为取时间序列内连续相邻两项之差；

步骤S23、根据BIC准则，确定格兰杰因果关系检验中的最大延迟。BIC可由下式求得：

BIC＝kln(n)-2ln(L) 公式(18)

式中，k为模型参数个数，n为样本数量，L为似然函数，取BIC最大时的个数作为格兰杰因果关系检验的最大延迟；

确定最大迟延后，采用多项式拟合的方法构建回归方程，多项式拟合的步骤如下；

多项式定义为：

其中，M是多项式的阶数，w0,…,wM是多项式的系数，记做W，使用下式来进行误差评估：

其中，E()代表误差，x

步骤S24、进行杜滨—沃特森检验，检测回归后的残差是否服从正态分布，如果没有通过检验，则参数之间不存在格兰杰因果关系。

所述杜滨—沃特森检验的步骤如下：

通过构建统计量d，依据d值所落的区间范围判定残差间是否存在自相关；d统计量可由下式求得：

式中e

为了改进无法判定的情况，提出新统计量d

d

式中a、b满足：

式中，E(d

P＝2(n-1)-tr[(X′AX)(X′X)

Q＝2(3n-4)-2tr[X′A

其中，n为样本容量，k为参数个数，tr为求矩阵的迹，X为样本的值，X′为X的转置。可判定：

若：d

若：4-d

具体判断逻辑如图3所示，其中，

步骤S25、在进行格兰杰因果关系检验时，给定两个时间序列{X

首先验证X是否是Y变化的格兰杰原因，该检验要求估计以下的无约束的回归模型：

以及有约束的回归模型：

其中，α

分别计算所述无约束的回归模型和所述有约束的回归模型的残差平方和为：RSS

构造F统计量，n为样本容量：

检验零假设，如果F≥F

交换X与Y的位置，以同样的方式检验Y是否是引起X变化的格兰杰原因；

对一对遥测参数进行格兰杰因果关系检验后，得到它们的相关关系，从而得到节点之间边的关系。

步骤S3、基于互相关函数判断参数之间的格兰杰因果关系的相关强度、相关方向、滞后时间，剔除相关强度低于阈值的相关关系。

该步骤中，所述遥测参数之间格兰杰因果关系强弱的确定包括以下步骤：

对遥测时间序列数据{x

其中，μ为均值，δ为标准差；对于两条归一化的曲线G＝{g

上式中，0的个数是|i|个，其中，-s＜i＜s；定义G

于是：

则最小值或者最大值的指标分别是：

sc

sc

令：

FCC(G,H)为最终得到的衡量遥测参数之间相关性的元组，FCC(G,H)[0]∈[-1,1],FCC(G,H)[0]的绝对值越接近于1表示参数之间的相关强度越强，正值的FCC(G,H)[0]表示参数之间是正相关，负值的FCC(G,H)[0]表示参数之间是负相关，,FCC(G,H)[1]代表参数之间影响的滞后时间，当|FCC(G,H)[0]|＜ε时，ε为阈值，代表参数之间的相关关系不成立，因此它们之间关系被剔除，这样就得到筛选之后的参数之间的相关关系；其中，FCC(G,H)[0]表示FCC(G,H)数组的第一个元素，FCC(G,H)[1]表示FCC(G,H)数组的第二个元素。

步骤S4、根据剔除之后剩余的遥测参数之间的关系构造邻接矩阵，画出参数之间的因果关系图。

步骤S5、根据遥测参数之间的邻接矩阵，应用解释结构模型，建立卫星遥测参数故障传播层级图。

该步骤中，所述卫星遥测参数故障传播层级图的确定如下：

根根筛选之后的卫星参数之间的格兰杰因果关系，然后将其转化成矩阵的形式，构建要素之间的直接关系矩阵，最后利用Warshall算法计算可达矩阵，并进一步计算求解，得到系统的多级递阶层次结构模型。

步骤S6、根据格兰杰因果关系模型，定量的构造步骤S3剔除之后剩余的卫星遥测参数之间的故障传播模型，根据故障传播模型来判断卫星故障发生的波及效应。

该步骤中，所述故障传播模型的构造如下：

根据FCC(G,H)[1]得到的参数间的滞后时间，利用格兰杰因果关系模型，估计如下的回归：

公式中的e和f为FCC(G,H)[1]得到的卫星遥测参数之间的滞后时间，α

实施例2

依托上述实施例中方法，本实施例公开一种卫星在轨故障传播和波及效应建模与预测方法的具体算例，包括以下步骤：

S11：收集卫星遥测参数历史数据；使用前值填充的方法对缺失值进行填充。

参数数据包括温度、电流、电压等遥测参数，因为卫星在轨运行过程中，需要对其部件的多种物理量进行测量，这些物理量就是遥测参数，所得到的数据就是遥测数据，遥测参数以表示形式为时间序列，具体表现为

S12：使用格兰杰因果关系模型判断卫星遥测参数之间是否具有格兰杰因果关系。遥测参数间的格兰杰因果关系进行两两判断，得出关于所有参数相关关系的邻接矩阵，画出其初步的格兰杰因果关系图，参见图4。

S13：基于互相关函数判断参数之间的格兰杰因果关系的相关强度、相关方向、滞后时间，剔除相关强度较弱的相关关系。

对已经判断出来具有相关关系的参数，使用FCC(G,H)对其进行相关强度、相关方向以及滞后时间的分析，参见表1，其中v,e代表遥测参数的名称。

表1：

S4：根据剔除之后的遥测参数之间的关系构造邻接矩阵，画出参数之间的因果关系图。

把相关强度过小的参数之间的相关关系剔除，得到新的邻接矩阵，画出筛选之后的格兰杰因果关系图如图5所示。

S5：根据遥测参数之间的邻接矩阵，应用解释结构模型，建立卫星遥测参数故障传播层级图。

对筛选之后的邻接矩阵，构建参数之间的直接关系矩阵，最后利用Warshall算法计算可达矩阵，并进一步计算求解，得到系统的多级递阶层次结构模型，得到的层次结构模型如表2所示。

表2：

S6：根据格兰杰因果关系模型，定量的构造剔除之后剩余的卫星遥测参数之间的故障传播模型，根据故障传播模型来判断卫星故障发生的波及效应，如表3所示。

表3：

实施例3

本实施例公开一种基于改进的格兰杰因果关系模型的卫星在轨故障传播和波及效应建模与预测系统，包括存储器、处理器以及存储在存储器上并可在处理器上运行的计算机程序，所述处理器执行所述计算机程序时实现上述两实施例中相对应方法的步骤。

综上，本发明上述各实施例所分别公开的卫星在轨故障传播和波及效应建模与预测方法及系统，至少具有以下有益效果：

本发明基于互相关函数判断参数之间的格兰杰因果关系的相关强度、相关方向、滞后时间，剔除相关强度低于阈值的相关关系，解决了传统的因果关系建模无法判断参数之间相关关系的强弱问题，并基于剔除相关强度低于阈值的相关关系之后剩余的遥测参数得到故障传播层级图，分析遥测参数之间的层级关系，解决了有关故障传播图无法得出参数的层次结构的问题，并确保了根据故障传播模型来判断卫星故障发生的波及效应的精度。

以上所述仅为本发明的优选实施例而已，并不用于限制本发明，对于本领域的技术人员来说，本发明可以有各种更改和变化。凡在本发明的精神和原则之内，所作的任何修改、等同替换、改进等，均应包含在本发明的保护范围之内。