[go: up one dir, main page]

CN109932908B - A multi-directional principal component analysis process monitoring method based on alarm reliability fusion - Google Patents

A multi-directional principal component analysis process monitoring method based on alarm reliability fusion Download PDF

Info

Publication number
CN109932908B
CN109932908B CN201910211439.1A CN201910211439A CN109932908B CN 109932908 B CN109932908 B CN 109932908B CN 201910211439 A CN201910211439 A CN 201910211439A CN 109932908 B CN109932908 B CN 109932908B
Authority
CN
China
Prior art keywords
spe
alarm
alarm reliability
reliability
matrix
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.)
Active
Application number
CN201910211439.1A
Other languages
Chinese (zh)
Other versions
CN109932908A (en
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.)
Hangzhou Dianzi University
Original Assignee
Hangzhou Dianzi University
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 Hangzhou Dianzi University filed Critical Hangzhou Dianzi University
Priority to CN201910211439.1A priority Critical patent/CN109932908B/en
Publication of CN109932908A publication Critical patent/CN109932908A/en
Application granted granted Critical
Publication of CN109932908B publication Critical patent/CN109932908B/en
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Landscapes

  • Management, Administration, Business Operations System, And Electronic Commerce (AREA)
  • Testing And Monitoring For Control Systems (AREA)

Abstract

本发明涉及一种基于报警信度融合的多向主元分析过程监测方法，包括模型数据采集、三维数据展开、二维矩阵标准化、PCA分解、主元个数选取、计算SPE和T2统计控制限，建立模糊隶属度函数，计算报警信度以及全局报警信度的权重，并进行报警决策。本发明克服了多向主元分析法模型统计控制限报警的局限性，利用历史报警信度和当前报警信度融合的方法来进行报警决策，大大降低多向主元分析在间歇过程监测过程应用中的误报率，为间歇过程的建模以及监测提供了新的途径。

Figure 201910211439

The invention relates to a multi-directional principal element analysis process monitoring method based on alarm reliability fusion, including model data acquisition, three-dimensional data expansion, two-dimensional matrix standardization, PCA decomposition, selection of the number of principal elements, calculation of SPE and T2 statistical control limits , establish a fuzzy membership function, calculate the alarm reliability and the weight of the global alarm reliability, and make an alarm decision. The invention overcomes the limitation of the multi-directional principal component analysis model statistical control limit alarm, uses the method of fusion of historical alarm reliability and current alarm reliability to make alarm decision, and greatly reduces the application of multi-directional principal component analysis in the intermittent process monitoring process The false alarm rate in this paper provides a new way for the modeling and monitoring of intermittent processes.

Figure 201910211439

Description

Multi-directional principal component analysis process monitoring method based on alarm reliability fusion
Technical Field
The invention relates to the field of intermittent process multivariable monitoring and fault monitoring, and relates to a multidirectional principal component analysis process monitoring method based on alarm reliability fusion.
Background
Although the multi-way principal component analysis (MPCA) is applied to the condition monitoring and fault diagnosis of the intermittent process, the method makes several assumptions in practical application: 1) the data of the intermittent process is linear, 2) the intermittent process is always in a stable running state, the process data conforms to uniform distribution, 3) the whole intermittent production process is single, and the like. The practical industrial production cannot meet the above assumptions, and the MPCA uses a single constant statistical control limit to easily cause alarm errors, for example, the 95% probability statistical control limit is too strict to cause too many statistical false alarm, and the 99.999% probability statistical control limit is too loose to cause abnormal process monitoring, so the probability statistical alarm control limit is greatly hindered in factory application.
Disclosure of Invention
The invention aims to overcome the defects of the prior art and provide a batch process monitoring alarm method based on information fusion, which can greatly improve the alarm accuracy of multi-directional principal component analysis in the application of the intermittent process monitoring process, and adopts the following technical scheme in order to achieve the aim:
the invention comprises the following steps:
1) model data collection
If an intermittent operation is provided with J measurement variables and K sampling points, each measurement batch can obtain a J multiplied by K matrix, and after the measurement steps of the I batch are repeated, the obtained data can be expressed as a three-dimensional matrixX(I.times.J.times.K), wherein the measured variables are temperature, speed, pressure, stroke, etc. state parameters that can be measured during the batch run.
2) Three-dimensional data expansion
Combining three-dimensional matricesX(I multiplied by J multiplied by K) is expanded according to the collection batch direction, namely variables on each sampling point in an operation batch are arranged according to the time sequence to obtain a two-dimensional matrix
Figure BDA0002000657950000012
Explicit matrix
Figure BDA0002000657950000013
Column KJ, line I.
3) Two dimensional matrix normalization
Setting a two-dimensional matrix
Figure BDA0002000657950000014
The variable at any point in the variable is
Figure BDA0002000657950000015
The variable is normalized by subtracting the mean and dividing the variance, and the calculation formula of the normalization is as follows:
Figure BDA0002000657950000011
wherein:
Figure BDA0002000657950000021
is that
Figure BDA0002000657950000027
Mean of any column of the matrix, sjkIs that
Figure BDA0002000657950000028
Variance of any column of the matrix
Figure BDA0002000657950000022
Figure BDA0002000657950000023
4) Multi-directional principal component analysis modeling
Standardizing the two-dimensional matrix of the last step
Figure BDA0002000657950000029
And (3) performing PCA decomposition to complete the modeling of the multi-directional principal component analysis method, wherein the PCA decomposition process is as follows:
Figure BDA0002000657950000024
S＝trace(TTT/(I-1))； (3)
wherein: t is tiIs orthogonal principal component directionAmount, piThe load vector is orthogonal normalization, S is the trace of the covariance matrix of the principal element and represents the interpretation degree of each principal element to the process;
in the formula (2), X is decomposed to obtain a score matrix T (I multiplied by JK) and a load matrix P (JK multiplied by JK).
5) Selecting the number R of principal elements
Restated equation (2) to the form:
Figure BDA0002000657950000025
wherein: t isr(I×R)、Pr(JK multiplied by R) are respectively a score matrix and a load matrix after R principal elements are reserved, and E is a residual error matrix;
through the above transformation, the multi-way principal component analysis model decomposes the original data space into principal component space and residual space, and the principal component space variables are highly correlated and generally sufficient to describe the variability of the data.
The number of principal elements R can be generally set according to the experience of the user or by adopting the Broken-Stick criterion, the content of the Broken-Stick is to keep the principal element when the interpretation degree S (R) of the R-th principal element accounts for a percentage of the total contribution sum (S) of all principal elements which is more than G (R), and otherwise, the result is terminated, wherein the calculation formula of G (R) is as follows:
Figure BDA0002000657950000026
wherein: s (r) is the degree of interpretation of the r-th principal element, and sum (S) is the sum of the contributions of all principal elements.
6) Two indices of the model data, T2 and SPE, and corresponding 95% and 99.999% control limits were calculated
a) 95% control limit of T2:
Figure BDA0002000657950000031
FR,n-R,αrepresenting a F distribution with a degree of freedom R, n-R, alpha being 5%;
99.999% control of T2Limiting:
Figure BDA0002000657950000032
FR,n-R,αrepresenting a F distribution with a degree of freedom R, n-R, alpha is 0.001%;
b) 95% control limit for SPE:
Figure BDA0002000657950000033
Figure BDA0002000657950000034
Cαis a critical value at the normal distribution calibration level alpha of 5%, Cα＝norminv(0.95,0,1)；
99.999% control limit for SPE:
Figure BDA0002000657950000035
wherein theta isiAnd h0Is calculated as above, where α is 0.001%, Cα＝norminv(0.99999,0,1)。
7) Index SPE and T2 for calculating two spaces of sample data
An index of the sample residual space, SPE, calculates equation (10), where I is the identity matrix,
SPE＝EET＝X(I-PrPr T)XT (10)
the index T2 of the sample pivot space calculates equation (11),
T2＝Tr TSr -1Tr (11)
8) obtaining historical alarm confidence based on fuzzy threshold
Respectively taking the 99.999 percent and 95 percent statistical control limits of T2 and SPE as fuzzy upper and lower limits to construct a fuzzy membership function
Figure BDA0002000657950000036
Converting the alarm of which the statistic exceeds the control threshold into the alarm reliability,
Figure BDA0002000657950000037
representing a 99.999% probability statistical control limit for SPE or T2,
Figure BDA0002000657950000038
a 95% probability statistical control limit representing SPE or T2;
Figure BDA0002000657950000039
respectively substituting the upper and lower limits of the two statistics in the SPE, the T2 and the step 6) calculated in the step 7) into a function
Figure BDA0002000657950000041
Obtaining two alarm reliability sequences Bspe{bspe(1),bspe(2)……bspe(t)}，BT2{bT2(1),bT2(2)……bT2(t), wherein t is the sampling time;
9) computing combining weights based on alarm confidence
Firstly, an alarm reliability sequence B is calculatedspe(t){bspe(1),bspe(2)……bspe(K) The combining weights of
(a) When t is 1, since there is no relevant alarm information before t is 1, the global alarm reliability when t is 1 is the alarm reliability b at the first time0:1(NA)＝bspe(1)；
(b) When t is 2, the global alarm reliability is weighted and fused by the global alarm reliability at the time when t is 1 and the alarm reliability at the time when t is 2, and the global alarm reliability b at the time when t is 2 is obtained0:2(NA), calculated as shown in formula:
b0:2(NA)＝α2b0:1(NA)+β2bspe(2) (13)
wherein alpha is20.5 denotes the pair b0:1Linear weighted value of (NA), beta2The weighted value of the alarm reliability at the current moment is expressed as 0.5;
(c) when t is more than or equal to 3, the global alarm reliability carries out weighted fusion by using the global alarm reliability of the first two moments and the alarm reliability of the current moment,
b0:t(NA)＝αt b0:t-1(NA)+βtbspe(t) (14)
wherein the weight coefficient alphatWeight factor, beta, representing the global confidence of the last momenttIs a weighted value of the alarm reliability at the current moment, and the weighted value satisfy the relation alphat＝1-βt. The formula is as follows (15) - (16)
Figure BDA0002000657950000042
Figure BDA0002000657950000043
Wherein Supt (b)i) The support degree of the alarm reliability is represented, and the calculation method is as follows:
(i) when t is more than or equal to 3, the alarm reliability b at the time t can be obtainedspeAnd global alarm evidence b at t-1 and t-20:t-1And b0:t-2Calculate b0:t-1And b0:t-2The distance between every two adopts an Euclidean distance formula:
Figure BDA0002000657950000044
then b0:t-1And b0:t-2The similarity between the alarm reliability degrees is determined,
Sim(b0:t-1,b0:t-2)＝1-dst(b0:t-1,b0:t-2) (18)
obtained in the same way, Sim (b)0:t-2,bspe)，Sim(b0:t-1,bspe)
(ii) Calculation of bspeAnd b0:t-1And b0:t-2The support degree of each alarm reliability relative to the other two alarm reliabilities, so that each evidenceThe degrees of support (19) to (21) are shown in the formulas:
Supt(bt)＝Sim(b0:t-1,bspe)+Sim(b0:t-2,bspe) (19)
Supt(b0:t-1)＝Sim(b0:t-1,bspe)+Sim(b0:t-2,b0:t-1) (20)
Supt(b0:t-2)＝Sim(b0:t-1,b0:t-2)+Sim(b0:t-2,bspe) (21)
in the same way, the sequence B of the available alarm credibilityT2Combining weights of
10) Alarm decision based on global alarm reliability
When b is0:tAnd (NA) is more than or equal to 0.5, monitoring and alarming.
Compared with the prior art, the invention has the beneficial effects that: the method integrates the monitoring information of the historical moment and the current moment to make alarm decision, improves the accuracy of the process monitoring fault diagnosis result, and provides a new practical method for process monitoring without process prior knowledge.
Drawings
FIG. 1 is a fuzzy membership function according to the present invention
Figure BDA0002000657950000051
FIG. 2 is a schematic diagram showing the development of three-dimensional data of a batch process based on multi-way principal component analysis according to the present invention.
FIG. 3 shows MPCA monitoring of a batch of an injection molding process in an embodiment of the present invention.
Detailed Description
The invention is further described with reference to the following drawings and specific embodiments:
injection molding is a typical batch process, one batch of which mainly comprises four stages of injection, pressure holding, plasticizing and cooling, specifically, in an injection section, a high pressure of a hydraulic cylinder pushes a screw forward to push molten plastic in a machine barrel to a mold cavity, when the mold cavity is completely or nearly filled, the process is switched to a pressure holding stage, in which the high pressure continues to fill a small amount of material into the mold cavity to supplement material shrinkage due to cooling and solidification; when the gate cools and the material in the mold cavity is no longer affected by the injection nozzle, the hold pressure stage ends. The screw is rotated and retracted to push a sufficient amount of molten plastic to the front of the screw. The screw retreats while the volume calculation is started. After the head melt reaches a certain injection quantity, the screw stops moving back and rotating, and the process state of the time is called a plasticizing stage. And (4) when the pressure maintaining section is finished and the plasticizing process is carried out, simultaneously carrying out the cooling section until the material in the die reaches the hardness capable of being ejected, and finishing the cooling section.
The injection molding process is continuously repeated to obtain stable and consistent products, and by taking the injection molding process as an example, the multivariate statistical monitoring method based on the alarm reliability comprises the following steps of:
(1) model data collection
If an intermittent operation is provided with J measurement variables and K sampling points, each measurement batch can obtain a J multiplied by K matrix, and after the measurement steps of the I batch are repeated, the obtained data can be expressed as a three-dimensional matrixX(I.times.JXK). In order to ensure that the detection data covers a working range of a long enough time, the value of a data batch I used for modeling in general industry is more than 100, and the measurement variables are state parameters which can be measured in the batch operation process such as temperature, speed, pressure, stroke and the like; based on the process time, the speed of process change and whether the model load is in a reasonable range, the number of sampling points K is generally less than 1000.
In this embodiment, the variables available for measuring the operating process of the variable laboratory injection molding machine are 8: the opening degree of a pressure valve, the opening degree of a flow valve, an injection stroke, an injection speed, an injection pressure and the temperature of a machine barrel (3 sections), 100 is taken for an operation batch I, and a sampling point K reserved for each batch is 488.
(2) Three-dimensional data expansion
See alsoFIG. 2, shows a three-dimensional matrixXExpanding according to the collection batch direction, namely arranging the variables on each sampling point in an operation batch according to the time sequence to obtain a two-dimensional matrix
Figure BDA0002000657950000065
Matrix array
Figure BDA0002000657950000066
Is line I, column KJ;
(3) two dimensional matrix normalization
Setting a two-dimensional matrix
Figure BDA0002000657950000067
The variable at any point in the variable is
Figure BDA0002000657950000068
The variable is normalized by subtracting the mean and dividing the variance, and the calculation formula of the normalization is as follows:
Figure BDA0002000657950000061
wherein:
Figure BDA0002000657950000062
is that
Figure BDA0002000657950000069
Mean of any column of the matrix, sjkIs that
Figure BDA00020006579500000610
The variance of any column of the matrix is,
Figure BDA0002000657950000063
Figure BDA0002000657950000064
the standardization processing of the step is equivalent to extracting the average running track of one operation in the intermittent process, and highlights normal random fluctuation among different operation batches in the intermittent process.
(4) MPCA modeling
The MPCA modeling is a method of firstly expanding a three-dimensional matrix into a large two-dimensional matrix and then performing conventional PCA decomposition, the step of performing PCA decomposition on the two-dimensional matrix X (I multiplied by JK) subjected to the normalization processing in the previous step comprises the following decomposition process:
Figure BDA0002000657950000071
S＝trace(TTT/(I-1))； (3)
wherein: t is tiBeing orthogonal principal component vectors, piFor an orthogonally normalized load vector, S is the trace of the covariance matrix of the principal elements, representing the magnitude of the interpretation of the process by each principal element.
The formula (2) X is decomposed to obtain a score matrix T (I multiplied by JK) and a load matrix P (JK multiplied by JK).
(5) Selecting the number of principal elements
Generally, the first few principal elements generally contain most of the variance information of the gap process, and the other principal elements may mainly contain noise information, so equation (2) can be restated as follows:
Figure BDA0002000657950000072
wherein: t isr(I×R)、Pr(JK multiplied by R) are respectively a score matrix and a load matrix after R principal elements are reserved, and E is a residual error matrix;
with the above transformation, the MPCA model decomposes the original data space into a principal component space and a residual space, the principal component space variables being highly correlated, generally sufficient to describe the variability of the data.
The number of principal elements R can be generally set according to the experience of the user or by adopting the Broken-Stick criterion, the content of the Broken-Stick is to keep the principal element when the interpretation degree S (R) of the R-th principal element accounts for a percentage of the total contribution sum (S) of all principal elements which is more than G (R), and otherwise, the result is terminated, wherein the calculation formula of G (R) is as follows:
Figure BDA0002000657950000073
wherein: s (R) is the interpretation degree of the R-th principal element, and sum (S) is the sum of all the principal elements, in this embodiment, the number of the principal elements R is selected to be 3, and the interpretation degree for the process is 86.73.
(6) Calculating model statistical control limits
a) 95% control limit of T2:
Figure BDA0002000657950000081
FR,n-R,αrepresenting a F distribution with a degree of freedom R, n-R, alpha is 5%.
99.999% control Limit for T2:
Figure BDA0002000657950000082
FR,n-R,αrepresenting a F distribution with a degree of freedom R, n-R, alpha is 0.001%.
After R is 3 and n is 100, the above formula is substituted, TsqCL95 is 8.3447 and TsqCL999 is 30.4307;
b) 95% control limit for SPE:
Figure BDA0002000657950000083
Figure BDA0002000657950000084
Cαis a critical value at the normal distribution calibration level alpha of 5%, Cα＝norminv(0.95,0,1)＝1.6449。
99.999% control limit for SPE:
Figure BDA0002000657950000085
wherein theta isiAnd h0Is calculated as above, where α is 0.001%, Cα＝norminv(0.99999,0,1)＝4.2649。
(7) Computing statistics of two spaces of sample data
SPE＝EET＝X(I-PrPr T)XT (10)
T2＝Tr TSr -1Tr (11)。
(8) Obtaining historical alarm confidence based on fuzzy threshold
Obtaining two alarm reliability sequences B of a residual error space and a principal component space by using formulas (10) to (11)spe{bspe(1),bspe(2)……bspe(t)}，BT2{bT2(1),bT2(2)……bT2(t), t is the sampling time, and the fuzzy membership function thereof
Figure BDA0002000657950000086
See figure 1.
(9) Separately calculating combining weights based on alarm confidence
Firstly, an alarm reliability sequence B is calculatedspe(t){bspe(1),bspe(2)……bspe(K) The combining weight of the { right } is greater than the first;
(a) when t is 1, since there is no relevant alarm information before t is 1, the global alarm reliability when t is 1 is the alarm reliability b at the first time0:1(NA)＝bspe(1)
(b) When t is 2, the global alarm reliability is weighted and fused by the global alarm reliability at the time when t is 1 and the alarm reliability at the time when t is 2, and the global alarm reliability b at the time when t is 2 is obtained0:2(NA), calculated as shown in formula:
b0:2(NA)＝α2b0:1(NA)+β2bspe(2) (12)
wherein alpha is20.5 denotes the pair b0:1Linearity of (NA)Weighted value, beta2The weighted value of the alarm reliability at the current moment is expressed as 0.5;
(c) when t is more than or equal to 3, the global alarm reliability carries out weighted fusion by using the global alarm reliability of the first two moments and the alarm reliability of the current moment,
b0:t(NA)＝αt b0:t-1(NA)+βtbspe(t) (13)
wherein the weight coefficient alphatWeight factor, beta, representing the global confidence of the last momenttThe weighted value of the alarm reliability at the current moment is calculated by the following method:
(i) when t is more than or equal to 3, the alarm reliability b at the time t can be obtainedspeAnd global alarm evidence b at t-1 and t-20:t-1And b0:t-2Calculate b0:t-1And b0:t-2The distance between every two adopts an Euclidean distance formula:
Figure BDA0002000657950000091
then b0:t-1And b0:t-2The similarity between the alarm reliability degrees is determined,
Sim(b0:t-1,b0:t-2)＝1-dst(b0:t-1,b0:t-2) (15)
obtained in the same way, Sim (b)0:t-2,bspe)，Sim(b0:t-1,bspe)
(ii) Calculation of bspeAnd b0:t-1And b0:t-2The support of each alarm confidence level relative to the other two alarm confidence levels, so the support of each evidence is shown by the equations (16) - (18):
Supt(bspe)＝Sim(b0:t-1,bspe)+Sim(b0:t-2,bspe) (16)
Supt(b0:t-1)＝Sim(b0:t-1,bspe)+Sim(b0:t-2,b0:t-1) (17)
Supt(b0:t-2)＝Sim(b0:t-1,b0:t-2)+Sim(b0:t-2,bspe) (18)
(iii) on the basis of the steps, the support degree of each alarm reliability relative to other two alarm reliability is obtained, and finally, a dynamically updated weight factor beta is obtainedtThe following formula (19);
Figure BDA0002000657950000092
αt＝1-βt
similarly, a principal component spatial belief sequence B can be computedT2{bT2(1),bT2(2)……bT2(t) } alarm confidence combining weights.
(10) Alarm decision based on global alarm reliability
When b is0:tAnd (NA) is more than or equal to 0.5, monitoring and alarming.
The results of this example are as follows, see FIG. 3:
Figure BDA0002000657950000101
according to the invention, the alarm reliability of the historical moment and the current moment is integrated to carry out alarm decision, so that the accuracy and effectiveness of the process monitoring fault diagnosis result are improved.
It should be understood that the present invention is not limited to the injection molding process of the above-described embodiments, and that equivalent modifications or substitutions can be made by those skilled in the art without departing from the spirit of the present invention, and are intended to be included within the scope of the appended claims.

Claims (5)

1.一种基于报警信度融合的多向主元分析过程监测方法，其特征在于，包括以下步骤：1. a multi-directional principal component analysis process monitoring method based on the fusion of alarm reliability, is characterized in that, comprises the following steps: 1)模型数据采集1) Model data collection 设一个间歇操作具有J个测量变量和K个采样点，则每一个测量批次可得到一个J×K的矩阵，重复I批次的测量步骤后，得到的数据表述为一个三维矩阵X(I×J×K)，所述的测量变量包括温度、速度、压力或行程；Assuming that an intermittent operation has J measurement variables and K sampling points, a J×K matrix can be obtained for each measurement batch. After repeating the measurement steps of I batches, the obtained data is expressed as a three-dimensional matrix X (I ×J×K), the measured variables include temperature, speed, pressure or stroke; 2)三维数据展开2) 3D data expansion 将三维矩阵
Figure FDA0003387896500000011
按照采集批次方向展开，即将一个操作批次内的各采样点上的变量按照时间顺序排开得到二维矩阵
Figure FDA0003387896500000012
3D matrix
Figure FDA0003387896500000011
Expand according to the direction of the collection batch, that is, arrange the variables on each sampling point in an operation batch in chronological order to obtain a two-dimensional matrix
Figure FDA0003387896500000012
 3)二维矩阵标准化3) Two-dimensional matrix normalization 设二维矩阵
Figure FDA0003387896500000013
内任意一点的变量为
Figure FDA0003387896500000014
对该变量进行减均值、除方差的标准化处理；Let a two-dimensional matrix
Figure FDA0003387896500000013
The variable at any point in the
Figure FDA0003387896500000014
Standardize the variable by subtracting the mean and dividing the variance; 4)多向主元分析法建模4) Multi-directional principal component analysis method modeling 对上一步经标准化后的二维矩阵
Figure FDA0003387896500000015
执行PCA分解，完成多向主元分析法的建模，其中PCA分解过程如下：For the two-dimensional matrix normalized in the previous step
Figure FDA0003387896500000015
Execute PCA decomposition to complete the modeling of multi-directional principal component analysis. The PCA decomposition process is as follows:
Figure FDA0003387896500000016
Figure FDA0003387896500000016
 S＝trace(TTT/(I-1))； (3)S=trace(T T T/(I-1)); (3) 其中：ti为正交的主元向量，pi为正交归一化的负载向量，S是主元的协方差矩阵的迹，代表各个主元对于过程的解释度大小；Among them: t i is the orthogonal pivot vector, pi is the orthonormalized load vector, S is the trace of the covariance matrix of the pivot, representing the degree of explanation of each pivot for the process; 公式(2)中X分解得到得分矩阵T(I×JK)及负载矩阵P(JK×JK)；In formula (2), X is decomposed to obtain the score matrix T(I×JK) and the load matrix P(JK×JK); 5)选取主元个数R5) Select the number of pivots R 将公式(2)重新表述成如下形式：Formula (2) can be reformulated as follows:
Figure FDA0003387896500000017
Figure FDA0003387896500000017
 其中：Tr(I×R)、Pr(JK×R)分别为保留R个主元后的得分矩阵和负载矩阵，E为残差矩阵；Among them: T r (I×R), P r (JK×R) are the score matrix and load matrix after retaining R pivot elements respectively, and E is the residual matrix; 6)计算模型数据的两个指标T2和SPE以及对应的95％和99.999％控制限6) Calculate the two indicators T2 and SPE of the model data and the corresponding 95% and 99.999% control limits a)T2的95％控制限：
Figure FDA0003387896500000018
a) 95% control limit for T2:
Figure FDA0003387896500000018
 FR,n-R,α表示自由度为R，n-R的F分布；F R, nR, α represents the F distribution with degrees of freedom R, nR; T2的99.999％控制限：
Figure FDA0003387896500000021
99.999% control limit for T2:
Figure FDA0003387896500000021
 FR,n-R,α表示自由度为R，n-R的F分布；F R, nR, α represents the F distribution with degrees of freedom R, nR; b)SPE的95％控制限：
Figure FDA0003387896500000022
b) 95% control limit for SPE:
Figure FDA0003387896500000022
 θi＝trace(Vi)，V＝ETE/(I-1)，i＝1，2，3，
Figure FDA0003387896500000023
Cα是正态分布校准水平α下的临界值；θ i =trace(V i ), V=E T E/(I-1), i=1, 2, 3,
Figure FDA0003387896500000023
C α is the critical value at the normal distribution calibration level α; SPE的99.999％控制限：
Figure FDA0003387896500000024
99.999% Control Limit for SPE:
Figure FDA0003387896500000024
 其中θi和h0的计算同上；The calculation of θ i and h 0 is the same as above; 7)计算样本数据两个空间的指标SPE和T27) Calculate the indicators SPE and T2 of the two spaces of the sample data 样本残差空间的指标SPE计算如下：The index SPE of the sample residual space is calculated as follows: SPE＝EET＝X(I′-PrPr T)XT (10)SPE=EE T =X(I'-P r Pr T )X T (10) 其中I′为单位矩阵；where I' is the identity matrix; 样本主元空间的指标T2计算如下：The index T2 of the sample pivot space is calculated as follows: T2＝Tr TSr -1Tr (11)T2=T r T S r -1 T r (11) 8)基于模糊阈获取历史报警信度8) Obtain historical alarm reliability based on fuzzy threshold 分别以T2和SPE的99.999％和95％统计控制限为模糊上下限，构造模糊隶属度函数
Figure FDA0003387896500000025
来将指标超过控制阈值的报警转化报警信度，
Figure FDA0003387896500000026
代表SPE或者T2的99.999％的概率统计控制限，
Figure FDA0003387896500000027
代表SPE或者T2的95％的概率统计控制限；Taking the 99.999% and 95% statistical control limits of T2 and SPE as the fuzzy upper and lower limits, respectively, construct the fuzzy membership function
Figure FDA0003387896500000025
to convert alarms whose indicators exceed the control threshold into alarm reliability,
Figure FDA0003387896500000026
represents the 99.999% probability statistical control limit of SPE or T2,
Figure FDA0003387896500000027
Represents the 95% probability statistical control limit of SPE or T2;
Figure FDA0003387896500000028
Figure FDA0003387896500000028
 将步骤7)计算的SPE、T2以及步骤6)中两个指标的上下限分别带入函数
Figure FDA0003387896500000029
得到两个报警信度序列Bspe{bspe(1),bspe(2)……bspe(t)}，BT2{bT2(1),bT2(2)……bT2(t)},t为采样的时刻；Bring the SPE, T2 calculated in step 7) and the upper and lower limits of the two indicators in step 6) into the function respectively
Figure FDA0003387896500000029
Obtain two alarm reliability sequences B spe {b spe (1),b spe (2)…b spe (t)}, B T2 {b T2 (1),b T2 (2)…b T2 (t )}, t is the sampling time; 9)基于报警信度计算组合权重9) Calculate the combined weight based on the alarm reliability 首先计算报警信度序列Bspe(t){bspe(1),bspe(2)……bspe(K)}的组合权重；First calculate the combined weight of the alarm reliability sequence B spe (t){b spe (1),b spe (2)...b spe (K)}; (a)当t＝1时，因在t＝1之前没有相关的报警信息，所以t＝1时全局的报警信度就是第一时刻的报警信度b0:1(NA)＝bspe(1)(a) When t=1, since there is no relevant alarm information before t=1, the global alarm reliability at t=1 is the alarm reliability at the first moment b 0:1 (NA)=b spe ( 1) (b)当t＝2时，全局的报警信度利用t＝1时刻的全局报警信度和t＝2时刻的报警信度进行加权融合，得到t＝2时刻的全局报警信度b0:2(NA)，其计算如式所示：(b) When t=2, the global alarm reliability uses the global alarm reliability at time t=1 and the alarm reliability at time t=2 to perform weighted fusion to obtain the global alarm reliability b at time t=2 : 2 (NA), which is calculated as: b0:2(NA)＝α2b0:1(NA)+β2bspe(2) (13)b 0:2 (NA)=α 2 b 0:1 (NA)+β 2 b spe (2) (13) 其中α2表示对b0:1(NA)的线性加权值，β2表示当前时刻报警信度的加权值；Wherein α 2 represents the linear weighted value of b 0:1 (NA), and β 2 represents the weighted value of the alarm reliability at the current moment; (c)当t≥3时,全局的报警信度利用前两时刻的全局报警信度和当前时刻的报警信度进行加权融合，(c) When t≥3, the global alarm reliability is weighted and fused using the global alarm reliability of the previous two moments and the alarm reliability of the current moment, b0:t(NA)＝αtb0:t-1(NA)+βtbspe(t) (14)b 0:t (NA)=α t b 0:t-1 (NA)+β t b spe (t) (14) 其中权重系数αt表示上一时刻全局信度的权重因子，βt是当前时刻报警信度的加权值，二者满足关系式αt＝1-βt；其公式如下(15)-(16)Among them, the weight coefficient α t represents the weight factor of the global reliability at the previous moment, and β t is the weighted value of the alarm reliability at the current moment . )
Figure FDA0003387896500000031
Figure FDA0003387896500000031
Figure FDA0003387896500000032
Figure FDA0003387896500000032
 同理，可得报警信度序列BT2的组合权重，其中Supt(bi)表示报警信度的支持度；Similarly, the combined weight of the alarm reliability sequence B T2 can be obtained, where Supt(b i ) represents the support degree of the alarm reliability; 10)基于全局报警信度的报警决策10) Alarm decision based on global alarm reliability 当b0:t(NA)≥0.5则监测报警。Monitor alarm when b 0:t (NA)≥0.5. 2.根据权利要求1所述的一种基于报警信度融合的多向主元分析过程监测方法，其特征在于，二维矩阵中的变量
Figure FDA0003387896500000033
进行减均值、除方差的标准化处理过程如下：2. a kind of multi-directional principal component analysis process monitoring method based on alarm reliability fusion according to claim 1 is characterized in that, the variable in the two-dimensional matrix
Figure FDA0003387896500000033
The standardization process of subtracting the mean and dividing the variance is as follows:
Figure FDA0003387896500000034
Figure FDA0003387896500000034
 其中：
Figure FDA0003387896500000035
是
Figure FDA0003387896500000036
矩阵任一列的均值，sjk是
Figure FDA0003387896500000037
矩阵任一列的方差，in:
Figure FDA0003387896500000035
Yes
Figure FDA0003387896500000036
the mean of any column of the matrix, s jk is
Figure FDA0003387896500000037
the variance of any column of the matrix,
Figure FDA0003387896500000041
Figure FDA0003387896500000041
Figure FDA0003387896500000042
Figure FDA0003387896500000042
 3.根据权利要求1所述的一种基于报警信度融合的多向主元分析过程监测方法，其特征在于，主元个数R根据用户的经验设定或者采用Broken-Stick准则。3 . The multi-directional principal component analysis process monitoring method based on the fusion of alarm reliability according to claim 1 , wherein the number of principal components R is set according to the user's experience or adopts the Broken-Stick criterion. 4 . 4.根据权利要求3所述的一种基于报警信度融合的多向主元分析过程监测方法，其特征在于，所述的Broken-Stick准则的内容是当第r个主元的解释度S(r)占所有主元总贡献sum(S)的百分比大于G(r)的时候保留该主元，否则终止，其中G(r)的计算公式如下：4. a kind of multi-directional principal component analysis process monitoring method based on alarm reliability fusion according to claim 3, is characterized in that, the content of described Broken-Stick criterion is when the interpretability S of the rth principal element (r) When the percentage of the total contribution sum(S) of all pivots is greater than G(r), keep the pivot, otherwise terminate, where G(r) is calculated as follows:
Figure FDA0003387896500000043
Figure FDA0003387896500000043
 其中：S(r)是第r个主元的解释度，sum(S)是所有主元的贡献和。where: S(r) is the explainability of the rth pivot, and sum(S) is the sum of the contributions of all pivots. 5.根据权利要求1所述的一种基于报警信度融合的多向主元分析过程监测方法，其特征在于，报警信度的支持度Supt(bi)计算如下：5. a kind of multi-directional principal component analysis process monitoring method based on alarm reliability fusion according to claim 1, is characterized in that, the support degree Supt (b i ) of alarm reliability is calculated as follows: (i)当t≥3之后，得到t时刻的报警信度bspe和t-1以及t-2时刻的全局报警证据b0:t-1和b0:t-2，计算b0:t-1以及b0:t-2两两之间的距离：(i) After t≥3, obtain the alarm reliability b spe and t-1 at time t and the global alarm evidence b 0:t-1 and b 0:t-2 at time t-2 , and calculate b 0:t The distance between -1 and b 0:t-2 pairwise:
Figure FDA0003387896500000044
Figure FDA0003387896500000044
 那么b0:t-1以及b0:t-2报警信度之间相似度，Then the similarity between the alarm reliability of b 0:t-1 and b 0:t-2 , Sim(b0:t-1,b0:t-2)＝1-dst(b0:t-1,b0:t-2) (18)Sim(b 0:t-1 ,b 0:t-2 )=1-dst(b 0:t-1 ,b 0:t-2 ) (18) 同理可得，Sim(b0:t-2,bspe)，Sim(b0:t-1,bspe)In the same way, Sim(b 0:t-2 ,b spe ), Sim(b 0:t-1 ,b spe ) (ii)计算bspe和b0:t-1以及b0:t-2每个报警信度相对于其他两个报警信度的支持程度，所以每个证据的支持度(19)-(21)式所示：(ii) Calculate the support degree of b spe and b 0:t-1 and b 0:t-2 of each alarm reliability relative to the other two alarm reliability, so the support degree of each evidence (19)-(21 ) is shown as: Supt(bt)＝Sim(b0:t-1,bspe)+Sim(b0:t-2,bspe) (19)Supt(b t )=Sim(b 0:t-1 ,b spe )+Sim(b 0:t-2 ,b spe ) (19) Supt(b0:t-1)＝Sim(b0:t-1,bspe)+Sim(b0:t-2,b0:t-1) (20)Supt(b 0:t-1 )=Sim(b 0:t-1 ,b spe )+Sim(b 0:t-2 ,b 0:t-1 ) (20) Supt(b0:t-2)＝Sim(b0:t-1,b0:t-2)+Sim(b0:t-2,bspe) (21)。Supt(b 0:t-2 )=Sim(b 0:t-1 ,b 0:t-2 )+Sim(b 0:t-2 ,b spe ) (21).
CN201910211439.1A 2019-03-20 2019-03-20 A multi-directional principal component analysis process monitoring method based on alarm reliability fusion Active CN109932908B (en)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201910211439.1A CN109932908B (en) 2019-03-20 2019-03-20 A multi-directional principal component analysis process monitoring method based on alarm reliability fusion

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201910211439.1A CN109932908B (en) 2019-03-20 2019-03-20 A multi-directional principal component analysis process monitoring method based on alarm reliability fusion

Publications (2)

Publication Number Publication Date
CN109932908A CN109932908A (en) 2019-06-25
CN109932908B true CN109932908B (en) 2022-03-01

Family

ID=66987806

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201910211439.1A Active CN109932908B (en) 2019-03-20 2019-03-20 A multi-directional principal component analysis process monitoring method based on alarm reliability fusion

Country Status (1)

Country Link
CN (1) CN109932908B (en)

Families Citing this family (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN116300774B (en) * 2023-05-23 2023-08-08 蓝星智云(山东)智能科技有限公司 Intermittent process visual monitoring method based on principal component analysis and nuclear density estimation

Citations (20)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
WO2003096130A1 (en) * 2001-10-23 2003-11-20 Brooks-Pri Automation, Inc. Semiconductor run-to-run control system with state and model parameter estimation
WO2005114338A1 (en) * 2004-05-20 2005-12-01 Mcmaster University Method for controlling the appearance of products and process performance by image analysis
CN101872182A (en) * 2010-05-21 2010-10-27 杭州电子科技大学 A Batch Process Monitoring Method Based on Recursive Nonlinear Partial Least Squares
CN102313578A (en) * 2011-08-04 2012-01-11 广州市香港科大霍英东研究院 Mechanical seal online monitoring system
CN102436524A (en) * 2011-10-19 2012-05-02 清华大学 Fuzzy reasoning method for soft fault diagnosis for analog circuit
CN102431136A (en) * 2011-09-16 2012-05-02 广州市香港科大霍英东研究院 Multi-stage batch process stage division method based on multidirectional principal component analysis method
CN102662321A (en) * 2012-03-23 2012-09-12 清华大学 Online updating method of principal component analysis monitoring model
CN102789676A (en) * 2012-08-10 2012-11-21 杭州电子科技大学 Method for designing industrial alarm on basis of alarm evidence fusion
CN103116306A (en) * 2013-02-05 2013-05-22 浙江大学 Automatic stepping type ordered time interval dividing method
CN103245759A (en) * 2013-03-28 2013-08-14 杭州电子科技大学 Product quality monitoring method based on autoregression total projection to latent structures (T-PLS)
CN103279123A (en) * 2013-05-21 2013-09-04 沈阳化工大学 Method of monitoring faults in sections for intermittent control system
CN105004542A (en) * 2015-07-15 2015-10-28 浙江中烟工业有限责任公司 Online monitoring and fault diagnosing method for mixing and flavouring process of cigarette filament production based on principal component analysis
CN105139086A (en) * 2015-08-13 2015-12-09 杭州电子科技大学 Track profile irregularity amplitude estimation method employing optimal belief rules based inference
CN105150482A (en) * 2015-08-19 2015-12-16 广州市香港科大霍英东研究院 Online fault detection method for non-return ring
CN106022366A (en) * 2016-07-04 2016-10-12 杭州电子科技大学 Rotary mechanical equipment fault diagnosis method based on neighbor evidence fusion
CN107346122A (en) * 2017-05-27 2017-11-14 四川用联信息技术有限公司 Improve the manufacturing process multivariate quality diagnostic classification device of fuzzy support vector machine
CN107490964A (en) * 2017-08-17 2017-12-19 杭州电子科技大学 A kind of rotating machinery fault feature reduction method of feature based evidence discretization
CN108257365A (en) * 2018-01-29 2018-07-06 杭州电子科技大学 A kind of industrial alarm designs method based on global nonspecific evidence dynamic fusion
CN108960332A (en) * 2018-07-11 2018-12-07 杭州电子科技大学 A kind of on-line monitoring method based on multidirectional the analysis of main elements
CN109145972A (en) * 2018-08-09 2019-01-04 杭州电子科技大学 A kind of watercraft electric propulsion system frequency converter alarm design method

Family Cites Families (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US8489360B2 (en) * 2006-09-29 2013-07-16 Fisher-Rosemount Systems, Inc. Multivariate monitoring and diagnostics of process variable data

Patent Citations (20)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
WO2003096130A1 (en) * 2001-10-23 2003-11-20 Brooks-Pri Automation, Inc. Semiconductor run-to-run control system with state and model parameter estimation
WO2005114338A1 (en) * 2004-05-20 2005-12-01 Mcmaster University Method for controlling the appearance of products and process performance by image analysis
CN101872182A (en) * 2010-05-21 2010-10-27 杭州电子科技大学 A Batch Process Monitoring Method Based on Recursive Nonlinear Partial Least Squares
CN102313578A (en) * 2011-08-04 2012-01-11 广州市香港科大霍英东研究院 Mechanical seal online monitoring system
CN102431136A (en) * 2011-09-16 2012-05-02 广州市香港科大霍英东研究院 Multi-stage batch process stage division method based on multidirectional principal component analysis method
CN102436524A (en) * 2011-10-19 2012-05-02 清华大学 Fuzzy reasoning method for soft fault diagnosis for analog circuit
CN102662321A (en) * 2012-03-23 2012-09-12 清华大学 Online updating method of principal component analysis monitoring model
CN102789676A (en) * 2012-08-10 2012-11-21 杭州电子科技大学 Method for designing industrial alarm on basis of alarm evidence fusion
CN103116306A (en) * 2013-02-05 2013-05-22 浙江大学 Automatic stepping type ordered time interval dividing method
CN103245759A (en) * 2013-03-28 2013-08-14 杭州电子科技大学 Product quality monitoring method based on autoregression total projection to latent structures (T-PLS)
CN103279123A (en) * 2013-05-21 2013-09-04 沈阳化工大学 Method of monitoring faults in sections for intermittent control system
CN105004542A (en) * 2015-07-15 2015-10-28 浙江中烟工业有限责任公司 Online monitoring and fault diagnosing method for mixing and flavouring process of cigarette filament production based on principal component analysis
CN105139086A (en) * 2015-08-13 2015-12-09 杭州电子科技大学 Track profile irregularity amplitude estimation method employing optimal belief rules based inference
CN105150482A (en) * 2015-08-19 2015-12-16 广州市香港科大霍英东研究院 Online fault detection method for non-return ring
CN106022366A (en) * 2016-07-04 2016-10-12 杭州电子科技大学 Rotary mechanical equipment fault diagnosis method based on neighbor evidence fusion
CN107346122A (en) * 2017-05-27 2017-11-14 四川用联信息技术有限公司 Improve the manufacturing process multivariate quality diagnostic classification device of fuzzy support vector machine
CN107490964A (en) * 2017-08-17 2017-12-19 杭州电子科技大学 A kind of rotating machinery fault feature reduction method of feature based evidence discretization
CN108257365A (en) * 2018-01-29 2018-07-06 杭州电子科技大学 A kind of industrial alarm designs method based on global nonspecific evidence dynamic fusion
CN108960332A (en) * 2018-07-11 2018-12-07 杭州电子科技大学 A kind of on-line monitoring method based on multidirectional the analysis of main elements
CN109145972A (en) * 2018-08-09 2019-01-04 杭州电子科技大学 A kind of watercraft electric propulsion system frequency converter alarm design method

Non-Patent Citations (3)

* Cited by examiner, † Cited by third party
Title
A belief rule-based evidence updating method for industrial alarm system design;Xiaobin Xu,Haiyang Xu,Chenglin Wen,Jianning Li,Pingzhi Hou;《Elsevier Science》;20181231;第81卷;73-84 *
基于置信规则库推理的证据滤波报警器设计;徐海洋，徐晓滨，文成林，李建宁;《山东科技大学学报》;20170831;第36卷(第4期);45-50 *
递推PLS方法及其在过程监控中的应用;苏金明，李春富，孙如田，侯平智，王建中;《机电工程》;20100331;第27卷(第3期);92-95 *

Also Published As

Publication number Publication date
CN109932908A (en) 2019-06-25

Similar Documents

Publication Publication Date Title
CN103336507A (en) Statistical modeling and on-line monitoring method based on multimodality collaboration time frame automatic division
CN103116306B (en) Automatic stepping type ordered time interval dividing method
CN109460890B (en) An intelligent self-healing method based on reinforcement learning and control performance monitoring
TW201715320A (en) Tool failure analysis using space-distorted similarity
CN118305980B (en) Control method and device of injection molding machine
CN117644625B (en) Intelligent injection molding method based on machine vision
CN109932908B (en) A multi-directional principal component analysis process monitoring method based on alarm reliability fusion
Sun et al. Modelling and prediction of injection molding process using copula entropy and multi-output SVR
CN103310095A (en) Intermittent process quality index soft measuring method
CN116401570B (en) Intelligent processing system for printing quality monitoring big data
EP3520987A1 (en) Method for monitoring and control of the injection process of plastic materials
Yassaie et al. Data-driven fault classification in large-scale industrial processes using reduced number of process variables
CN116834244A (en) Image monitoring and alarm system for injection mold and method thereof
Zhou et al. Objectives, challenges, and prospects of batch processes: Arising from injection molding applications
US20240131765A1 (en) Reinforcement Learning Method, Non-Transitory Computer Readable Recording Medium, Reinforcement Learning Device and Molding Machine
US12135536B2 (en) Method, system and computer program product for monitoring a shaping process
CN108537249B (en) Industrial process data clustering method for density peak clustering
Yan et al. Automated process monitoring in injection molding via representation learning and setpoint regression
US20220402183A1 (en) Learning Model Generation Method, Non-Transitory Computer Readable Recording Medium, Set Value Determination Device, Molding Machine, and Molding Apparatus System
Ding et al. Novel deep learning based soft sensor feature extraction for part weight prediction in injection molding processes
CN114364503A (en) Artificial intelligence-based injection molding system and method for generating molding conditions in injection molding system
CN109177101A (en) A kind of injection molding machine batch process fault detection method
JPH068298A (en) Method for characterizing performance of injection molding process
CN113762609B (en) Product quality prediction method and device
US20240326306A1 (en) Dataset Creation Method, Learning Model Generation Method, Non-Transitory Computer Readable Recording Medium, and Dataset Creation Device

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