CN111474920A - 一种基于ica与acva的故障检测方法 - Google Patents

一种基于ica与acva的故障检测方法 Download PDF

Info

Publication number
CN111474920A
CN111474920A CN202010352906.5A CN202010352906A CN111474920A CN 111474920 A CN111474920 A CN 111474920A CN 202010352906 A CN202010352906 A CN 202010352906A CN 111474920 A CN111474920 A CN 111474920A
Authority
CN
China
Prior art keywords
matrix
ica
obtaining
acva
variable
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Granted
Application number
CN202010352906.5A
Other languages
English (en)
Other versions
CN111474920B (zh
Inventor
凡时财
张建
邹见效
徐红兵
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
University of Electronic Science and Technology of China
Original Assignee
University of Electronic Science and Technology of China
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by University of Electronic Science and Technology of China filed Critical University of Electronic Science and Technology of China
Priority to CN202010352906.5A priority Critical patent/CN111474920B/zh
Publication of CN111474920A publication Critical patent/CN111474920A/zh
Application granted granted Critical
Publication of CN111474920B publication Critical patent/CN111474920B/zh
Expired - Fee Related legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • G—PHYSICS
    • G05—CONTROLLING; REGULATING
    • G05B—CONTROL OR REGULATING SYSTEMS IN GENERAL; FUNCTIONAL ELEMENTS OF SUCH SYSTEMS; MONITORING OR TESTING ARRANGEMENTS FOR SUCH SYSTEMS OR ELEMENTS
    • G05B23/00—Testing or monitoring of control systems or parts thereof
    • G05B23/02—Electric testing or monitoring
    • G05B23/0205—Electric testing or monitoring by means of a monitoring system capable of detecting and responding to faults
    • G05B23/0259—Electric testing or monitoring by means of a monitoring system capable of detecting and responding to faults characterized by the response to fault detection
    • G05B23/0262—Confirmation of fault detection, e.g. extra checks to confirm that a failure has indeed occurred
    • G—PHYSICS
    • G05—CONTROLLING; REGULATING
    • G05B—CONTROL OR REGULATING SYSTEMS IN GENERAL; FUNCTIONAL ELEMENTS OF SUCH SYSTEMS; MONITORING OR TESTING ARRANGEMENTS FOR SUCH SYSTEMS OR ELEMENTS
    • G05B2219/00—Program-control systems
    • G05B2219/20—Pc systems
    • G05B2219/24—Pc safety
    • G05B2219/24065—Real time diagnostics

Landscapes

  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Engineering & Computer Science (AREA)
  • Automation & Control Theory (AREA)
  • Complex Calculations (AREA)

Abstract

本发明公开了一种基于ICA与ACVA的故障检测方法，能够从时序的独立成分特征中提取出时序特征，有效解决了传统的ICA方法难以提取时序特征的问题，构建处理时序故障的故障检测方法。本发明通过提取故障的时序特征，能够得到更好的检测效果，提升了故障检测能力，适用于工业生产过程的在线监控，是一种相比现有技术更优的故障检测方法。本发明采用ICA方法求取了系统的独立成分特征，然后采用ACVA方法构建了时序模型，构建了ICA‑ACVA的故障检测方法，能够对工业生产过程中产生的故障进行精准的检测。

Description

一种基于ICA与ACVA的故障检测方法
技术领域
本发明属于工业过程的故障检测技术领域，具体涉及一种基于ICA与ACVA的故障检测方法。
背景技术
现代工业生产已越来越规模化、复杂化，若生产过程出现故障，不仅影响产品质量，更容易对人的生命安全造成威胁。因此常常采用故障检测技术来监控工业生产状况。由于复杂的工业过程难以采用精确的机理建模方法，因而使用数据驱动技术进行故障检测的方法得到了广泛关注。数据驱动技术的思想是使用历史过程数据中的内在信息来判断当前生产过程的状况，从而实现故障的检测与诊断。由于过程数据中常常包含大量非高斯信息，因而基于独立成分分析(Independent Component Analysis，ICA)方法的数据驱动方法得到了很多研究。ICA方法主要利用数据的高阶信息实现独立成分特征的提取，然后基于这些独立特征计算相应的统计量并估计控制限，便得到故障检测模型，并通过判断新数据的统计量是否超过控制限来实现故障检测。但是在实际工业过程中，故障在随时间演变的过程，具有与正常工况不同的时序特征，所以故障均是时序故障。然而，ICA方法虽利用了多变量的相关性，却忽略了变量的自相关性，故对于多数时序故障的检测效果难以提高，尤其是对一些变换缓慢的时序故障，其检测效果比较差。本发明目的在于提高对时序故障的检测能力，提供一种基于独立成分分析(Independent Component Analysis，ICA)方法和自适应规范变量分析(Adaptive Canonical Variate Analysis，ACVA)方法结合的ICA-ACVA方法的故障检测方法，用于工业过程的故障检测。
发明内容
针对现有技术中的上述不足，本发明提供的一种基于ICA与ACVA的故障检测方法，包括以下步骤：
S1、采集若干工业生产系统的工作信号作为样本信号，并对样本信号进行标准化，得到标准化样本矩阵；
S2、根据标准化后的工作信号，使用快速ICA算法求解ICA模型，得到解混矩阵；
S3、使用解混矩阵对ACVA模型求解，分别得到规范空间和规范残差空间的右奇异矩阵；
S4、根据规范空间和规范残差空间的右奇异矩阵，获取规范变量和规范残差变量；
S5、分别利用规范变量和规范残差变量获取第一统计量，并根据第一统计量预估控制限；
S6、采集tnew时刻工业生产系统的工作信号，使用步骤S1-S4的方法获取tnew时刻及其过去时刻工作信号的第二统计量；
S7、判断第二统计量是否超出其控制限，若是，则判定为故障状态，否则判定为正常状态。
进一步地，所述步骤S1包括以下分步骤：
S1.1、采集n个时刻的工作信号作为样本信号，得到n个m维的测量向量；
S1.2、将n个测量向量按照采样时刻进行升序排列，得到数据矩阵X＝[x(1),x(2),...,x(r),...,x(n)]，x(r)表示第r时刻得到的测量向量，r＝1,2,...,n；
S1.3、对数据矩阵X中测量向量进行标准化，得到标准化的样本矩阵
Figure BDA0002472465110000021
Figure BDA0002472465110000022
表示第r时刻标准化后的测量向量。
进一步地，所述步骤S2包括以下分步骤：
S2.1、通过PAC白化算法对标准化样本矩阵
Figure BDA0002472465110000031
进行白化处理，得到白化变换矩阵Q；
S2.2、根据白化变换矩阵Q，获取标准化样本矩阵
Figure BDA0002472465110000032
的白化矩阵Z，所述白化矩阵Z具体为：
Figure BDA0002472465110000033
S2.3、根据快速ICA算法对ICA模型Sz＝BZ进行求解，得到白化矩阵Z的第一解混矩阵B，Sz表示将Z经过解混后的独立向量矩阵；
S2.4、根据第一解混矩阵B，获取标准化样本矩阵
Figure BDA0002472465110000034
的第二解混矩阵W为：
W＝BTQ
其中，T表示转置。
进一步地，所述步骤S3包括以下分步骤：
S3.1、设定分量数d，将第二解混矩阵W的行向量按照范数从大到小的顺序依次选取d个行向量组成矩阵，得到组合矩阵Wd；
S3.2、使用组合矩阵Wd获取标准化样本矩阵
Figure BDA0002472465110000035
中第r时刻的测量向量
Figure BDA0002472465110000036
的独立向量s(r)，得到独立向量矩阵S＝[s(1),s(2),...,s(r),...,s(n)]，所述独立向量s(r)为：
Figure BDA0002472465110000037
S3.3、设定时滞长度q，根据独立向量矩阵S获取第t时刻的过去样本向量sp(t)和未来样本向量sf(t)为：
Figure BDA0002472465110000038
Figure BDA0002472465110000041
其中，q+1≤t≤n-q+1，Rdp表示d×q行1列的实数矩；
S3.4、根据过去样本向量sp(t)和未来样本向量sf(t)，获取过去Hankel矩阵Sp和未来Hankel矩阵Sf为：
Sp＝[sp(q+1),sp(q+2),...,sp(t),...,sp(q+N)]∈Rdq×N
Sf＝[sf(q+1),sf(q+2),...,sf(t),...,sf(q+N)]∈Rdq×N
其中，N＝n-2q+1，Rdq×N表示dq行N列的实数矩阵；
S3.5、分别对过去Hankel矩阵Sp和未来Hankel矩阵Sf进行标准化，分别得到标准化过去Hankel矩阵
Figure BDA0002472465110000042
和标准化未来Hankel矩阵
Figure BDA0002472465110000043
为：
Figure BDA0002472465110000044
Figure BDA0002472465110000045
S3.6、根据标准化过去Hankel矩阵
Figure BDA0002472465110000046
和标准化未来Hankel矩阵
Figure BDA0002472465110000047
获取过去Hankel矩阵协方差∑pp、未来Hankel矩阵协方差∑ff和互协方差∑fp为：
Figure BDA0002472465110000048
Figure BDA0002472465110000049
Figure BDA00024724651100000410
其中，T表示转置；
S3.7、根据过去Hankel矩阵协方差∑pp、未来Hankel矩阵协方差∑ff和互协方差∑fp获取ACVA模型H为：
Figure BDA00024724651100000411
其中，
Figure BDA0002472465110000051
表示对∑ff求取逆矩阵后再进行矩阵开方计算，
Figure BDA0002472465110000052
表示对∑pp求取逆矩阵后再进行矩阵开方计算；
S3.8、对ACVA模型H取奇异值分解，得到H＝UΛVT，U＝[u1,u2,...,udq]，V＝[v1,v2,...,vdq]，Λ＝diag(σ1,σ2,...,σdq)，U表示左奇异矩阵，u1,u2,...,udq均表示左奇异矩阵U的列向量，Λ表示奇异值矩阵，diag表示对角矩阵，σ1,σ2,...,σdq均表示奇异值，V表示右奇异矩阵，v1,v2,...,vdq均表示右奇异矩阵的列向量；
S3.9、设定阈值threshold，按奇异值从大到小的顺序依次选择奇异值，并判断选择奇异值的和是否大于或等于threshold，若是，则停止选择，记录选择的奇异值个数g，并进入步骤S3.10，否则，重复步骤S3.9；
S3.10、选取右奇异矩阵V中前g个列向量，得到规范空间的右奇异矩阵Vd，并将右奇异矩阵V中剩余的列向量组成规范残差空间的右奇异矩阵Ve。
进一步地，所述步骤S4中获取规范变量和规范残差变量的具体方法为：通过规范空间的右奇异矩阵Vd和规范残差空间的右奇异矩阵Ve获取第t时刻的规范变量cd(t)和规范残差变量ce(t)为：
Figure BDA0002472465110000053
Figure BDA0002472465110000054
进一步地，所述步骤S5中分别利用规范变量和规范残差变量获取第一统计量的具体方法为：根据规范变量cd(t)和规范残差变量ce(t)，获取规范变量cd(t)和规范残差变量ce(t)的统计量为：
Figure BDA0002472465110000055
Figure BDA0002472465110000056
其中，
Figure BDA0002472465110000057
表示与规范变量cd(t)对应的统计量，
Figure BDA0002472465110000058
表示与规范残差变量ce(t)对应的统计量，将
Figure BDA0002472465110000059
和
Figure BDA00024724651100000510
作为第一统计量。
进一步地，所述步骤S5中根据第一统计量预估控制限的具体方法为：
B1、使用核密度估计法KDE分别估计统计量
Figure BDA0002472465110000061
的概率密度
Figure BDA0002472465110000062
和统计量
Figure BDA0002472465110000063
的概率密度
Figure BDA0002472465110000064
B2、根据统计量
Figure BDA0002472465110000065
和统计量
Figure BDA0002472465110000066
获取数据跨度diff为diff＝max(data)-min(data)
其中，data取值为
Figure BDA0002472465110000067
或
Figure BDA0002472465110000068
B3、根据数据跨度diff选取控制限求解区间[start,end]，start＝min(data)-diff/2，end＝max(data)+diff/2；
B4、设定小区间分割个数为c和控制限的置信度α，将控制限求解区间[start,end]分割为c个小区间；
B5、根据概率密度
Figure BDA0002472465110000069
和概率密度
Figure BDA00024724651100000610
在控制限求解区间[start,end]上从左至右叠加小区间的面积值，并判断叠加的面积值是否大于或等于c个小区间总面积的α，若是，则得到对应的叠加区间数k，并进入步骤B6，否则重复步骤B5；
B6、通过叠加区间数k获取控制限UCLα为start+k△，△表示每个小区间的区间长度，△＝(end-start)/c。
进一步地，所述步骤B5中k值的求取公式为：
Figure BDA00024724651100000611
其中，j取值为d或e，
Figure BDA00024724651100000612
表示统计量
Figure BDA00024724651100000613
在(start+i△)处的概率密度，将其作为在第i个小区间上的高度；
Figure BDA00024724651100000614
表示统计量
Figure BDA00024724651100000615
在(start+i△)处概率密度，将其作为在第i个小区间上的高度，△表示小区间的宽度。
本发明的有益效果为：
(1)本发明能够从时序的独立成分特征中提取出时序特征，有效解决了传统的ICA方法难以提取时序特征的问题，构建处理时序故障的故障检测方法。
(2)本发明通过提取故障的时序特征，能够得到更好的检测效果，提升了故障检测能力，适用于工业生产过程的在线监控，是一种相比现有技术更优的故障检测方法。
(3)本发明采用ICA方法求取了系统的独立成分特征，然后采用ACVA方法构建了时序模型，构建了ICA-ACVA的故障检测方法，能够对工业生产过程中产生的故障进行精准的检测。
附图说明
图1为本发明提出的一种基于ICA与ACVA的故障检测方法流程图。
具体实施方式
下面对本发明的具体实施方式进行描述，以便于本技术领域的技术人员理解本发明，但应该清楚，本发明不限于具体实施方式的范围，对本技术领域的普通技术人员来讲，只要各种变化在所附的权利要求限定和确定的本发明的精神和范围内，这些变化是显而易见的，一切利用本发明构思的发明创造均在保护之列。
下面结合附图详细说明本发明的实施例。
如图1所示，一种基于ICA与ACVA的故障检测方法，包括以下步骤：
S1、采集若干工业生产系统的工作信号作为样本信号，并对样本信号进行标准化，得到标准化样本矩阵；
S2、根据标准化后的工作信号，使用快速ICA算法求解ICA模型，得到解混矩阵；
S3、使用解混矩阵对ACVA模型求解，分别得到规范空间和规范残差空间的右奇异矩阵；
S4、根据规范空间和规范残差空间的右奇异矩阵，获取规范变量和规范残差变量；
S5、分别利用规范变量和规范残差变量获取第一统计量，并根据第一统计量预估控制限；
S6、采集tnew时刻工业生产系统的工作信号，使用步骤S1-S4的方法获取tnew时刻及其过去时刻工作信号的第二统计量；
S7、判断第二统计量是否超出其控制限，若是，则判定为故障状态，否则判定为正常状态。
所述步骤S1包括以下分步骤：
S1.1、采集n个时刻的工作信号作为样本信号，得到n个m维的测量向量；
S1.2、将n个测量向量按照采样时刻进行升序排列，得到数据矩阵X＝[x(1),x(2),...,x(r),...,x(n)]，x(r)表示第r时刻得到的测量向量，r＝1,2,...,n；
S1.3、对数据矩阵X中测量向量进行标准化，得到标准化的样本矩阵
Figure BDA0002472465110000081
Figure BDA0002472465110000082
表示第r时刻标准化后的测量向量。
在本实施例中，标准化测量向量
Figure BDA0002472465110000083
中的第l个元素
Figure BDA0002472465110000084
的计算公式如下：
Figure BDA0002472465110000085
其中，x(r)(l)表示测量向量x(r)的第l个元素，mean(X(l))表示矩阵X中第l行的样本均值，std(X(l))表示矩阵X中第l行的样本标准差，l＝1,2,...,m。
所述步骤S2包括以下分步骤：
S2.1、通过PAC白化算法对标准化样本矩阵
Figure BDA0002472465110000086
进行白化处理，得到白化变换矩阵Q；
S2.2、根据白化变换矩阵Q，获取标准化样本矩阵
Figure BDA0002472465110000087
的白化矩阵Z，所述白化矩阵Z具体为：
Figure BDA0002472465110000091
S2.3、根据快速ICA算法对ICA模型Sz＝BZ进行求解，得到白化矩阵Z的第一解混矩阵B，Sz表示将Z经过解混后的独立向量矩阵；
S2.4、根据第一解混矩阵B，获取标准化样本矩阵
Figure BDA0002472465110000092
的第二解混矩阵W为：
W＝BTQ
其中，T表示转置。
所述步骤S3包括以下分步骤：
S3.1、设定分量数d，将第二解混矩阵W的行向量按照范数从大到小的顺序依次选取d个行向量组成矩阵，得到组合矩阵Wd；
在本实施例中，Wd∈Rd×m，Rd×m表示d行m列的实数矩。
S3.2、使用组合矩阵Wd获取标准化样本矩阵
Figure BDA0002472465110000093
中第r时刻的测量向量
Figure BDA0002472465110000094
的独立向量s(r)，得到独立向量矩阵S＝[s(1),s(2),...,s(r),...,s(n)]，所述独立向量s(r)为：
Figure BDA0002472465110000095
S3.3、设定时滞长度q，根据独立向量矩阵S获取第t时刻的过去样本向量sp(t)和未来样本向量sf(t)为：
Figure BDA0002472465110000096
Figure BDA0002472465110000097
其中，q+1≤t≤n-q+1，Rdp表示d×q行1列的实数矩阵；
S3.4、根据过去样本向量sp(t)和未来样本向量sf(t)，获取过去Hankel矩阵Sp和未来Hankel矩阵Sf为：
Sp＝[sp(q+1),sp(q+2),...,sp(t),...,sp(q+N)]∈Rdq×N
Sf＝[sf(q+1),sf(q+2),...,sf(t),...,sf(q+N)]∈Rdq×N
其中，N＝n-2q+1，Rdq×N表示dq行N列的实数矩阵；
S3.5、分别对过去Hankel矩阵Sp和未来Hankel矩阵Sf进行标准化，分别得到标准化过去Hankel矩阵
Figure BDA0002472465110000101
和标准化未来Hankel矩阵
Figure BDA0002472465110000102
为：
Figure BDA0002472465110000103
在本实施例中，标准化的过去样本向量
Figure BDA0002472465110000104
中的第lp个元素
Figure BDA0002472465110000105
的计算公式如下：
Figure BDA0002472465110000106
其中，mean(Sp(lp))表示过去Hankel矩阵Sp中第lp行的样本均值，sp(t)(lp)表示过去样本向量sp(t)中第lp个元素，lp＝1,2,...,dp。
标准化的过去样本向量
Figure BDA0002472465110000107
中的第lf个元素
Figure BDA0002472465110000108
的计算公式如下：
Figure BDA0002472465110000109
其中，mean(Sf(lf))表示未来Hankel矩阵Sf中第lf行的样本均值，sf(t)(lf)表示未来样本向量sf(t)中第lf个元素，lf＝1,2,...,dp。
S3.6、根据标准化过去Hankel矩阵
Figure BDA00024724651100001010
和标准化未来Hankel矩阵
Figure BDA00024724651100001011
获取过去Hankel矩阵协方差∑pp、未来Hankel矩阵协方差∑ff和互协方差∑fp为：
Figure BDA0002472465110000111
Figure BDA0002472465110000112
Figure BDA0002472465110000113
其中，T表示转置；
S3.7、根据过去Hankel矩阵协方差∑pp、未来Hankel矩阵协方差∑ff和互协方差∑fp获取ACVA模型H为：
Figure BDA0002472465110000114
其中，
Figure BDA0002472465110000115
表示对∑ff求取逆矩阵后再进行矩阵开方计算，
Figure BDA0002472465110000116
表示对∑pp求取逆矩阵后再进行矩阵开方计算；
S3.8、对ACVA模型H取奇异值分解，得到H＝UΛVT，U＝[u1,u2,...,udq]，V＝[v1,v2,...,vdq]，Λ＝diag(σ1,σ2,...,σdq)，U表示左奇异矩阵，u1,u2,...,udq均表示左奇异矩阵U的列向量，Λ表示奇异值矩阵，diag表示对角矩阵，σ1,σ2,...,σdq均表示奇异值，V表示右奇异矩阵，v1,v2,...,vdq均表示右奇异矩阵的列向量；
S3.9、设定阈值threshold，按奇异值从大到小的顺序依次选择奇异值，并判断选择奇异值的和是否大于或等于threshold，若是，则停止选择，记录选择的奇异值个数g，并进入步骤S3.10，否则，重复步骤S3.9；
在本实施例中，g值符合以下条件：
Figure BDA0002472465110000117
其中，y∈[1,2,...,dp]。
S3.10、选取右奇异矩阵V中前g个列向量，得到规范空间的右奇异矩阵Vd，并将右奇异矩阵V中剩余的列向量组成规范残差空间的右奇异矩阵Ve。
在本实施例中，V＝[Vd Ve]。
所述步骤S5中分别利用规范变量和规范残差变量获取第一统计量的具体方法为：
A1、通过规范空间的右奇异矩阵Vd和规范残差空间的右奇异矩阵Ve获取第t时刻的规范变量cd(t)和规范残差变量ce(t)为：
Figure BDA0002472465110000121
Figure BDA0002472465110000122
A2、根据规范变量cd(t)和规范残差变量ce(t)，获取规范变量cd(t)和规范残差变量ce(t)的统计量为：
Figure BDA0002472465110000123
Figure BDA0002472465110000124
其中，
Figure BDA0002472465110000125
表示与规范变量cd(t)对应的统计量，
Figure BDA0002472465110000126
表示与规范残差变量ce(t)对应的统计量，将
Figure BDA0002472465110000127
和
Figure BDA0002472465110000128
作为第一统计量。
S3.10、选取右奇异矩阵V中前g个列向量，得到规范空间的右奇异矩阵Vd，并将右奇异矩阵V中剩余的列向量组成规范残差空间的右奇异矩阵Ve。
所述步骤S4中获取规范变量和规范残差变量的具体方法为：通过规范空间的右奇异矩阵Vd和规范残差空间的右奇异矩阵Ve获取第t时刻的规范变量cd(t)和规范残差变量ce(t)为：
Figure BDA0002472465110000129
Figure BDA00024724651100001210
所述步骤S5中分别利用规范变量和规范残差变量获取第一统计量的具体方法为：根据规范变量cd(t)和规范残差变量ce(t)，获取规范变量cd(t)和规范残差变量ce(t)的统计量为：
Figure BDA0002472465110000131
Figure BDA0002472465110000132
其中，
Figure BDA0002472465110000133
表示与规范变量cd(t)对应的统计量，
Figure BDA0002472465110000134
表示与规范残差变量ce(t)对应的统计量，将
Figure BDA0002472465110000135
和
Figure BDA0002472465110000136
作为第一统计量。
所述步骤S5中根据第一统计量预估控制限的具体方法为：
B1、使用核密度估计法KDE分别估计统计量
Figure BDA0002472465110000137
的概率密度
Figure BDA0002472465110000138
和统计量
Figure BDA0002472465110000139
的概率密度
Figure BDA00024724651100001310
B2、根据统计量
Figure BDA00024724651100001311
和统计量
Figure BDA00024724651100001312
获取数据跨度diff为diff＝max(data)-min(data)
其中，data取值为
Figure BDA00024724651100001313
或
Figure BDA00024724651100001314
B3、根据数据跨度diff选取控制限求解区间[start,end]，start＝min(data)-diff/2，end＝max(data)+diff/2；
B4、设定小区间分割个数为c和控制限的置信度α，将控制限求解区间[start,end]分割为c个小区间；
B5、根据概率密度
Figure BDA00024724651100001315
和概率密度
Figure BDA00024724651100001316
在控制限求解区间[start,end]上从左至右叠加小区间的面积值，并判断叠加的面积值是否大于或等于c个小区间总面积的α，若是，则得到对应的叠加区间数k，并进入步骤B6，否则重复步骤B5；
B6、通过叠加区间数k获取控制限UCLα为start+k△，△表示每个小区间的区间长度，△＝(end-start)/c。
所述步骤B5中k值的求取公式为：
Figure BDA00024724651100001317
其中，j取值为d或e，
Figure BDA0002472465110000141
表示统计量
Figure BDA0002472465110000142
在(start+i△)处的概率密度，将其作为在第i个小区间上的高度；
Figure BDA0002472465110000143
表示统计量
Figure BDA0002472465110000144
在(start+i△)处概率密度，将其作为在第i个小区间上的高度，△表示小区间的宽度。
在本实施例中，进行故障检测的具体步骤为：
在实时工业生产系统运行过程中，对于第tnew时刻的测量数据向量xnew(tnew)，首先对测量数据向量xnew(tnew)进行标准化，得到
Figure BDA0002472465110000145
然后应用ICA变换，得到独立向量
Figure BDA0002472465110000146
将其于前q-1个时刻所计算的独立向量snew(tnew-1),....,snew(tnew-q+1)一起组成实时过去样本向量。如果前q-1个时刻的独立向量数目不够，则当前样本不执行故障检测步骤。该实时过去样本向量公式如下：
Figure BDA0002472465110000147
然后按照过去样本向量的标准化步骤，得到标准化的实时过去样本向量
Figure BDA0002472465110000148
按照下式计算
Figure BDA0002472465110000149
的实时规范变量cd,new(tnew)和实时规范残差变量ce,new(tnew)：
Figure BDA00024724651100001410
Figure BDA00024724651100001411
其中，Vd(tnew-1)和Ve(tnew-1)分别是第tnew-1时刻的规范空间和规范残差空间的右奇异矩阵。
分别通过实时规范变量cd,new(tnew)和实时规范残差变量ce,new(tnew)获取统计量
Figure BDA00024724651100001412
和统计量
Figure BDA00024724651100001413
并判断统计量
Figure BDA00024724651100001414
是否超过其控制限
Figure BDA0002472465110000151
若是，则判定为故障状态，否则下一判断步骤；下一判断步骤为：判断统计量
Figure BDA0002472465110000152
是否超过其控制限
Figure BDA0002472465110000153
若是，则判定为故障状态，否则为正常状态。统计量
Figure BDA0002472465110000154
和统计量
Figure BDA0002472465110000155
为第二统计量。
在本实施例中，在实时工业生产系统运行过程中，应对ACVA模型进行更新。ACVA模型在满足下列更新条件时才进行更新：
Figure BDA0002472465110000156
其中，Thupdate表示更新控制倍数，显然，若Thupdate＝1，则该条件含义即为只有在两个统计量均没有报告故障的时候，才更新模型。若不满足上述条件，则从该时刻算起，共Tdelay个时刻内不执行模型更新步骤。其中Tdelay表示延时更新时刻数，1≤Tdelay≤2q，若Tdelay＝1，则只有当前时刻不更新，若Tdelay＝2q，则当前时刻的采样样本不会构成在更新过去或未来样本向量中，可防止离群值样本对模型更新产生不利影响。
若满足更新条件。则将第tnew时刻计算出的独立向量snew(tnew)，与前2q-1个时刻所计算的独立向量snew(tnew-1),....,snew(tnew-2q+1)一起组成更新过去样本向量和更新未来样本向量。如果前2q-1个时刻的独立向量数目不够，则当前时刻不执行更新步骤，该更新过去样本向量sp,update(tnew)和更新未来样本向量sf,update公式如下：
Figure BDA0002472465110000157
Figure BDA0002472465110000158
然后对更新过去样本向量sp,update(tnew)和更新未来样本向量sf,update进行标准化，得到标准化的更新过去样本向量
Figure BDA0002472465110000161
和更新未来样本向量
Figure BDA0002472465110000162
计算如下中间变量：
Figure BDA0002472465110000163
Figure BDA0002472465110000164
Figure BDA0002472465110000165
Figure BDA0002472465110000166
Figure BDA0002472465110000167
Figure BDA0002472465110000168
Figure BDA0002472465110000169
Figure BDA00024724651100001610
其中，yo表示第一更新向量，yh表示第二更新向量，
Figure BDA00024724651100001611
表示第一更新参数，
Figure BDA00024724651100001612
表示第二更新参数，
Figure BDA00024724651100001613
表示第三更新参数，
Figure BDA00024724651100001614
表示第四更新参数，b＝1,2,...,dq，
Figure BDA00024724651100001615
yo,b表示第一更新向量中下标为b的元素，yf,b表示第二更新向量中下标为b的元素。β表示遗忘因子，其值被设置在0至1之间。VT(tnew-1)表示tnew-1时刻的矩阵V的转置，vb(tnew-1)表示(tnew-1)对应的右奇异矩阵V的列向量，UT(tnew-1)表示tnew-1时刻的矩阵U的转置，ub(tnew-1)表示(tnew-1)对应的左奇异矩阵U的列向量，σb-1(tnew-1)表示(tnew-1)对应的奇异值矩阵中第b-1个奇异值。
对ACVA模型中左奇异列向量、右奇异列向量和奇异值进行更新，且更新公式如下：
Figure BDA0002472465110000171
Figure BDA0002472465110000172
Figure BDA0002472465110000173
其中，
Figure BDA0002472465110000174
表示更新后的左奇异列向量，
Figure BDA0002472465110000175
表示更新后的右奇异列向量，
Figure BDA0002472465110000176
表示更新后的奇异值，T表示转置。
左奇异列向量、右奇异列向量和奇异值的更新后的标准化公式如下：
Figure BDA0002472465110000177
Figure BDA0002472465110000178
Figure BDA0002472465110000179
其中，ub(tnew)表示标准化后的左奇异列向量，vb(tnew)表示标准化后的右奇异列向量，σb(tnew)表示标准化后的奇异值，
Figure BDA00024724651100001710
表示
Figure BDA00024724651100001711
的模长，
Figure BDA00024724651100001712
表示
Figure BDA00024724651100001713
的模长，
Figure BDA00024724651100001714
表示
Figure BDA00024724651100001715
的绝对值。
将ub(tnew)、vb(tnew)和σb(tnew)代入tnew时刻的左奇异矩阵、右奇异矩阵和奇异值矩阵，得到tnew时刻的左奇异矩阵、右奇异矩阵和奇异值矩阵如下：
U(tnew)＝[u1(tnew),u2(tnew),...,ub(tnew),...,umq(tnew)]
V(tnew)＝[v1(tnew),v2(tnew),...,vb(tnew),...,vmq(tnew)]
Λ(tnew)＝diag(σ1(tnew),σ2(tnew),...,σb(tnew),...,σdq(tnew))
此时ACVA模型完成了更新。
在本实施里中，采用基于ICA的故障检测方法作为对比方法，表1是本实施例中ICA-ACVA模型的
Figure BDA00024724651100001716
统计量和
Figure BDA00024724651100001717
统计量与ICA模型的I2统计量和SPE统计量的故障检测结果召回率统计表。
表1
Figure BDA0002472465110000181
故障召回率值越高，表示对故障的检测效果越好，因此可以从表1看出，本发明方法在工业生产系统的故障检测上具有优异的检测效果。
表2是本实施例中ICA-ACVA模型的
Figure BDA0002472465110000182
统计量和
Figure BDA0002472465110000183
统计量与ICA模型的I2统计量和SPE统计量对正常工况样本的故障检测结果召回率统计表。
表2
Figure BDA0002472465110000191
正常样本的召回率值越高，表明越不容易对正常工况造成误报。从表2可以看出，本发明方法对工业生产系统的正常工况具有很好的识别效果。从正常情况的召回率和故障情况的召回率两方面看，本发明方法均具备良好的效果。
本发明的有益效果为：能够从时序的独立成分特征中提取出时序特征，有效解决了传统的ICA方法难以提取时序特征的问题，构建处理时序故障的故障检测方法。本发明通过提取故障的时序特征，能够得到更好的检测效果，提升了故障检测能力，适用于工业生产过程的在线监控，是一种相比现有技术更优的故障检测方法，本发明采用ICA方法求取了系统的独立成分特征，然后采用ACVA方法构建了时序模型，构建了ICA-ACVA的故障检测方法，能够对工业生产过程中产生的故障进行精准的检测。

Claims (8)

1.一种基于ICA与ACVA的故障检测方法，其特征在于，包括以下步骤：
S1、采集若干工业生产系统的工作信号作为样本信号，并对样本信号进行标准化，得到标准化样本矩阵；
S2、根据标准化后的工作信号，使用快速ICA算法求解ICA模型，得到解混矩阵；
S3、使用解混矩阵对ACVA模型求解，分别得到规范空间和规范残差空间的右奇异矩阵；
S4、根据规范空间和规范残差空间的右奇异矩阵，获取规范变量和规范残差变量；
S5、分别利用规范变量和规范残差变量获取第一统计量，并根据第一统计量预估控制限；
S6、采集tnew时刻工业生产系统的工作信号，使用步骤S1-S4的方法获取tnew时刻及其过去时刻工作信号的第二统计量；
S7、判断第二统计量是否超出其控制限，若是，则判定为故障状态，否则判定为正常状态。
2.根据权利要求1所述的基于ICA与ACVA的故障检测方法，其特征在于，所述步骤S1包括以下分步骤：
S1.1、采集n个时刻的工作信号作为样本信号，得到n个m维的测量向量；
S1.2、将n个测量向量按照采样时刻进行升序排列，得到数据矩阵X＝[x(1),x(2),...,x(r),...,x(n)]，x(r)表示第r时刻得到的测量向量，r＝1,2,...,n；
S1.3、对数据矩阵X中测量向量进行标准化，得到标准化的样本矩阵
Figure FDA0002472465100000011
Figure FDA0002472465100000012
表示第r时刻标准化后的测量向量。
3.根据权利要求2所述的基于ICA与ACVA的故障检测方法，其特征在于，所述步骤S2包括以下分步骤：
S2.1、通过PAC白化算法对标准化样本矩阵
Figure FDA0002472465100000021
进行白化处理，得到白化变换矩阵Q；
S2.2、根据白化变换矩阵Q，获取标准化样本矩阵
Figure FDA0002472465100000022
的白化矩阵Z，所述白化矩阵Z具体为：
Figure FDA0002472465100000023
S2.3、根据快速ICA算法对ICA模型Sz＝BZ进行求解，得到白化矩阵Z的第一解混矩阵B，Sz表示将Z经过解混后的独立向量矩阵；
S2.4、根据第一解混矩阵B，获取标准化样本矩阵
Figure FDA0002472465100000024
的第二解混矩阵W为：
W＝BTQ
其中，T表示转置。
4.根据权利要求3所述的基于ICA与ACVA的故障检测方法，其特征在于，所述步骤S3包括以下分步骤：
S3.1、设定分量数d，将第二解混矩阵W的行向量按照范数从大到小的顺序依次选取d个行向量组成矩阵，得到组合矩阵Wd；
S3.2、使用组合矩阵Wd获取标准化样本矩阵
Figure FDA0002472465100000025
中第r时刻的测量向量
Figure FDA0002472465100000026
的独立向量s(r)，得到独立向量矩阵S＝[s(1),s(2),...,s(r),...,s(n)]，所述独立向量s(r)为：
Figure FDA0002472465100000027
S3.3、设定时滞长度q，根据独立向量矩阵S获取第t时刻的过去样本向量sp(t)和未来样本向量sf(t)为：
Figure FDA0002472465100000031
Figure FDA0002472465100000032
其中，q+1≤t≤n-q+1，Rdp表示d×q行1列的实数矩；
S3.4、根据过去样本向量sp(t)和未来样本向量sf(t)，获取过去Hankel矩阵Sp和未来Hankel矩阵Sf为：
Sp＝[sp(q+1),sp(q+2),...,sp(t),...,sp(q+N)]∈Rdq×N
Sf＝[sf(q+1),sf(q+2),...,sf(t),...,sf(q+N)]∈Rdq×N
其中，N＝n-2q+1，Rdq×N表示dq行N列的实数矩阵；
S3.5、分别对过去Hankel矩阵Sp和未来Hankel矩阵Sf进行标准化，分别得到标准化过去Hankel矩阵
Figure FDA0002472465100000033
和标准化未来Hankel矩阵
Figure FDA0002472465100000034
为：
Figure FDA0002472465100000035
Figure FDA0002472465100000036
S3.6、根据标准化过去Hankel矩阵
Figure FDA0002472465100000037
和标准化未来Hankel矩阵
Figure FDA0002472465100000038
获取过去Hankel矩阵协方差∑pp、未来Hankel矩阵协方差∑ff和互协方差∑fp为：
Figure FDA0002472465100000039
Figure FDA00024724651000000310
Figure FDA00024724651000000311
其中，T表示转置；
S3.7、根据过去Hankel矩阵协方差∑pp、未来Hankel矩阵协方差∑ff和互协方差∑fp获取ACVA模型H为：
Figure FDA0002472465100000041
其中，
Figure FDA0002472465100000042
表示对∑ff求取逆矩阵后再进行矩阵开方计算，
Figure FDA0002472465100000043
表示对∑pp求取逆矩阵后再进行矩阵开方计算；
S3.8、对ACVA模型H取奇异值分解，得到H＝UΛVT，U＝[u1,u2,...,udq]，V＝[v1,v2,...,vdq]，Λ＝diag(σ1,σ2,...,σdq)，U表示左奇异矩阵，u1,u2,...,udq均表示左奇异矩阵U的列向量，Λ表示奇异值矩阵，diag表示对角矩阵，σ1,σ2,...,σdq均表示奇异值，V表示右奇异矩阵，v1,v2,...,vdq均表示右奇异矩阵的列向量；
S3.9、设定阈值threshold，按奇异值从大到小的顺序依次选择奇异值，并判断选择奇异值的和是否大于或等于threshold，若是，则停止选择，记录选择的奇异值个数g，并进入步骤S3.10，否则，重复步骤S3.9；
S3.10、选取右奇异矩阵V中前g个列向量，得到规范空间的右奇异矩阵Vd，并将右奇异矩阵V中剩余的列向量组成规范残差空间的右奇异矩阵Ve。
5.根据权利要求4所述的基于ICA与ACVA的故障检测方法，其特征在于，所述步骤S4中获取规范变量和规范残差变量的具体方法为：通过规范空间的右奇异矩阵Vd和规范残差空间的右奇异矩阵Ve获取第t时刻的规范变量cd(t)和规范残差变量ce(t)为：
Figure FDA0002472465100000044
Figure FDA0002472465100000045
。
6.根据权利要求5所述的基于ICA与ACVA的故障检测方法，其特征在于，所述步骤S5中分别利用规范变量和规范残差变量获取第一统计量的具体方法为：根据规范变量cd(t)和规范残差变量ce(t)，获取规范变量cd(t)和规范残差变量ce(t)的统计量为：
Figure FDA0002472465100000051
Figure FDA0002472465100000052
其中，
Figure FDA0002472465100000053
表示与规范变量cd(t)对应的统计量，
Figure FDA0002472465100000054
表示与规范残差变量ce(t)对应的统计量，将
Figure FDA0002472465100000055
和
Figure FDA0002472465100000056
作为第一统计量。
7.根据权利要求6所述的基于ICA与ACVA的故障检测方法，其特征在于，所述步骤S5中根据第一统计量预估控制限的具体方法为：
B1、使用核密度估计法KDE分别估计统计量
Figure FDA0002472465100000057
的概率密度
Figure FDA0002472465100000058
和统计量
Figure FDA0002472465100000059
的概率密度
Figure FDA00024724651000000510
B2、根据统计量
Figure FDA00024724651000000511
和统计量
Figure FDA00024724651000000512
获取数据跨度diff为
diff＝max(data)-min(data)
其中，data取值为
Figure FDA00024724651000000513
或
Figure FDA00024724651000000514
B3、根据数据跨度diff选取控制限求解区间[start,end]，start＝min(data)-diff/2，end＝max(data)+diff/2；
B4、设定小区间分割个数为c和控制限的置信度α，将控制限求解区间[start,end]分割为c个小区间；
B5、根据概率密度
Figure FDA00024724651000000515
和概率密度
Figure FDA00024724651000000516
在控制限求解区间[start,end]上从左至右叠加小区间的面积值，并判断叠加的面积值是否大于或等于c个小区间总面积的α，若是，则得到对应的叠加区间数k，并进入步骤B6，否则重复步骤B5；
B6、通过叠加区间数k获取控制限UCLα为start+k△，△表示每个小区间的区间长度，△＝(end-start)/c。
8.根据权利要求7所述的基于ICA与ACVA的故障检测方法，其特征在于，所述步骤B5中k值的求取公式为：
Figure FDA0002472465100000061
其中，j取值为d或e，
Figure FDA0002472465100000062
表示统计量
Figure FDA0002472465100000065
在(start+i△)处的概率密度，将其作为在第i个小区间上的高度；
Figure FDA0002472465100000063
表示统计量
Figure FDA0002472465100000064
在(start+i△)处概率密度，将其作为在第i个小区间上的高度，△表示小区间的宽度。
CN202010352906.5A 2020-04-29 2020-04-29 一种基于ica与acva的故障检测方法 Expired - Fee Related CN111474920B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202010352906.5A CN111474920B (zh) 2020-04-29 2020-04-29 一种基于ica与acva的故障检测方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202010352906.5A CN111474920B (zh) 2020-04-29 2020-04-29 一种基于ica与acva的故障检测方法

Publications (2)

Publication Number Publication Date
CN111474920A true CN111474920A (zh) 2020-07-31
CN111474920B CN111474920B (zh) 2021-04-06

Family

ID=71761997

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202010352906.5A Expired - Fee Related CN111474920B (zh) 2020-04-29 2020-04-29 一种基于ica与acva的故障检测方法

Country Status (1)

Country Link
CN (1) CN111474920B (zh)

Cited By (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN111898313A (zh) * 2020-06-30 2020-11-06 电子科技大学 一种基于ica与svm集成学习的故障检测方法
CN112229435A (zh) * 2020-10-28 2021-01-15 震兑工业智能科技有限公司 基于规范变量分析的船舶主机偶发性故障检测方法及系统

Citations (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN109272154A (zh) * 2018-09-11 2019-01-25 浙江大学 一种基于典型变量分析与隐马尔可夫的智能电厂送风机故障退化状态预测方法
CN109538311A (zh) * 2018-09-21 2019-03-29 浙江大学 面向高端发电装备中汽轮机的控制性能实时监测方法

Patent Citations (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN109272154A (zh) * 2018-09-11 2019-01-25 浙江大学 一种基于典型变量分析与隐马尔可夫的智能电厂送风机故障退化状态预测方法
CN109538311A (zh) * 2018-09-21 2019-03-29 浙江大学 面向高端发电装备中汽轮机的控制性能实时监测方法

Non-Patent Citations (4)

* Cited by examiner, † Cited by third party
Title
LI X等: "Adaptive Canonical Variate Analysis for Machine Fault Detection and Identification", 《2018 PROGNOSTICS AND SYSTEM HEALTH MANAGEMENT CONFERENCE (PHM-CHONGQING)》 *
ZHANG JIAN等: "A Novel Deep DPCA-SVM Method for Fault Detection in Industrial Processes", 《2019 IEEE 58TH CONFERENCE ON DECISIONANDCONTROL(CDC)》 *
李志农 等: "基于CVA-ICA的机械故障源动态盲分离方法", 《机械工程学报》 *
赵小强 等: "基于滑动窗的CVA故障诊断算法", 《兰州理工大学学报》 *

Cited By (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN111898313A (zh) * 2020-06-30 2020-11-06 电子科技大学 一种基于ica与svm集成学习的故障检测方法
CN111898313B (zh) * 2020-06-30 2022-05-20 电子科技大学 一种基于ica与svm集成学习的故障检测方法
CN112229435A (zh) * 2020-10-28 2021-01-15 震兑工业智能科技有限公司 基于规范变量分析的船舶主机偶发性故障检测方法及系统

Also Published As

Publication number Publication date
CN111474920B (zh) 2021-04-06

Similar Documents

Publication Publication Date Title
CN108062565B (zh) 基于化工te过程的双主元-动态核主元分析故障诊断方法
CN108803520B (zh) 一种基于变量非线性自相关性剔除的动态过程监测方法
CN109675935B (zh) 一种变控制限的ipca轧制过程在线故障诊断方法
CN105259895B (zh) 一种工业过程微小故障的检测和分离方法及其监测系统
CN108762228A (zh) 一种基于分布式pca的多工况故障监测方法
CN111474920B (zh) 一种基于ica与acva的故障检测方法
CN111949012B (zh) 一种基于双权重多邻域保持嵌入算法的间歇过程故障检测方法
CN104699077B (zh) 一种基于嵌套迭代费舍尔判别分析的故障变量隔离方法
CN108664002A (zh) 一种面向质量的非线性动态过程监控方法
CN111368428B (zh) 一种基于监控二阶统计量的传感器精度下降故障检测方法
US20220317672A1 (en) A Visualization Method for Process Monitoring Based on Bi-kernel T-distributed Stochastic Neighbor Embedding
CN110738259B (zh) 一种基于Deep DPCA-SVM的故障检测方法
CN111324110A (zh) 一种基于多块收缩自动编码器的发酵过程故障监测方法
CN103926919A (zh) 基于小波变换和Lasso函数的工业过程故障检测方法
JP3718765B2 (ja) プラント診断装置
Jia et al. Dynamic higher-order cumulants analysis for state monitoring based on a novel lag selection
CN112947649B (zh) 一种基于互信息矩阵投影的多变量过程监控方法
CN109325065B (zh) 基于动态隐变量模型的多采样率软测量方法
CN112149054B (zh) 基于时序扩展的正交邻域保持嵌入模型的构建与应用
CN113341926A (zh) 一种基于稀疏加权邻域保持嵌入的多阶段间歇过程故障检测方法
CN103995985B (zh) 基于Daubechies小波变换和弹性网的故障检测方法
CN111188761B (zh) 一种基于Fourier-CVA模型面向机泵设备的监测方法
CN116304823A (zh) 一种垃圾焚烧过程在线诊断方法
CN115047853A (zh) 基于递推规范变量残差和核主元分析的微小故障检测方法
CN114200914A (zh) 一种基于mw-occa的质量相关早期故障检测方法

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant
CF01 Termination of patent right due to non-payment of annual fee
CF01 Termination of patent right due to non-payment of annual fee

Granted publication date: 20210406