一种基于KPCA的工业过程故障诊断方法

摘要

本发明涉及一种基于KPCA的工业过程故障诊断方法，对工业生产过程的历史正常数据利用KPCA方法提取出主元子空间的负载方向和残差子空间的负载方向，对已知故障的历史故障数据利用KPCA方法提取出主元子空间的负载方向和残差子空间的负载方向，对已知故障类型的历史故障数据进行基于T2统计量重构和基于SPE统计量重构，提取出历史故障数据针对T2统计量重构的故障特征方向和针对SPE统计量重构的故障特征方向，得到重构故障特征方向集合，实时采集工业生产过程的新数据，利用KPCA方法计算新数据的T2统计量和SPE统计量，判断实时采集的工业生产过程是否发生故障，利用重构故障特征方向集合对新数据进行故障方向重构，判断当前的工业生产过程的故障类型。

说明书

技术领域

本发明属于过程控制领域，具体涉及一种基于KPCA的工业过程故障诊断方法。

背景技术

运用工业生产过程中得到的数据进行建模并对生产过程中出现的故障进行检测与诊断是一个很有挑战性的问题，近些年来受到广泛关注。许多学者研究了利用PCA、PLS等多变量统计方法进行生产过程中故障的检测与诊断。这些方法能够提取测量数据的潜在特征，并且根据这些特征利用统计学原理定义检测统计量以及其在正常生产情况下的控制限。进行在线监测时，通过新的采样数据计算相应的统计量，若结果超限报警，则认为有故障发生。故障被检测出来之后，还需快速诊断出故障发生的原因，以便能够在较短时间内将生产过程恢复到正常情况。

故障重构是故障诊断的一种重要手段，该方法在生产过程发生故障时，结合已知的故障特征方向从故障数据中恢复其正常值。其目的在于，通过利用不同类型的故障特征方向对故障数据进行恢复，找出能够使故障数据成功恢复至正常值的故障特征方向，最终经过故障特征方向与实际故障的匹配找出故障发生的原因。

传统的基于PCA方法的故障重构技术将故障数据空间分解为两个相互垂直的子空间，即主元子空间和残差子空间，并将这两个子空间的负载方向分别作为重构统计量T²和SPE的故障特征方向。通过PCA方法得到的负载方向能够反映故障数据的分布方向，将其作为故障特征方向有一定的合理性。但是，传统的故障重构方法属于线性建模方法，不能反映出数据的非线性特征；并且该方法只关注故障数据的内部关系，不能有效地区别数据中的故障信息与正常信息，直接使用该方法提取的故障特征方向进行故障重构会导致重构过度的情况。

发明内容

针对现有技术的不足，本发明提出一种基于KPCA的工业过程故障诊断方法。

本发明技术方案如下：

一种基于KPCA的工业过程故障诊断方法，包括以下步骤：

步骤1：对工业生产过程的历史正常数据利用KPCA方法进行子空间划分，得到历史正常数据的主元子空间的负载方向P和历史正常数据的残差子空间的负载方向P^*；

步骤1.1：对工业生产过程的n个历史正常数据进行m次采样，获取工业生产过程的历史正常数据集X；

步骤1.2：对历史正常数据集X进行中心化处理后将其映射到高维特征空间，得到高维特征空间的历史正常数据Φ(X)；

步骤1.3：对高维特征空间的历史正常数据Φ(X)进行PCA分解，分解为历史正常数据的主元子空间和历史正常数据的残差子空间，得到历史正常数据的主元部分和历史正常数据的残差部分

步骤1.4：计算历史正常数据的主元子空间的负载方向P和历史正常数据的残差子空间的负载方向P^*。

步骤2：对工业生产过程已知故障的历史故障数据利用KPCA方法进行子空间划分，得到历史故障数据的主元子空间的负载方向P_f和历史故障数据的残差子空间的负载方向

步骤2.1：对工业生产过程的n个历史故障数据进行m次采样，获取工业生产过程的已知故障的历史故障数据集X_f；

步骤2.2：对历史故障数据集X_f进行中心化处理后将其映射到高维特征空间，得到高维特征空间的历史故障数据Φ(X_f)；

步骤2.3：对高维特征空间的历史故障数据Φ(X_f)进行PCA分解，分解为历史故障数据的主元子空间和历史故障数据的残差子空间，得到历史故障数据的主元部分和历史正常数据的残差部分

步骤2.4：计算历史故障数据的主元子空间的负载方向P_f和历史故障数据的残差子空间的负载方向

步骤3：对历史故障数据进行基于T²统计量重构，提取出历史故障数据针对T²统计量重构的故障特征方向

步骤3.1：将历史正常数据的主元部分映射到历史故障数据的主元子空间，获得针对T²统计量重构的历史故障数据的主元子空间中与故障相关的特征方向P_fr；

步骤3.1.1：将历史正常数据的主元部分映射到历史故障数据的主元子空间，得到历史正常数据的主元部分映射到历史故障数据的主元子空间的映射矩阵

步骤3.1.2：对历史正常数据的主元部分映射到历史故障数据的主元子空间的映射矩阵进行PCA分解，得到映射矩阵在历史故障数据的主元子空间的负载方向P_r；

步骤3.1.3：计算映射矩阵在其历史故障数据的主元子空间的负载方向P_r上的得分T_r，计算历史故障数据的主元部分在映射矩阵的历史故障数据的主元子空间的负载方向P_r上的得分T_fr；

步骤3.1.4：定义在映射矩阵的历史故障数据的主元子空间的负载方向P_r上故障数据与正常数据的得分比矩阵RT，其中，μ＝1，2，…，d，d∈1…n，RT_μ为得分比矩阵RT的第μ个得分比，var(·)表示方差运算，T_fr(：，μ)为得分T_fr的第μ列，T_r(：，μ)为得分T_r的第μ列；

步骤3.1.5：设定得分比率阂值η，将得分比矩阵RT中的得分比与得分比率阂值η对比，在映射矩阵在历史故障数据的主元子空间的负载方向P_r中提取出得分比大于得分比率阂值η的负载方向，得到针对T²统计量重构的历史故障数据的主元子空间中与故障相关的特征方向P_fr。

步骤3.2：将历史正常数据的主元部分映射到历史故障数据的残差子空间，获得针对T²统计量重构的历史故障数据的残差子空间中与故障相关的特征方向

步骤3.2.1：将历史正常数据的主元部分映射到历史故障数据的残差子空间，得到历史正常数据的主元部分映射到历史故障数据的残差子空间的映射矩阵

步骤3.2.2：对历史正常数据的主元部分映射到历史故障数据的残差子空间的映射矩阵进行PCA分解，得到映射矩阵在历史故障数据的残差子空间的负载方向

步骤3.2.3：计算映射矩阵在其历史故障数据的残差子空间的负载方向上的得分计算历史故障数据的残差部分在映射矩阵的历史故障数据的残差子空间的负载方向上的得分

步骤3.2.4：定义在映射矩阵的历史故障数据的残差子空间的负载方向上故障数据与正常数据的得分比矩阵RT^*，其中，为得分比矩阵RT^*的第μ个得分比，为得分的第μ列，为得分的第μ列；

步骤3.2.5：将得分比矩阵RT^*中的得分比与得分比率阈值η对比，在映射矩阵在历史故障数据的残差子空间的负载方向中提取出得分比大于得分比率阈值η的负载方向，得到针对T²统计量重构的历史故障数据的残差子空间中与故障相关的特征方向

步骤3.3：获得历史故障数据针对T²统计量重构的故障特征方向

步骤4：对历史故障数据进行基于SPE统计量重构，提取出历史故障数据针对SPE统计量重构的故障特征方向

步骤4.1：将历史正常数据的残差部分映射到历史故障数据的主元子空间，获得针对SPE统计量重构的历史故障数据的主元子空间中与故障相关的特征方向V_fr；

步骤4.1.1：将历史正常数据的残差部分映射到历史故障数据的主元子空间，得到历史正常数据的残差部分映射到历史故障数据的主元子空间的映射矩阵

步骤4.1.2：对历史正常数据的残差部分映射到历史故障数据的主元子空间的映射矩阵进行PCA分解，得到映射矩阵在历史故障数据的主元子空间的负载方向V_r；

步骤4.1.3：定义在映射矩阵的历史故障数据的主元子空间的负载方向V_r上故障数据与正常数据的平方误差比矩阵RT_e，其中，k＝1，2，…，d，RT_e，k为平方误差比矩阵RT_e的第k个平方误差比，||·||表示欧氏距离，V_r，k为负载方向V_r的第k列；

步骤4.1.4：设定平方误差比率阈值θ，将平方误差比矩阵RT_e中的平方误差比与平方误差比率阈值θ对比，在映射矩阵在历史故障数据的主元子空间的负载方向V_r中提取出平方误差比大于平方误差比率阈值θ的负载方向，得到针对SPE统计量重构的历史故障数据的主元子空间中与故障相关的特征方向V_fr。

步骤4.2：将历史正常数据的残差部分映射到历史故障数据的残差子空间，获得针对SPE统计量重构的历史故障数据的残差子空间中与故障相关的特征方向

步骤4.2.1：将历史正常数据的残差部分映射到历史故障数据的残差子空间，得到历史正常数据的残差部分映射到历史故障数据的残差子空间的映射矩阵

步骤4.2.2：对历史正常数据的残差部分映射到历史故障数据的残差子空间的映射矩阵进行PCA分解，得到映射矩阵在历史故障数据的残差子空间的负载方向

步骤4.2.3：定义在映射矩阵的历史故障数据的残差子空间的负载方向上故障数据与正常数据的平方误差比矩阵其中，为平方误差比矩阵的第k个平方误差比，为负载方向的第k列；

步骤4.2.4：将平方误差比矩阵中的平方误差比与平方误差比率阈值θ对比，在映射矩阵在历史故障数据的残差子空间的负载方向中提取出平方误差比大于平方误差比率阈值θ的负载方向，得到针对SPE统计量重构的历史故障数据的残差子空间中与故障相关的特征方向

步骤4.3：获得历史故障数据针对SPE统计量重构的故障特征方向

步骤5：重复步骤1-步骤4，对c类已知故障类型的历史故障数据进行基于T²统计量重构和基于SPE统计量重构，得到重构故障特征方向集合；

步骤6：实时采集工业生产过程的新数据，利用KPCA方法进行子空间划分，计算新数据的T²统计量和SPE统计量，并根据新数据的T²统计量和SPE统计量判断实时采集的工业生产过程是否发生故障；

步骤6.1：实时采集工业生产过程的新数据集X_new，对其进行中心化处理后将其映射到高维特征空间，得到高维特征空间的新数据Φ(X_new)；

步骤6.2：对高维特征空间的新数据Φ(X_new)进行PCA分解，得到新数据的得分t_new＝<P，Φ(X_new)>和新数据的残差e_new＝(I-PP^T)Φ(X_new)，其中，I为单位矩阵；

步骤6.3：计算新数据的T²统计量和新数据的SPE统计量SPE_new；

步骤6.4：判断当前的工业生产过程是否发生故障：若新数据的T²统计量或SPE统计量偏离其对应的正常值域，执行步骤7，否则返回步骤6.1；

步骤7：利用重构故障特征方向集合对实时采集的工业生产过程的新数据进行故障方向重构，判断当前的工业生产过程的故障类型；

步骤7.1：选取故障类型l，l＝1，2，…，c，利用第l类历史故障数据针对T²统计量重构的故障特征方向对新数据进行故障方向重构：其中，为针对第l类故障T²统计量重构后的新数据，P_rec，l为第l类历史故障数据针对T²统计量重构的故障特征方向；

步骤7.2：利用第l类历史故障数据针对SPE统计量重构的故障特征方向对新数据进行故障方向重构：其中，Φ_SPE，f，l(X_new)为针对第l类故障SPE统计量重构后的新数据，V_rec，l为第l类历史故障数据针对SPE统计量重构的故障特征方向；

步骤7.3：计算针对第l类故障T²统计量重构后的新数据的得分：计算针对第l类故障SPE统计量重构后的新数据Φ_SPE，f，l(X_new)的残差e_f，l＝(I-PP^T)Φ_SPE，f，l(X_new)；

步骤7.4：计算针对第l类故障T²统计量重构后的新数据的T²统计量和针对第l类故障SPE统计量重构后的新数据Φ_SPE，f，l(X_new)的SPE统计量$>{\mathrm{SPE}}_{f,l}=e_{f,l}^T{e}_{f,l};$>

步骤7.5：判断当前的工业生产过程的故障类型：若针对第l类故障T²统计量重构后的新数据的T²统计量和针对第l类故障SPE统计量重构后的新数据Φ_SPE，f，l(X_new)的SPE统计量SPE_f，l均未偏离其对应的正常值域，则当前的工业生产过程的故障类型为第l类故障，否则，重新选择故障类型l，返回步骤7.1；

本发明的有益效果是：

本发明提出的一种基于KPCA的工业过程故障诊断方法，解决了非线性过程的故障分离问题，采用KPCA方法将故障数据空间分解为主元子空间和残差子空间，利用所得的负载方向对正常数据进行投影，利用PCA方法对投影的数据进行分析，通过比较各个方向上故障数据与正常数据的得分提取出引起统计量超限的故障方向，利用重构故障特征方向集合依次对检测到的故障数据进行重构，只有当前故障所对应的重构故障特征方向能够正确去除数据中的故障信息，消除检测统计量超限报警现象，据此可确定故障类别，达到故障分离的目的。

附图说明

图1为本发明具体实施方式中的电熔镁炉的基本工作原理示意图；

其中，1-变压器，2-电路短网，3-电机升降装置，4-电级，5-电熔镁炉炉壳，6-车体，7-电弧，8-炉料；

图2为本发明具体实施方式中的基于KPCA的工业过程故障诊断方法的流程图；

图3为本发明具体实施方式中的对历史故障数据进行基于T²统计量重构的流程图；

图4为本发明具体实施方式中的对历史故障数据进行基于SPE统计量重构的流程图；

图5为本发明具体实施方式中的对电熔镁炉冶炼过程的新数据进行故障诊断的流程图；

图6为本发明具体实施方式中的新数据1的T²统计量检测图和SPE统计量检测图；

其中，(a)为新数据1的T²统计量检测图，(b)为新数据1的SPE统计量检测图；

图7为本发明具体实施方式中的针对第l类故障T²统计量重构后的新数据1的T²统计量检测图和针对第l类故障SPE统计量重构后的新数据1的SPE统计量检测图；

其中，(a)为第l＝1类故障T²统计量重构后的新数据1的T²统计量检测图；

(b)为针对第l＝1类故障SPE统计量重构后的新数据1的SPE统计量检测图；

(c)为针对第l＝2类故障T²统计量重构后的新数据1的T²统计量检测图；

(d)为针对第l＝2类故障SPE统计量重构后的新数据1的SPE统计量检测图；

图8为本发明具体实施方式中的新数据1去除与故障无关的第l＝2类故障信息后的T²统计量检测图和SPE统计量检测图；

其中，(a)为新数据1去除与故障无关的第l＝2类故障信息后的T²统计量检测图；

(b)为新数据1去除与故障无关的第l＝2类故障信息后的SPE统计量检测图；

图9为本发明具体实施方式中的新数据2的T²统计量检测图和SPE统计量检测图；

其中，(a)为新数据2的T²统计量检测图，(b)为新数据2的SPE统计量检测图；

图10为本发明具体实施方式中的针对第l类故障T²统计量重构后的新数据2的T²统计量检测图和针对第l类故障SPE统计量重构后的新数据2的SPE统计量检测图；

其中，(a)为第l＝1类故障T²统计量重构后的新数据2的T²统计量检测图；

(b)为针对第l＝1类故障SPE统计量重构后的新数据2的SPE统计量检测图；

(c)为针对第l＝2类故障T²统计量重构后的新数据2的T²统计量检测图；

(d)为针对第l＝2类故障SPE统计量重构后的新数据2的SPE统计量检测图；

图11为本发明具体实施方式中的新数据2去除与故障无关的第l＝1类故障信息后的T²统计量检测图和SPE统计量检测图；

其中，(a)为新数据2去除与故障无关的第l＝1类故障信息后的T²统计量检测图；

(b)为新数据2去除与故障无关的第l＝1类故障信息后的SPE统计量检测图。

具体实施方式

下面结合附图对本发明具体实施方式加以详细的说明。

考虑到实际工业生产过程中变量之间通常呈现非线性关系，本发明在传统方法的基础上，进一步分析了故障数据与正常数据的关系，进而提取出与故障相关的故障特征方向。提出了基于KPCA的故障重构方法，该方法首先采用KPCA方法将故障数据空间分解为主元子空间和残差子空间，利用所得的负载方向对正常数据进行投影，利用PCA方法在高维特征空间中对投影的数据进行分析，并且通过比较各个方向上故障数据与正常数据的得分提取出引起统计量超限的故障方向，通过得到的故障特征方向对故障测试数据进行重构，消除统计量超限的现象，利用本文提出的方法对电熔镁熔炼过程的故障进行检测及故障诊断可以验证所提方法的有效性。

电熔镁砂是一种广泛应用于化学、航天、冶金等领域的重要耐火材料，电熔镁炉是用于生产电熔镁砂的主要设备之一。电熔镁炉是一种以电弧为热源的熔炼炉，它的热量集中，有利于熔炼电熔镁砂。电熔镁炉的整体设备组成主要包括：变压器、电路短网、电极、电极升降装置以及炉体等。炉边设有控制室，可控制电极升降。电熔镁炉的基本工作原理示意图如图1所示。

电熔镁炉通过电极引入大电流形成弧光产生高温来完成熔炼过程。目前我国多数电熔镁炉冶炼过程自动化程度还比较低，往往导致故障频繁和异常情况时有发生，其中由于电极执行器故障等原因导致电极距离电熔镁炉的炉壁过近，使得炉温异常，可以导致电熔镁炉的炉体熔化，熔炉一旦发生将会导致大量的财产损失以及危害人身安全。另外，由于炉体固定，执行器异常等原因导致电极长时间位置不变造成炉温不均，造成距离电极附近温度高，而距离电极远的区域温度低，一旦电极附近区域温度过高，容易造成“烧飞”炉料；而远离电极的区域温度过低形成死料区，这将严重影响产品产量和质量。这就需要及时地检测过程中异常和故障，因此，对电熔镁炉工作过程进行过程监测是十分必要和有意义的。

本实施方式中，采用本发明基于KPCA的工业过程故障诊断方法对电熔镁炉冶炼过程进行故障诊断，如图2所示，包括以下步骤：

步骤1：对电熔镁炉冶炼过程的历史正常数据利用KPCA方法进行子空间划分，得到历史正常数据的主元子空间的负载方向P和历史正常数据的残差子空间的负载方向P^*。

步骤1.1：对电熔镁炉冶炼过程的n个历史正常数据进行m次采样，获取电熔镁炉冶炼过程的历史正常数据集X。

本实施方式中，历史正常数据类型包括电流、电压、电极之间的相对位子和电熔镁炉的温度4种，即n为4，采样次数m为200。

步骤1.2：对历史正常数据集X进行中心化处理后将其映射到高维特征空间，得到高维特征空间的历史正常数据Φ(X)＝[Φ(x₁)，Φ(x₂)，…，Φ(x_n)]。

步骤1.3：对高维特征空间的历史正常数据Φ(X)进行PCA分解，分解为历史正常数据的主元子空间和历史正常数据的残差子空间，得到历史正常数据的主元部分和历史正常数据的残差部分

本实施方式中，定义核矩阵K∈R^n×n，[K]_i，j＝<Φ(x_i)，Φ(x_j)>，其中，i，j∈n，对高维特征空间的历史正常数据Φ(X)进行PCA分解，如式(1)所示：

$>\Phi \left(X\right)=\hat{\Phi }\left(X\right)+\tilde{\Phi }\left(X\right)---\left(1\right)$>

其中，为历史正常数据的主元部分，为历史正常数据的残差部分。

步骤1.4：计算历史正常数据的主元子空间的负载方向P和历史正常数据的残差子空间的负载方向P^*。

本实施方式中，历史正常数据的主元子空间的负载方向P，如式(2)所示：

P＝Φ(X)A    (2)

其中，A∈R^n×d为历史正常数据的主元子空间对应的核矩阵K的的负载矩阵，即前d个主元负载，d为常数，d∈1…n，本实施方式中，d取值为2。

历史正常数据的残差子空间的负载方向P^*，如式(3)所示：

P^*＝Φ(X)A^*    (3)

其中，A^*∈R^n×(n-d)为历史正常数据的残差子空间对应的核矩阵K的负载矩阵。

步骤2：对电熔镁炉冶炼过程的已知故障的历史故障数据利用KPCA方法进行子空间划分，得到历史故障数据的主元子空间的负载方向P_f和历史故障数据的残差子空间的负载方向$>P_f^{*}.$>

步骤2.1：对工业生产过程的n个历史故障数据进行m次采样，获取工业生产过程的已知故障的历史故障数据集X_f。

步骤2.2：对历史故障数据集X_f进行中心化处理后将其映射到高维特征空间，得到高维特征空间的历史故障数据Φ(X_f)＝[Φ(x_f，1)，Φ(x_f，2)，…，Φ(x_f，n)]。

步骤2.3：对高维特征空间的历史故障数据Φ(X_f)进行PCA分解，分解为历史故障数据的主元子空间和历史故障数据的残差子空间，得到历史故障数据的主元部分和历史正常数据的残差部分

本实施方式中，定义核矩阵K_f∈R^n×n，[K_f]_i，j＝<Φ(x_f，i)，Φ(x_f，i)>，其中，i，j∈n，对高维特征空间的历史故障数据Φ(X_f)进行PCA分解，如式(4)所示：

$>\Phi \left(X_f\right)=\hat{\Phi }\left(X_f\right)+\tilde{\Phi }\left(X_f\right)---\left(4\right)$>

其中，为历史故障数据的主元部分，历史正常数据的残差部分。

步骤2.4：计算历史故障数据的主元子空间的负载方向P_f和历史故障数据的残差子空间的负载方向

本实施方式中，历史故障数据的主元子空间的负载方向P_f，如式(5)所示：

P_f＝Φ(X_f)A_f    (5)

其中，A_f为历史故障数据的主元子空间对应的核矩阵K_f的负载矩阵。

历史故障数据的残差子空间的负载方向如式(6)所示：

$>P_f^{*}=\Phi \left(X_f\right)A_f^{*}---\left(6\right)$>

其中，为历史故障数据的残差子空间对应的核矩阵K_f的负载矩阵。

步骤3：对历史故障数据进行基于T²统计量重构，提取出历史故障数据针对T²统计量重构的故障特征方向如图3所示。

步骤3.1：将历史正常数据的主元部分映射到历史故障数据的主元子空间，获得针对T²统计量重构的历史故障数据的主元子空间中与故障相关的特征方向P_fr。

步骤3.1.1：将历史正常数据的主元部分映射到历史故障数据的主元子空间，得到历史正常数据的主元部分映射到历史故障数据的主元子空间的映射矩阵

本实施方式中，历史正常数据的主元部分映射到历史故障数据的主元子空间的映射矩阵如式(7)所示：

$>\left(\begin{array}{c}{\hat{\Phi }}_f^T\left(X\right)={\hat{\Phi }}^T\left(X\right)P_f{P}_f^T\\ ={\Phi }^T\left(X\right)\Phi \left(X\right){\mathrm{AA}}^T{\Phi }^T\left(X\right)\Phi \left(X_f\right)A_f{A}_f^T{\Phi }^T\left(X_f\right)\\ ={\mathrm{KAA}}^T{K}_{\alpha }{A}_f{A}_f^T{\Phi }^T\left(X_f\right)\end{array}\right)---\left(7\right)$>

其中，核矩阵K_α定义为[K_α]_i，j＝<Φ(x_i)，Φ(x_f，j)>。

步骤3.1.2：对历史正常数据的主元部分映射到历史故障数据的主元子空间的映射矩阵进行PCA分解，得到映射矩阵在历史故障数据的主元子空间的负载方向P_r。

本实施方式中，映射矩阵在历史故障数据的主元子空间的负载方向P_r，如式(8)所示：

$>P_r={\hat{\Phi }}_f\left(X\right)A_r---\left(8\right)$>

其中，A_r为映射矩阵对应核矩阵K_α的负载矩阵。

步骤3.1.3：计算映射矩阵在其历史故障数据的主元子空间的负载方向P_r上的得分T_r，计算历史故障数据的主元部分在映射矩阵的历史故障数据的主元子空间的负载方向P_r上的得分T_fr。

本实施方式中，映射矩阵在其历史故障数据的主元子空间的负载方向P_r上的得分T_r，如式(9)所示：

$>T_r={\hat{\Phi }}_f^T\left(X\right)P_r---\left(9\right)$>

历史故障数据的主元部分在映射矩阵的历史故障数据的主元子空间的负载方向P_r上的得分T_fr，如式(10)所示：

$>T_{\mathrm{fr}}={\hat{\Phi }}^T\left(X_f\right)P_r---\left(10\right)$>

步骤3.1.4：定义在映射矩阵的历史故障数据的主元子空间的负载方向P_r上故障数据与正常数据的得分比矩阵RT，其中，μ＝1，2，…，d，RT_μ为得分比矩阵RT的第μ个得分比，var(·)表示方差运算，T_fr(：，μ)为得分T_fr的第μ列，T_r(：，μ)为得分T_r的第μ列。

步骤3.1.5：设定得分比率阈值η，将得分比矩阵RT中的得分比与得分比率阈值η对比，在映射矩阵在历史故障数据的主元子空间的负载方向P_r中提取出得分比大于得分比率阈值η的负载方向，得到针对T²统计量重构的历史故障数据的主元子空间中与故障相关的特征方向P_fr。

本实施方式中，设定得分比率阈值η≥1，针对T²统计量重构的历史故障数据的主元子空间中与故障相关的特征方向P_fr，表示如式(11)所示：

$>P_{\mathrm{fr}}={\hat{\Phi }}_f\left(X\right)A_{\mathrm{fr}}---\left(11\right)$>

其中，A_fr为负载矩阵A_r的列向量中得分比大于得分比率阈值η对应的列向量矩阵。

步骤3.2：将历史正常数据的主元部分映射到历史故障数据的残差子空间，获得针对T²统计量重构的历史故障数据的残差子空间中与故障相关的特征方向

步骤3.2.1：将历史正常数据的主元部分映射到历史故障数据的残差子空间，得到历史正常数据的主元部分映射到历史故障数据的残差子空间的映射矩阵

本实施方式中，历史正常数据的主元部分映射到历史故障数据的残差子空间的映射矩阵如式(12)所示：

$>{\hat{\Phi }}_f^{*T}\left(X\right)={\hat{\Phi }}^T\left(X\right)P_f^{*}{P}_f^{*T}---\left(12\right)$>

步骤3.2.2：对历史正常数据的主元部分映射到历史故障数据的残差子空间的映射矩阵进行PCA分解，得到映射矩阵在历史故障数据的残差子空间的负载方向

本实施方式中，映射矩阵在历史故障数据的残差子空间的负载方向如式(13)所示：

$>P_r^{*}={\hat{\Phi }}_f^{*}\left(X\right)A_r^{*}---\left(13\right)$>

其中，为映射矩阵对应核矩阵K_α的负载矩阵。

步骤3.2.3：计算映射矩阵在其历史故障数据的残差子空间的负载方向上的得分计算历史故障数据的残差部分在映射矩阵的历史故障数据的残差子空间的负载方向上的得分

本实施方式中，映射矩阵在其历史故障数据的残差子空间的负载方向上的得分如式(14)所示：

$>T_r^{*}={\hat{\Phi }}_f^{*T}\left(X\right)P_r^{*}---\left(14\right)$>

历史故障数据的残差部分在映射矩阵的历史故障数据的残差子空间的负载方向上的得分如式(15)所示：

$>T_{\mathrm{fr}}^{*}={\hat{\Phi }}_f^{*T}\left(X\right)P_r^{*}---\left(15\right)$>

步骤3.2.4：定义在映射矩阵的历史故障数据的残差子空间的负载方向上故障数据与正常数据的得分比矩阵RT^*，其中，为得分比矩阵RT^*的第μ个得分比，为得分的第μ列，为得分的第μ列。

步骤3.2.5：将得分比矩阵RT^*中的得分比与得分比率阈值η对比，在映射矩阵在历史故障数据的残差子空间的负载方向中提取出得分比大于得分比率阈值η的负载方向，得到针对T²统计量重构的历史故障数据的残差子空间中与故障相关的特征方向

本实施方式中，针对T²统计量重构的历史故障数据的残差子空间中与故障相关的特征方向如式(16)所示：

$>P_{\mathrm{fr}}^{*}={\hat{\Phi }}_f^{*}\left(X\right)A_{\mathrm{fr}}^{*}---\left(16\right)$>

其中，为负载矩阵的列向量中得分比大于得分比率阈值η对应的列向量矩阵。

步骤3.3：获得历史故障数据针对T²统计量重构的故障特征方向P_rec。

本实施方式中，历史故障数据针对T²统计量重构的故障特征方向P_rec如式(17)所示：

$>P_{\mathrm{rec}}=[P_{\mathrm{fr}},P_{\mathrm{fr}}^{*}]=[{\hat{\Phi }}_f\left(X\right)A_{\mathrm{fr}},{\hat{\Phi }}_f^{*}\left(X\right)A_{\mathrm{fr}}^{*}]---\left(17\right)$>

步骤4：对历史故障数据进行基于SPE统计量重构，提取出历史故障数据针对SPE统计量重构的故障特征方向如图4所示。

步骤4.1：将历史正常数据的残差部分映射到历史故障数据的主元子空间，获得针对SPE统计量重构的历史故障数据的主元子空间中与故障相关的特征方向V_fr。

步骤4.1.1：将历史正常数据的残差部分映射到历史故障数据的主元子空间，得到历史正常数据的残差部分映射到历史故障数据的主元子空间的映射矩阵

本实施方式中，历史正常数据的残差部分映射到历史故障数据的主元子空间的映射矩阵如式(18)所示：

$>\left(\begin{array}{c}{\tilde{\Phi }}_f^T\left(X\right)={\tilde{\Phi }}^T\left(X\right)P_f{P}_f^T\\ ={\Phi }^T\left(X\right)\Phi \left(X\right){A^{*}A}^{*T}{\Phi }^T\left(X\right)\Phi \left(X_f\right)A_f{A}_f^T{\Phi }^T\left(X_f\right)\\ ={K{A}^{*}A}^{*T}{K}_{\alpha }{A}_f{A}_f^T{\Phi }^T\left(X_f\right)\end{array}\right)---\left(18\right)$>

步骤4.1.2：对历史正常数据的残差部分映射到历史故障数据的主元子空间的映射矩阵进行PCA分解，得到映射矩阵在历史故障数据的主元子空间的负载方向V_r。

本实施方式中，映射矩阵在历史故障数据的主元子空间的负载方向V_r，如式(19)所示：

$>V_r={\tilde{\Phi }}_f\left(X\right)A_e---\left(19\right)$>

其中，A_e为映射矩阵对应核矩阵K_α的负载矩阵。

步骤4.1.3：定义在映射矩阵的历史故障数据的主元子空间的负载方向V_r上故障数据与正常数据的平方误差比矩阵RT_e，其中，k＝1，2，…，d，RT_e，k为平方误差比矩阵RT_e的第k个平方误差比，||·||表示欧氏距离，V_r，k为负载方向V_r的第k列。

步骤4.1.4：设定平方误差比率阈值θ，将平方误差比矩阵RT_e中的平方误差比与平方误差比率阈值θ对比，在映射矩阵在历史故障数据的主元子空间的负载方向V_r中提取出平方误差比大于平方误差比率阈值θ的负载方向，得到针对SPE统计量重构的历史故障数据的主元子空间中与故障相关的特征方向V_fr。

本实施方式中，设定平方误差比率阈值θ≥1，针对SPE统计量重构的历史故障数据的主元子空间中与故障相关的特征方向V_fr，如式(20)所示：

$>V_{\mathrm{fr}}={\tilde{\Phi }}_f\left(X\right)A_{e,\mathrm{fr}}---\left(20\right)$>

其中，A_e，fr为负载矩阵A_e的列向量中平方误差比大于平方误差比率阈值θ对应的列向量矩阵。

步骤4.2：将历史正常数据的残差部分映射到历史故障数据的残差子空间，获得针对SPE统计量重构的历史故障数据的残差子空间中与故障相关的特征方向

步骤4.2.1：将历史正常数据的残差部分映射到历史故障数据的残差子空间，得到历史正常数据的残差部分映射到历史故障数据的残差子空间的映射矩阵

本实施方式中，历史正常数据的残差部分映射到历史故障数据的残差子空间的映射矩阵如式(21)所示：

$>{\tilde{\Phi }}_f^{*T}\left(X\right)={\tilde{\Phi }}^T\left(X\right)P_f^{*}{P}_f^{*T}---\left(21\right)$>

步骤4.2.2：对历史正常数据的残差部分映射到历史故障数据的残差子空间的映射矩阵进行PCA分解，得到映射矩阵在历史故障数据的残差子空间的负载方向

本实施方式中，映射矩阵在历史故障数据的残差子空间的负载方向如式(22)所示：

$>V_r^{*}={\tilde{\Phi }}_f\left(X\right)A_e^{*}---\left(22\right)$>

其中，为映射矩阵对应核矩阵K_α的负载矩阵。

步骤4.2.3：定义在映射矩阵的历史故障数据的残差子空间的负载方向上故障数据与正常数据的平方误差比矩阵其中，为平方误差比矩阵的第k个平方误差比，为负载方向的第k列。

步骤4.2.4：将平方误差比矩阵中的平方误差比与平方误差比率阈值θ对比，在映射矩阵在历史故障数据的残差子空间的负载方向中提取出平方误差比大于平方误差比率阈值θ的负载方向，得到针对SPE统计量重构的历史故障数据的残差子空间中与故障相关的特征方向

本实施方式中，针对SPE统计量重构的历史故障数据的残差子空间中与故障相关的特征方向如式(23)所示：

$>V_{\mathrm{fr}}^{*}={\tilde{\Phi }}_f\left(X\right)A_{e,\mathrm{fr}}^{*}---\left(23\right)$>

其中，为负载矩阵的列向量中平方误差比大于平方误差比率阈值θ对应的列向量矩阵。

步骤4.3：获得历史故障数据针对SPE统计量重构的故障特征方向V_rec。

本实施方式中，历史故障数据针对SPE统计量重构的故障特征方向V_rec如式(24)所示：

$>V_{\mathrm{rec}}=[V_{\mathrm{fr}},V_{\mathrm{fr}}^{*}]=[{\tilde{\Phi }}_f\left(X\right)A_{e,\mathrm{fr}},{\tilde{\Phi }}_f^{*}\left(X\right)A_{e,\mathrm{fr}}^{*}]---\left(24\right)$>

步骤5：重复步骤1-步骤4，对c类已知故障类型的历史故障数据进行基于T²统计量重构和基于SPE统计量重构，得到重构故障特征方向集合。

本实施方式中，选取c＝2类已知故障类型的历史故障数据，包括阶跃故障和斜坡故障。

对电熔镁炉冶炼过程的新数据进行故障诊断的过程，如图5所示。

步骤6：实时采集电熔镁炉冶炼过程的新数据，利用KPCA方法进行子空间划分，计算新数据的T²统计量和SPE统计量，并根据新数据的T²统计量和SPE统计量判断实时采集的电熔镁炉冶炼过程是否发生故障。

步骤6.1：实时采集电熔镁炉冶炼过程的新数据集X_new，对其进行中心化处理后将其映射到高维特征空间，得到高维特征空间的新数据Φ(X_new)。

本实施方式中，采集的新数据为400个。

步骤6.2：对高维特征空间的新数据Φ(X_new)进行PCA分解，得到新数据的得分t_new＝<P，Φ(X_new)>和新数据的残差e_new＝(I-PP^T)Φ(X_new)，其中，I为单位矩阵。

步骤6.3：计算新数据的T²统计量和新数据的SPE统计量SPE_new。

本实施方式中，计算的新数据的T²统计量如式(25)所示：

$>T_{\mathrm{new}}^2=t_{\mathrm{new}}{\Lambda }^{-1}{t}_{\mathrm{new}}^T---\left(25\right)$>

新数据的SPE统计量SPE_new，如式(26)所示：

$>{\mathrm{SPE}}_{\mathrm{new}}=e_{\mathrm{new}}^T{e}_{\mathrm{new}}---\left(26\right)$>

步骤6.4：判断当前电熔镁炉冶炼过程是否发生故障：若新数据的T²统计量或SPE统计量偏离其对应的正常值域，执行步骤7，否则返回步骤6.1。

步骤7：利用重构故障特征方向集合对实时采集的电熔镁炉冶炼过程的新数据进行故障方向重构，判断当前电熔镁炉冶炼过程的故障类型。

步骤7.1：选取故障类型l，利用第l类历史故障数据针对T²统计量重构的故障特征方向对新数据进行故障方向重构。

利用第l类历史故障数据针对T²统计量重构的故障特征方向对新数据进行故障方向重构公式如式(27)所示：

$>{\Phi }_{T^2,f,l}^T\left(X_{\mathrm{new}}\right)={\Phi }^T\left(X_{\mathrm{new}}\right)(I-P_{\mathrm{rec},l}{P}_{\mathrm{rec},l}^T)---\left(27\right)$>

其中，为针对第l类故障T²统计量重构后的新数据，P_rec，l为第l类历史故障数据针对T²统计量重构的故障特征方向，l＝1，2。

步骤7.2：利用第l类历史故障数据针对SPE统计量重构的故障特征方向对新数据进行故障方向重构。

利用第l类历史故障数据针对SPE统计量重构的故障特征方向对新数据进行故障方向重构公式如式(28)所示：

$>{\Phi }_{\mathrm{SPE},f,l}^T\left(X_{\mathrm{new}}\right)={\Phi }^T\left(X_{\mathrm{new}}\right)(I-V_{\mathrm{rec},l}{V}_{\mathrm{rec},l}^T)---\left(28\right)$>

其中，Φ_SPE，f，l(X_new)为针对第l类故障SPE统计量重构后的新数据，V_rec，l为第l类历史故障数据针对SPE统计量重构的故障特征方向。

步骤7.3：计算针对第l类故障T²统计量重构后的新数据的得分：计算针对第l类故障SPE统计量重构后的新数据Φ_SPE，f，l(X_new)的残差e_f，l＝(I-PP^T)Φ_SPE，f，l(X_new)。

步骤7.4：计算针对第l类故障T²统计量重构后的新数据的T²统计量和针对第l类故障SPE统计量重构后的新数据Φ_SPE，f，l(X_new)的SPE统计量$>{\mathrm{SPE}}_{f,l}=e_{f,l}^T{e}_{f,l}.$>

步骤7.5：判断当前电熔镁炉冶炼过程的故障类型：若针对第l类故障T²统计量重构后的新数据的T²统计量和针对第l类故障SPE统计量重构后的新数据Φ_SPE，f，l(X_new)的SPE统计量SPE_f，l均未偏离其对应的正常值域，则当前电熔镁炉冶炼过程的故障类型为第l类故障，否则，重新选择故障类型l，返回步骤7.1。

实施例1：

实时采集电熔镁炉冶炼过程的400个新数据1，新数据1的T²统计量检测图和SPE统计量检测图如图6所示，其中，图6(a)为新数据的T²统计量检测图，图6(b)为新数据的SPE统计量检测图，从图中可见，新数据1的T²统计量和SPE统计量都大约从第101个采样点开始出现了超限现象，并形成稳定的报警，提示故障发生。

对新数据1利用第l类历史故障数据针对T²统计量重构的故障特征方向对新数据进行故障方向重构，利用第l类历史故障数据针对SPE统计量重构的故障特征方向对新数据进行故障方向重构。

针对第l类故障T²统计量重构后的新数据1的T²统计量检测图和针对第l类故障SPE统计量重构后的新数据1的SPE统计量检测图如图7所示，其中，图7(a)为第l＝1类故障T²统计量重构后的新数据1的T²统计量检测图；图7(b)为针对第l＝1类故障SPE统计量重构后的新数据1的SPE统计量检测图；图7(c)为针对第l＝2类故障T²统计量重构后的新数据1的T²统计量检测图；图7(d)为针对第l＝2类故障SPE统计量重构后的新数据1的SPE统计量检测图。

由上述实验结果可知，第l＝1类故障的重构能够去除新数据1中的故障信息，消除统计量中的超限现象。而第l＝2类故障的重构模型无法完全去除新数据1中的故障信息，无法消除统计量超限的现象。据此，可判断该新数据1中含有第l＝1类故障的信息，当前发生的故障为第l＝1类故障。

对该结果进行检验，新数据1去除与故障无关的第l＝2类故障信息后的T²统计量检测图和SPE统计量检测图如图8所示，其中，图8(a)为新数据1去除与故障无关的第l＝2类故障信息后的T²统计量检测图；图8(b)为新数据1去除与故障无关的第l＝2类故障信息后的SPE统计量检测图；由图可知，均为能消除统计量中的超限现象，说明该结果正确。

实施例2：

实时采集电熔镁炉冶炼过程的400个新数据2，新数据2的T²统计量检测图和SPE统计量检测图如图9所示，其中，图9(a)为新数据的T²统计量检测图，图9(b)为新数据的SPE统计量检测图，从图中可见，新数据2的T²统计量和SPE统计量均从第150个采样附近开始出现了超限现象，并形成稳定的报警，提示故障发生。

对新数据2利用第l类历史故障数据针对T²统计量重构的故障特征方向对新数据进行故障方向重构，利用第l类历史故障数据针对SPE统计量重构的故障特征方向对新数据进行故障方向重构。

针对第l类故障T²统计量重构后的新数据2的T²统计量检测图和针对第l类故障SPE统计量重构后的新数据2的SPE统计量检测图如图10所示，其中，图10(a)为第l＝1类故障T²统计量重构后的新数据2的T²统计量检测图；图10(b)为针对第l＝1类故障SPE统计量重构后的新数据2的SPE统计量检测图；图10(c)为针对第l＝2类故障T²统计量重构后的新数据2的T²统计量检测图；图10(d)为针对第l＝2类故障SPE统计量重构后的新数据2的SPE统计量检测图。

由上述实验结果可知，第l＝1类故障的重构能够去除新数据2中的故障信息，仍存在超限报警现象。而第l＝2类故障的重构模型无法完全去除新数据2中的故障信息，能够完全消除超限现象。据此，可判断该新数据2中含有第l＝2类故障的信息，当前发生的故障为第l＝2类故障。

对该结果进行检验，新数据2去除与故障无关的第l＝1类故障信息后的T²统计量检测图和SPE统计量检测图如图11所示，其中，图11(a)为新数据2去除与故障无关的第l＝1类故障信息后的T²统计量检测图，图11(b)为新数据2去除与故障无关的第l＝1类故障信息后的SPE统计量检测图，由图可知，均为能消除统计量中的超限现象，说明该结果正确。