CN111259730A - State monitoring method and system based on multivariate state estimation - Google Patents

State monitoring method and system based on multivariate state estimation Download PDF

Info

Publication number
CN111259730A
CN111259730A CN201911413592.9A CN201911413592A CN111259730A CN 111259730 A CN111259730 A CN 111259730A CN 201911413592 A CN201911413592 A CN 201911413592A CN 111259730 A CN111259730 A CN 111259730A
Authority
CN
China
Prior art keywords
equipment
importance
feature
calculating
state estimation
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
CN201911413592.9A
Other languages
Chinese (zh)
Other versions
CN111259730B (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 AIMS Intelligent Technology Co Ltd
Original Assignee
Hangzhou AIMS Intelligent Technology Co Ltd
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 AIMS Intelligent Technology Co Ltd filed Critical Hangzhou AIMS Intelligent Technology Co Ltd
Priority to CN201911413592.9A priority Critical patent/CN111259730B/en
Publication of CN111259730A publication Critical patent/CN111259730A/en
Application granted granted Critical
Publication of CN111259730B publication Critical patent/CN111259730B/en
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F2218/00—Aspects of pattern recognition specially adapted for signal processing
    • G06F2218/08—Feature extraction
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01H—MEASUREMENT OF MECHANICAL VIBRATIONS OR ULTRASONIC, SONIC OR INFRASONIC WAVES
    • G01H17/00—Measuring mechanical vibrations or ultrasonic, sonic or infrasonic waves, not provided for in the preceding groups
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01M—TESTING STATIC OR DYNAMIC BALANCE OF MACHINES OR STRUCTURES; TESTING OF STRUCTURES OR APPARATUS, NOT OTHERWISE PROVIDED FOR
    • G01M99/00—Subject matter not provided for in other groups of this subclass
    • Y—GENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
    • Y02—TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
    • Y02P—CLIMATE CHANGE MITIGATION TECHNOLOGIES IN THE PRODUCTION OR PROCESSING OF GOODS
    • Y02P90/00—Enabling technologies with a potential contribution to greenhouse gas [GHG] emissions mitigation
    • Y02P90/02—Total factory control, e.g. smart factories, flexible manufacturing systems [FMS] or integrated manufacturing systems [IMS]

Abstract

The invention relates to the field of fault diagnosis of rotating equipment, and discloses a state monitoring method and a state monitoring system based on multivariate state estimation, which comprise the following steps: A) acquiring equipment working condition parameter signals and vibration signals to obtain an m x n dimensional characteristic parameter matrix; B) calculating the importance of each feature by adopting a GBDT method, and determining the feature selection number q according to the importance; C) calculating the state estimation value of the equipment by adopting a Mahalanobis distance method; D) and establishing an equipment abnormity judgment mechanism by adopting an exponential weighted average method, and judging the state condition of the equipment. The method can effectively select key parameters capable of reflecting the equipment state, reduce the characteristic dimension, reduce the state estimation time, avoid the interference of invalid characteristics, adopt the Mahalanobis distance method to fuse multi-dimensional characteristic parameters into an estimation value, carry out state estimation through the exponential weighted average method, reduce the fluctuation caused by random errors, absorb the capacity of instantaneous burst, and have strong timeliness and high accuracy.

Description

State monitoring method and system based on multivariate state estimation
Technical Field
The invention relates to the field of fault diagnosis of rotating equipment, in particular to a state monitoring method and system based on multivariate state estimation.
Background
In recent years, the domestic aerospace, rail transit, electric power and other industrial production is developed towards systematization, digitization and automation, a mechanical system is gradually complicated, and requirements on technological parameter requirements and equipment reliability are increasingly improved. The rotating equipment bears long-time running work, and is very easy to be abnormal, so that the production efficiency and the equipment performance are further influenced. Therefore, the rotating machine is monitored, the abnormity of the rotating equipment is monitored in time, the early abnormity judgment of the rotating equipment is realized, the direct economic loss caused by the accidental shutdown of the automatic production line of a factory due to faults and maintenance is avoided, the near-zero out-of-plan shutdown and carefree production are realized, and the production efficiency is improved.
When the rotating equipment is in an abnormal condition, the temperature is increased, the vibration index is increased, and a frequency domain generates a fault frequency band. The conventional rotating equipment abnormity judgment method is a method for judging a threshold value by a single characteristic index, an effective value of a vibration signal is compared and calculated according to an ISO2372 standard threshold value, and if the effective value exceeds the threshold value, the equipment is judged to be in an abnormal condition. ISO2372 is a general threshold judgment standard, and is used for judging that the false alarm rate and the missing report rate are high due to different vibration conditions in a scene without application. Under general conditions, the general threshold is larger, and when the vibration exceeds the general threshold, the rotary equipment has irreversible damage; in addition, the existing method can only judge the equipment state at the current moment, and cannot realize prediction, and in the existing state monitoring system, the false alarm rate and the missing alarm rate of the early-stage abnormality alarm mechanism of the equipment are higher. Therefore, a feature capable of reflecting the operation condition of the equipment needs to be constructed according to the details of the equipment, further, the state evaluation of the equipment is realized based on multi-feature parameters, the abnormal early warning of the equipment can be realized, and the false alarm rate of the abnormal judgment of the rotating equipment is reduced.
Disclosure of Invention
The method aims to solve the problems that in the existing abnormal judgment method for the rotating equipment, the threshold value is fixed, the size of the threshold value cannot be adjusted in a self-adaptive mode aiming at the working condition of the equipment, the characteristic parameters are selected according to general experience, and the proper characteristic parameters are not selected aiming at the equipment; the invention can effectively select key parameters capable of reflecting the equipment state, reduce the characteristic dimension, reduce the state estimation time, avoid the interference of invalid characteristics, improve the accuracy of state estimation, adopt the Mahalanobis distance method to fuse the multidimensional characteristic parameters into an estimation value, and carry out state estimation by the exponential weighted average method, can reduce the fluctuation caused by random errors, and can absorb the capacity of instantaneous burst, and has strong timeliness.
In order to achieve the purpose, the invention adopts the following technical scheme:
a state monitoring method based on multivariate state estimation comprises the following steps:
A) acquiring equipment working condition parameter signals and vibration signals to obtain m groups of equipment signal samples, respectively extracting features of each group of equipment signal samples to obtain n1 features of each group of equipment vibration signal samples and n2 features of the working condition parameter signals to obtain an m multiplied by n dimensional feature parameter matrix, wherein n is n1+ n 2;
B) calculating the importance of each feature by adopting a GBDT method, determining the feature selection number q according to the importance, and obtaining a feature parameter matrix after feature selection;
C) calculating the state estimation value of the equipment by adopting a Mahalanobis distance method;
D) and establishing an equipment abnormity judgment mechanism by adopting an exponential weighted average method, and judging the state condition of the equipment.
The invention can extract the multidimensional characteristics reflecting the running state of the equipment to the maximum extent by extracting the characteristics of the equipment vibration signal in time domain, frequency domain and time-frequency domain. By the GBDT method, characteristic parameters which can truly reflect the running state of the equipment can be extracted. And the multi-dimensional characteristic parameter estimation of the running state of the equipment is realized by a Mahalanobis distance method, and the abnormity early warning of the equipment is realized. An abnormity alarm mechanism is established according to the state estimation value, and the abnormity early warning false alarm rate can be reduced and the accuracy rate of judging the state condition of the equipment is high through an exponential weighted average mode.
Further, in the step a), feature extraction is performed on each group of device signal samples, so as to obtain n1 features of each group of device vibration signal samples and n2 features of the working condition parameter signals, and obtain an m × n-dimensional feature parameter matrix, including the steps of:
A1) f1 time domain features of each group of device vibration signal samples are extracted;
A2) f2 frequency domain features of each group of device vibration signal samples are extracted;
A3) wavelet packet transformation is carried out on each group of vibration signals to obtain P component signals, and the frequency band energy and the energy entropy of the P component signals are respectively calculated; the band energy is calculated by the formula
Figure RE-GDA0002417421820000021
The calculation formula of the energy entropy is
Figure RE-GDA0002417421820000024
Wherein i represents the ith component signal, N represents the number of component signal data points, j represents the number of component signal data points,
Figure RE-GDA0002417421820000022
a band energy representing an i-th component signal of the k-th group of the set vibration signal samples;
A4) obtaining the energy sum of the vibration signals of the kth group of devices
Figure RE-GDA0002417421820000023
A5) Obtaining n1 features of each group of equipment vibration signal samples, wherein the n1 features comprise time domain features, frequency band energy, energy entropy and energy sum, the n2 working condition parameter signal features comprise information such as current, voltage and rotating speed, a characteristic parameter matrix is obtained, and a characteristic parameter matrix is obtained
Figure RE-GDA0002417421820000031
According to the invention, through carrying out wavelet packet conversion on the equipment vibration signal, not only can equal frequency band division be realized, but also energy in each frequency band can be extracted, and multidimensional characteristics reflecting the operation state of the equipment can be extracted to the maximum extent.
Further, in the step B), the importance of each feature is calculated by using a gbdt (gradient Boosting decision tree) method, and the number q of feature choices is determined according to the importance, so as to obtain a feature parameter matrix after feature selection, including:
B1) the importance degree of the n characteristics on each tree is calculated respectively, and the calculation formula of the importance degree of the jth characteristic in a single tree is
Figure RE-GDA0002417421820000032
Wherein L is the number of leaf nodes of the tree, L-1 is the number of non-leaf nodes of the tree, vtIs a feature associated with the node t,
Figure RE-GDA0002417421820000033
is the reduction of the squared loss after splitting of node t;
B2) the average value of the importance of the n characteristics on each tree is calculated, and the average value of the importance of the jth characteristic on each tree is calculated according to the formula
Figure RE-GDA0002417421820000034
Wherein M is the number of trees;
B3) calculating the average value of the importance of n characteristics divided by the maximum importance value to obtain the importance value of each characteristic, wherein the calculation formula of the importance value of the jth characteristic is
Figure RE-GDA0002417421820000035
Represents the maximum importance value of the jth feature;
B4) the importance values are sequentially arranged from large to small to obtain an arranged importance arrangement set { I }1,I2,…,InDetermining the number q of feature choices according to the importance;
B5) according to the characteristic number q, obtaining a selected characteristic parameter matrix
Figure RE-GDA0002417421820000036
The invention adopts a GBDT method to calculate the importance of each feature, determines the number q of feature selection according to the importance, reduces the dimension of a feature parameter matrix, obtains the first q features with the maximum importance, and reconstructs the feature parameter matrix according to the first q features with the maximum importance to obtain a selected feature parameter matrix V. By adopting the GBDT importance degree calculation method, the key parameters capable of reflecting the equipment state can be effectively selected, and the parameter selection method can reduce the characteristic dimension and reduce the state estimation time on one hand; and on the other hand, the interference of invalid features is avoided, and the accuracy of state estimation can be improved.
Further, the step C) of calculating the state estimation value of the device by using the mahalanobis distance method includes the steps of:
C1) calculating the expected value of each dimensional matrix according to the selected characteristic parameter matrix V
μ＝[μ1,μ2,…,μq]T；
C2) Calculating a covariance matrix sigma of a characteristic parameter matrix according to the expected value mu;
C3) and calculating the Mahalanobis distance D according to the covariance matrix sigma of the characteristic parameter matrix, and taking the Mahalanobis distance D as the state estimation value of the equipment.
Mahalanobis distance (Mahalanobis distance) represents the covariance distance of data, and is an effective method for calculating the similarity between two unknown sample sets, and unlike euclidean distance, which considers the relationship between various characteristics, Mahalanobis distance is not affected by dimension, Mahalanobis distance between two points is independent of the unit of measurement of the original data, and Mahalanobis distance between two points calculated from normalized data and centralized data (i.e., the difference between the original data and the mean) is the same. Mahalanobis distance can also exclude variablesInterference of correlation between them. For a mean of mu, the covariance matrix is sigma-1Of a multivariate vector x having a Mahalanobis distance of
Figure BDA0002350610880000042
∑-1An inverse matrix corresponding to the covariance matrix is represented.
Further, in step D), an equipment abnormality judgment mechanism is established by using an exponential weighted moving average method to judge the equipment state condition, including the steps of:
D1) acquiring a Mahalanobis distance set of normal sample data in a period of time, and calculating a mean value mu of the Mahalanobis distance set1Sum variance σ1Determining the space range of the normal sample as [ mu ] according to the 3 sigma rule1-3σ，μ1+3σ]；
D2) Obtaining a new incoming observation sample xtCalculating an observation sample xtMahalanobis distance D oft；
D3) Calculating an observation sample estimation value D after weighted moving averaget＝w0Dt+w1Dt-1+W3Dt-2+…+WxDt-nWherein D ist-nRepresenting an observed sample xt-nMahalanobis distance, wnWeight, w, representing the mahalanobis distance at time t-n0+w1+…wn＝1；
D4) Judging the device estimation value DtWhether in normal sample space [ mu ]1-3σ，μ1+3σ]If yes, indicating that the equipment is normal; if not, the equipment is abnormal.
The exponential weighted average estimation method adopted by the invention can track the ability of the actual data to change suddenly, and has the ability of reducing short-term fluctuation and retaining long-term development trend.
A state monitoring system based on multivariate state estimation comprises a data acquisition module, a feature extraction module, a feature selection module, a state estimation module and an abnormality judgment module;
the data acquisition module is used for acquiring equipment working condition parameter signals and vibration signals;
the characteristic extraction module is used for extracting the characteristics of the equipment working condition parameter signals and the vibration signals acquired by the data acquisition module, wherein the characteristics comprise time domain characteristics, frequency band energy, energy entropy and energy sum;
the characteristic selection module is used for calculating the importance of each characteristic by adopting a GBDT method, calculating the importance of the equipment signal characteristics extracted from the characteristic extraction module, and selecting the equipment signal characteristics according to the importance to obtain a characteristic parameter matrix after the characteristics are selected;
the state estimation module is used for calculating a state estimation value of the equipment by adopting a Mahalanobis distance method to obtain the state estimation value of the equipment;
and the abnormality judgment module is used for establishing an equipment abnormality judgment mechanism by adopting an exponential weighted moving average method and judging the equipment state condition.
The invention has the following beneficial effects: the invention adopts a GBDT importance degree calculation method, can effectively select key parameters capable of reflecting the equipment state, and can reduce the feature dimension and the state estimation time through the parameter selection method; on the other hand, interference of invalid features is avoided, and the accuracy of state estimation can be improved. By adopting the Mahalanobis distance method, the multi-dimensional characteristic parameters are fused into an estimated value, and the state estimation is carried out by an exponential weighted average method, so that the fluctuation caused by random errors can be reduced, the instantaneous burst capacity can be absorbed, and the timeliness is strong.
Drawings
FIG. 1 is a schematic flow chart of an embodiment.
FIG. 2 is a graph of Mahalanobis distance calculations for the first 100 sample groups of the example.
FIG. 3 is a graph of the results of evaluation of the last 300 groups of samples according to one embodiment.
FIG. 4 is a diagram of the evaluation results using the 37-dimensional features according to one embodiment.
Detailed Description
The invention is further described with reference to the following detailed description and accompanying drawings.
In a first embodiment, as shown in fig. 1, a state monitoring method based on multivariate state estimation, taking a certain type of motor as an example, the rotation speed of the motor is 2500rpm, includes the steps of:
A) collecting working condition parameter signals and vibration signals of equipment, wherein the sampling frequency is 25600Hz, m groups of vibration signal samples are obtained, m is 400,
respectively extracting the characteristics of each group of equipment signal samples to obtain n1 characteristics of each group of equipment vibration signal samples and n2 characteristics of working condition parameter signals to obtain an m × n-dimensional characteristic parameter matrix, wherein n is n1+ n2, and the method comprises the following steps of:
A1) extracting 15 time domain characteristics of each group of equipment vibration signal samples, wherein the 15 time domain characteristics comprise a mean value, an absolute mean value, a variance, a standard deviation, a maximum value, a minimum value, a peak-to-peak value, an effective value, a skewness, a peak index, a pulse index, a margin index, a waveform index, a skewness index, a kurtosis index and a kurtosis index;
A2) extracting 4 frequency domain characteristics of each group of equipment vibration signal samples, wherein the 4 frequency domain characteristics comprise a spectrum effective value, a spectrum variance, a spectrum mean value and a spectrum center;
A3) carrying out 3-layer wavelet packet transformation on each group of vibration signals to obtain 8 component signals, and respectively calculating the frequency band energy and the energy entropy of the 8 component signals; the band energy is calculated by the formula
Figure RE-GDA0002417421820000061
The calculation formula of the energy entropy is
Figure RE-GDA0002417421820000069
Wherein i represents the ith component signal, N represents the number of component signal data points, j represents the number of component signal data points,
Figure RE-GDA0002417421820000062
a band energy representing an i-th component signal of the k-th group of the set vibration signal samples;
A4) obtaining the energy sum of the vibration signals of the kth group of devices
Figure RE-GDA0002417421820000063
A5) Obtaining n1 features of each group of device vibration signal samples and n2 features of the working condition parameter signals, wherein n is 37, the features comprise 15 time domain features, 4 frequency domain features, 8 frequency band energies, 8 energy entropies and 1 energy sum, and obtaining a 400 x 37-dimensional feature parameter matrix
Figure RE-GDA0002417421820000064
B) The method adopts a GBDT method to calculate the importance of each feature, determines the number q of feature choices according to the importance, and comprises the following steps:
B1) the importance of each tree is calculated for 37 characteristics, and the importance of the jth characteristic in a single tree is calculated by the formula
Figure RE-GDA0002417421820000065
Wherein L is the number of leaf nodes of the tree, L-1 is the number of non-leaf nodes of the tree, vtIs a feature associated with the node t,
Figure RE-GDA0002417421820000066
is the reduction of the squared loss after splitting of node t;
B2) the average value of the importance of the 37 th feature on each tree is calculated by the formula
Figure RE-GDA0002417421820000067
Wherein M is the number of trees;
B3) calculating the average value of the importance of 37 features divided by the maximum importance value to obtain the importance value of each feature, wherein the calculation formula of the importance value of the jth feature is
Figure RE-GDA0002417421820000068
Represents the maximum importance value of the jth feature;
B4) the importance values are sequentially arranged from large to small to obtain an arranged importance arrangement set { I }1,I2,…,InAnd determining the feature selection number q according to the importance. And (4) selecting the features with the importance degree exceeding 50%, namely selecting a waveform index, a spectrum gravity center, energy entropy values of all frequency bands and energy entropy of all frequency bands according to the importance degree sequence, wherein 11 features are recorded in total, namely q is 11.
B5) Obtaining the selected characteristic parameter matrix according to the characteristic number
Figure BDA0002350610880000077
C) The method for calculating the state estimation value of the equipment by adopting the Mahalanobis distance method comprises the following steps of:
C1) selecting the first 100 groups of samples from the characteristic parameter matrix V, and calculating the expected value mu of each dimensional matrix as [ mu ]1，μ2，…，μq]TThe expected value results are shown in table 1.
TABLE 1 corresponding average value table of characteristic parameters of the first 100 groups of samples
Characteristic parameter Mean value Characteristic parameter Mean value Characteristic parameter Mean value
Center of gravity of music score 6154.1 Waveform index 1.3328 H 0.8953
h1 0.1070 h2 0.0985 h3 0.1221
h4 0.1172 h5 0.0982 h6 0.1023
h7 0.1260 h8 0.1241
TABLE 1
C2) Calculating a covariance matrix sigma of a characteristic parameter matrix according to the expected value mu;
C3) and calculating the Mahalanobis distance D according to the covariance matrix sigma of the characteristic parameter matrix, and taking the Mahalanobis distance D as the state estimation value of the equipment. Figure 2 is the mahalanobis distance results for the first 100 sets of signals.
D) An exponential weighted average method is adopted to establish an equipment abnormity judgment mechanism and judge the equipment state condition, and the method comprises the following steps:
D1) acquiring the mahalanobis distance set of the first 100 normal sample data, and calculating the mean value mu of the mahalanobis distance of the first 100 groups of data1Sum variance σ1Determining the space range of the normal sample as [ mu ] according to the 3 sigma rule1-3σ，μ1+3σ](ii) a According to the calculation, mu is obtained1＝7.315，σ1At 5.096, the normal sample spatial range is [0,22.6037 ]]。
D2) Calculating the Mahalanobis distance D from the 101 th group to the 400 th group of observation samplestObtaining a new incoming observation sample xtCalculating an observation sample xtMahalanobis distance D oft；
D3) Calculating an observation sample estimation value D after weighted moving averaget＝w0Dt+w1Dt-1+W2Dt-2+…+WxDt-nWherein D ist-nRepresenting an observed sample xt-nMahalanobis distance, wnWeight, w, representing the mahalanobis distance at time t-n0+w1+…wn＝1；
D4) Judging the device estimation value DtWhether in normal sample space [0,22.6037]If yes, indicating that the equipment is normal; if notThen the device is abnormal. Figure 3 is a graph showing the mahalanobis distance estimate calculated for 300 samples, which exceeds the alarm line at 160 samples, indicating an abnormal condition of the equipment.
Fig. 4 is a graph of the results of feature evaluation with all 37-dimensional features selected for state estimation, with the device estimate exceeding the alarm line at group 178 samples later than the alarm time at which the device was found to be abnormal according to the present invention at group 160 samples. Therefore, the invention can find the abnormal condition of the equipment in advance and realize the advance alarm.
According to the invention, through wavelet packet transformation of the equipment vibration signal, not only can self-adaptive frequency band division be realized, but also energy in each frequency band can be extracted, and through characteristic extraction on time domain, frequency domain and time-frequency domain of the equipment vibration signal, multidimensional characteristics reflecting the operation state of the equipment can be extracted to the maximum extent. Key parameters capable of reflecting the equipment state can be effectively selected by the GBDT method, and the characteristic dimension can be reduced and the state estimation time can be reduced by the parameter selection method; on the other hand, interference of invalid features is avoided, and the accuracy of state estimation can be improved. And the multi-dimensional characteristic parameter estimation of the running state of the equipment is realized by a Mahalanobis distance method, and the abnormity early warning of the equipment is realized. An abnormity alarm mechanism is established according to the state estimation value, and the abnormity early warning and false alarm rate can be reduced in an exponential weighted average mode, so that the accuracy for judging the state condition of the equipment is high. The invention adopts the Mahalanobis distance method, integrates the multidimensional characteristic parameters into an estimated value, and carries out state estimation by an exponential weighted average method, thereby reducing the fluctuation caused by random errors, absorbing the capacity of instantaneous burst and having stronger timeliness.
Although the embodiments of the present invention have been described with reference to the accompanying drawings, it is not intended to limit the scope of the present invention, and it should be understood by those skilled in the art that various modifications and variations can be made without inventive efforts by those skilled in the art based on the technical solution of the present invention.

Claims (6)

1. A state monitoring method based on multivariate state estimation is characterized by comprising the following steps:
A) acquiring equipment working condition parameter signals and vibration signals to obtain m groups of equipment signal samples, respectively extracting features of each group of equipment signal samples to obtain n1 features of each group of equipment vibration signal samples and n2 features of the working condition parameter signals to obtain an m multiplied by n dimensional feature parameter matrix, wherein n is n1+ n 2;
B) calculating the importance of each feature by adopting a GBDT method, determining the feature selection number q according to the importance, and obtaining a feature parameter matrix after feature selection;
C) calculating the state estimation value of the equipment by adopting a Mahalanobis distance method;
D) and establishing an equipment abnormity judgment mechanism by adopting an exponential weighted average method, and judging the state condition of the equipment.
2. The multivariate state estimation-based state monitoring method as claimed in claim 1, wherein in the step a), feature extraction is respectively performed on each group of equipment signal samples, so as to obtain n1 features of each group of equipment vibration signal samples and n2 features of the working condition parameter signals, and obtain a dimensional feature parameter matrix, and the method comprises the following steps:
A1) f1 time domain features of each group of device vibration signal samples are extracted;
A2) f2 frequency domain features of each group of device vibration signal samples are extracted;
A3) wavelet packet transformation is carried out on each group of vibration signals to obtain P component signals, and the frequency band energy and the energy entropy of the P component signals are respectively calculated; the band energy is calculated by the formula
Figure RE-FDA0002417421810000011
The calculation formula of the energy entropy is
Figure RE-FDA0002417421810000012
Wherein i represents the ith component signal, N represents the number of component signal data points, j represents the number of component signal data points,
Figure RE-FDA0002417421810000013
a band energy representing an i-th component signal of the k-th group of the set vibration signal samples;
A4) obtaining the energy sum of the vibration signals of the kth group of devices
Figure RE-FDA0002417421810000014
A5) Obtaining n1 features of each group of equipment vibration signal samples, wherein the n1 features comprise time domain features, frequency band energy, energy entropy and energy sum, the n2 working condition parameter signal features comprise information such as current, voltage and rotating speed, and a feature parameter matrix is obtained
Figure RE-FDA0002417421810000015
3. The multivariate state estimation-based state monitoring method as claimed in claim 1 or 2, wherein the step B) of calculating the importance of each feature by using the GBDT method, determining the number q of feature choices according to the importance, and obtaining the feature parameter matrix after feature choice comprises:
B1) the importance degree of the n characteristics on each tree is calculated respectively, and the calculation formula of the importance degree of the jth characteristic in a single tree is
Figure RE-FDA0002417421810000021
Wherein L is the number of leaf nodes of the tree, L-1 is the number of non-leaf nodes of the tree, vtIs a feature associated with the node t,
Figure RE-FDA0002417421810000022
is the reduction of the square loss after splitting of node t;
B2) the average value of the importance of the n characteristics on each tree is calculated, and the average value of the importance of the jth characteristic on each tree is calculated according to the formula
Figure RE-FDA0002417421810000023
Wherein M is the number of trees;
B3) calculating the average value of the importance of n characteristics divided by the maximum importance value to obtain the importance value of each characteristic, wherein the calculation formula of the importance value of the jth characteristic is
Figure RE-FDA0002417421810000024
Figure RE-FDA0002417421810000025
Represents the maximum importance value of the jth feature;
B4) sequentially arranging the importance values from large to small to obtain an arranged importance arrangement set I1,I2,…,InDetermining the number q of feature choices according to the importance;
B5) according to the characteristic number q, obtaining a selected characteristic parameter matrix
Figure RE-FDA0002417421810000026
4. A multivariate state estimation based state monitoring method as defined in claim 3, wherein the mahalanobis distance method is used in step C) to calculate the state estimation value of the device, comprising the steps of:
C1) calculating the expected value mu [ mu ] of each dimensional matrix according to the selected characteristic parameter matrix V1,μ2,…,μq]T；
C2) Calculating a covariance matrix sigma of a characteristic parameter matrix according to the expected value mu;
C3) and calculating the Mahalanobis distance D according to the covariance matrix sigma of the characteristic parameter matrix, and taking the Mahalanobis distance D as the state estimation value of the equipment.
5. The multivariate state estimation-based state monitoring method as claimed in claim 1 or 4, wherein the step D) adopts an exponential weighted moving average method to establish a device abnormality judgment mechanism to judge the state of the device, and comprises the steps of:
D1) acquiring a Mahalanobis distance set of normal sample data in a period of time, and calculating a mean value mu of the Mahalanobis distance set1Sum variance σ1Determining the space range of the normal sample as [ mu ] according to the 3 sigma rule1-3σ，μ1+3σ]；
D2) Obtaining a new incoming observation sample xtCalculating an observation sample xtMahalanobis distance D oft；
D3) Calculating an observation sample estimation value D after weighted moving averaget＝w0Dt+w1Dt-1+w2Dt-2+…+wnDt-nWherein D ist-nRepresenting an observed sample xt-nMahalanobis distance, wnWeight, w, representing the mahalanobis distance at time t-n0+w1+…wn＝1；
D4) Judging the device estimation value DtWhether in normal sample space [ mu ]1-3σ，μ1+3σ]If yes, indicating that the equipment is normal; if not, the equipment is abnormal.
6. A multivariate state estimation-based state monitoring system, which is suitable for the multivariate state estimation-based state monitoring method according to any one of claims 1 to 5, and is characterized by comprising a data acquisition module, a feature extraction module, a feature selection module, a state estimation module and an abnormality judgment module;
the data acquisition module is used for acquiring equipment working condition parameter signals and vibration signals;
the characteristic extraction module is used for extracting characteristics of the equipment working condition parameter signals and the vibration signals acquired by the data acquisition module, wherein the characteristics comprise time domain characteristics, frequency band energy, energy entropy and energy sum;
the characteristic selection module is used for calculating the importance of each characteristic by adopting a GBDT method, calculating the importance of the equipment signal characteristics extracted from the characteristic extraction module, and selecting the equipment signal characteristics according to the importance to obtain a characteristic parameter matrix after the characteristics are selected;
the state estimation module is used for calculating a state estimation value of the equipment by adopting a Mahalanobis distance method to obtain the state estimation value of the equipment;
and the abnormality judgment module is used for establishing an equipment abnormality judgment mechanism by adopting an exponential weighted moving average method and judging the equipment state condition.
CN201911413592.9A 2019-12-31 2019-12-31 State monitoring method and system based on multivariate state estimation Active CN111259730B (en)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201911413592.9A CN111259730B (en) 2019-12-31 2019-12-31 State monitoring method and system based on multivariate state estimation

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201911413592.9A CN111259730B (en) 2019-12-31 2019-12-31 State monitoring method and system based on multivariate state estimation

Publications (2)

Publication Number Publication Date
CN111259730A true CN111259730A (en) 2020-06-09
CN111259730B CN111259730B (en) 2022-08-23

Family

ID=70948519

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201911413592.9A Active CN111259730B (en) 2019-12-31 2019-12-31 State monitoring method and system based on multivariate state estimation

Country Status (1)

Country Link
CN (1) CN111259730B (en)

Cited By (9)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN112015153A (en) * 2020-09-09 2020-12-01 江南大学 System and method for detecting abnormity of sterile filling production line
CN112763056A (en) * 2020-12-29 2021-05-07 上海交大智邦科技有限公司 Method and system for online real-time monitoring and evaluation of state of numerical control machine tool system
CN113723245A (en) * 2021-08-20 2021-11-30 西安交通大学 Method, system, equipment and storage medium for monitoring running state of reciprocating compressor
CN113746701A (en) * 2021-09-03 2021-12-03 四川英得赛克科技有限公司 Data acquisition method, system, storage medium and electronic equipment
CN115062678A (en) * 2022-08-19 2022-09-16 山东能源数智云科技有限公司 Training method of equipment fault detection model, fault detection method and device
CN115500829A (en) * 2022-11-24 2022-12-23 广东美赛尔细胞生物科技有限公司 Depression detection and analysis system applied to neurology
CN115829422A (en) * 2023-02-21 2023-03-21 创银科技(南通)有限公司 Industrial equipment operation abnormal state identification method based on big data
CN117268535A (en) * 2023-11-22 2023-12-22 四川中测仪器科技有限公司 Motor rotating shaft state monitoring method based on vibration data
CN117290670A (en) * 2023-11-27 2023-12-26 南京中鑫智电科技有限公司 Transformer bushing insulation state estimation method based on enhancement filter algorithm

Citations (6)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP2005241089A (en) * 2004-02-25 2005-09-08 Mitsubishi Electric Corp Apparatus diagnosing device, refrigeration cycle device, apparatus diagnosing method, apparatus monitoring system and refrigeration cycle monitoring system
US20060210141A1 (en) * 2005-03-16 2006-09-21 Omron Corporation Inspection method and inspection apparatus
CN103674511A (en) * 2013-03-18 2014-03-26 北京航空航天大学 Mechanical wearing part performance assessment and prediction method based on EMD (empirical mode decomposition)-SVD (singular value decomposition) and MTS (Mahalanobis-Taguchi system)
CN105718876A (en) * 2016-01-18 2016-06-29 上海交通大学 Evaluation method of health states of ball screw
CN107146004A (en) * 2017-04-20 2017-09-08 浙江大学 A kind of slag milling system health status identifying system and method based on data mining
US20190035170A1 (en) * 2017-07-27 2019-01-31 Toyota Motor Engineering & Manufacturing North America, Inc. Servicing schedule method based on prediction of degradation in electrified vehicles

Patent Citations (6)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP2005241089A (en) * 2004-02-25 2005-09-08 Mitsubishi Electric Corp Apparatus diagnosing device, refrigeration cycle device, apparatus diagnosing method, apparatus monitoring system and refrigeration cycle monitoring system
US20060210141A1 (en) * 2005-03-16 2006-09-21 Omron Corporation Inspection method and inspection apparatus
CN103674511A (en) * 2013-03-18 2014-03-26 北京航空航天大学 Mechanical wearing part performance assessment and prediction method based on EMD (empirical mode decomposition)-SVD (singular value decomposition) and MTS (Mahalanobis-Taguchi system)
CN105718876A (en) * 2016-01-18 2016-06-29 上海交通大学 Evaluation method of health states of ball screw
CN107146004A (en) * 2017-04-20 2017-09-08 浙江大学 A kind of slag milling system health status identifying system and method based on data mining
US20190035170A1 (en) * 2017-07-27 2019-01-31 Toyota Motor Engineering & Manufacturing North America, Inc. Servicing schedule method based on prediction of degradation in electrified vehicles

Cited By (16)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN112015153A (en) * 2020-09-09 2020-12-01 江南大学 System and method for detecting abnormity of sterile filling production line
CN112015153B (en) * 2020-09-09 2021-06-22 江南大学 System and method for detecting abnormity of sterile filling production line
WO2022052510A1 (en) * 2020-09-09 2022-03-17 江南大学 Anomaly detection system and method for sterile filling production line
CN112763056A (en) * 2020-12-29 2021-05-07 上海交大智邦科技有限公司 Method and system for online real-time monitoring and evaluation of state of numerical control machine tool system
CN113723245A (en) * 2021-08-20 2021-11-30 西安交通大学 Method, system, equipment and storage medium for monitoring running state of reciprocating compressor
CN113723245B (en) * 2021-08-20 2023-08-18 西安交通大学 Method, system, equipment and storage medium for monitoring running state of reciprocating compressor
CN113746701B (en) * 2021-09-03 2023-01-06 四川英得赛克科技有限公司 Data acquisition method, system, storage medium and electronic equipment
CN113746701A (en) * 2021-09-03 2021-12-03 四川英得赛克科技有限公司 Data acquisition method, system, storage medium and electronic equipment
CN115062678A (en) * 2022-08-19 2022-09-16 山东能源数智云科技有限公司 Training method of equipment fault detection model, fault detection method and device
CN115500829A (en) * 2022-11-24 2022-12-23 广东美赛尔细胞生物科技有限公司 Depression detection and analysis system applied to neurology
CN115829422A (en) * 2023-02-21 2023-03-21 创银科技(南通)有限公司 Industrial equipment operation abnormal state identification method based on big data
CN115829422B (en) * 2023-02-21 2024-01-02 北京瀚海蓝山智能科技有限公司 Industrial equipment operation abnormal state identification method based on big data
CN117268535A (en) * 2023-11-22 2023-12-22 四川中测仪器科技有限公司 Motor rotating shaft state monitoring method based on vibration data
CN117268535B (en) * 2023-11-22 2024-01-26 四川中测仪器科技有限公司 Motor rotating shaft state monitoring method based on vibration data
CN117290670A (en) * 2023-11-27 2023-12-26 南京中鑫智电科技有限公司 Transformer bushing insulation state estimation method based on enhancement filter algorithm
CN117290670B (en) * 2023-11-27 2024-01-26 南京中鑫智电科技有限公司 Transformer bushing insulation state estimation method based on enhancement filter algorithm

Also Published As

Publication number Publication date
CN111259730B (en) 2022-08-23

Similar Documents

Publication Publication Date Title
CN111259730B (en) State monitoring method and system based on multivariate state estimation
CN111089726B (en) Rolling bearing fault diagnosis method based on optimal dimension singular spectrum decomposition
CN106404399B (en) Method for Bearing Fault Diagnosis based on self-adaptive redundant Lifting Wavelet packet decomposition tree
CN111044902B (en) Motor fault diagnosis method based on current and voltage signals
CN109883706A (en) A kind of rolling bearing local damage Weak fault feature extracting method
CN113569990B (en) Strong noise interference environment-oriented performance equipment fault diagnosis model construction method
CN117235557B (en) Electrical equipment fault rapid diagnosis method based on big data analysis
CN115481657A (en) Wind generating set communication slip ring fault diagnosis method based on electric signals
CN111881594A (en) Non-stationary signal state monitoring method and system for nuclear power equipment
CN111639852B (en) Real-time evaluation method and system for vibration state of hydroelectric generating set based on wavelet singular value
CN111766513B (en) Capsule network-based variable-working-condition multi-fault diagnosis method for three-phase induction motor
Jiang et al. A novel rolling-element bearing faults classification method combines lower-order moment spectra and support vector machine
Rehab et al. Bearings fault detection using hidden Markov models and principal component analysis enhanced features
CN108828419B (en) Switch room partial discharge orientation method based on likelihood estimation
CN108594156B (en) Improved current transformer saturation characteristic identification method
CN115717993A (en) Multi-channel signal self-adaptive decomposition method
CN115561575A (en) Method for distinguishing electrical abnormal state of offshore wind farm based on multi-dimensional matrix profile
Yaqub et al. Machine health monitoring based on stationary wavelet transform and fourth-order cumulants
Mezni et al. A comparative study for ball bearing fault classification using kernel-SVM with Kullback Leibler divergence selected features
Zheng et al. Wavelet packet decomposition and neural network based fault diagnosis for elevator excessive vibration
Pan et al. Fast fault diagnosis method of rolling bearings based on compression features in multi-sensor redundant observation environment
Pan et al. Fast fault diagnosis method of rolling bearings in multi-sensor measurement enviroment
Yong Bearings fault diagnosis based on HMM and fractal dimensions spectrum
KR102418118B1 (en) Apparatus and method of deep learning-based facility diagnosis using frequency synthesis
CN113469467B (en) Wind power ultra-short term prediction method and device based on band-pass filtering

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
CB03 Change of inventor or designer information
CB03 Change of inventor or designer information

Inventor after: Li Qian

Inventor after: Ni Jun

Inventor after: Liu Shulin

Inventor after: Cai Yibiao

Inventor after: Yang Haojie

Inventor after: Sun Fengcheng

Inventor before: Li Qian

Inventor before: Liu Shulin

Inventor before: Cai Yibiao

Inventor before: Yang Haojie

Inventor before: Sun Fengcheng

GR01 Patent grant
GR01 Patent grant