WO1997030329A1 - Procede de detection d'anomalies et systeme de detection d'anomalies - Google Patents

Procede de detection d'anomalies et systeme de detection d'anomalies Download PDF

Info

Publication number
WO1997030329A1
WO1997030329A1 PCT/JP1997/000367 JP9700367W WO9730329A1 WO 1997030329 A1 WO1997030329 A1 WO 1997030329A1 JP 9700367 W JP9700367 W JP 9700367W WO 9730329 A1 WO9730329 A1 WO 9730329A1
Authority
WO
WIPO (PCT)
Prior art keywords
vector
time
calculated
generating
abnormality
Prior art date
Application number
PCT/JP1997/000367
Other languages
English (en)
French (fr)
Inventor
Takaya Miyano
Original Assignee
Sumitomo Metal Industries Limited
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 Sumitomo Metal Industries Limited filed Critical Sumitomo Metal Industries Limited
Priority to AU16716/97A priority Critical patent/AU1671697A/en
Priority to EP97902662A priority patent/EP0881467A4/en
Publication of WO1997030329A1 publication Critical patent/WO1997030329A1/ja

Links

Classifications

    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01D—MEASURING NOT SPECIALLY ADAPTED FOR A SPECIFIC VARIABLE; ARRANGEMENTS FOR MEASURING TWO OR MORE VARIABLES NOT COVERED IN A SINGLE OTHER SUBCLASS; TARIFF METERING APPARATUS; MEASURING OR TESTING NOT OTHERWISE PROVIDED FOR
    • G01D1/00—Measuring arrangements giving results other than momentary value of variable, of general application
    • G01D1/18—Measuring arrangements giving results other than momentary value of variable, of general application with arrangements for signalling that a predetermined value of an unspecified parameter has been exceeded
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01D—MEASURING NOT SPECIALLY ADAPTED FOR A SPECIFIC VARIABLE; ARRANGEMENTS FOR MEASURING TWO OR MORE VARIABLES NOT COVERED IN A SINGLE OTHER SUBCLASS; TARIFF METERING APPARATUS; MEASURING OR TESTING NOT OTHERWISE PROVIDED FOR
    • G01D1/00—Measuring arrangements giving results other than momentary value of variable, of general application
    • 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

Definitions

  • the present invention relates to a method for detecting an abnormality in a manufacturing apparatus, a manufacturing plant, a product, or the like, and a system used for implementing the method.
  • a main object of the present invention is to provide an abnormality detection method capable of detecting various abnormalities with high accuracy by using a time-series data and a prediction model thereof, and a system used for implementing the method. is there. Disclosure of the invention
  • the present inventor has constructed a prediction model based on a plurality of past time-series signals, provided the obtained time-series signals to the prediction model, predicted the trends of the time-series signals, and obtained the prediction results.
  • linear autoregression for example, GEP ⁇ . ⁇ and GM Jenkins, Time Series Analysis: Forecasting and Control Holden-Day Inc., 1976
  • neural network for example, No. LA-UR 87-2662 (1987). No. LA-U 88-418 (1988)
  • Dynamic Replacement Sheet Rule 26
  • Radial basis function network for example, M. Casdagli, Physica D, Vol. 35, pp. 335-356 (1989)
  • simplex projection for example, Sugihara and RM May, Nature, Vol. 344) , pp.734-741 (1990)
  • An abnormality detection method is a method for measuring a physical quantity that changes with time emitted by an object at predetermined time intervals, and detecting an abnormality of the object based on obtained time-series data.
  • a plurality of first vectors of appropriate dimensions are generated from the time-series data, a second vector is generated by translating each first vector by an appropriate time, and the generated second vector and first vector are generated. Calculate the deviation from the torque, compare the calculated value with the set threshold, and judge whether or not there is an abnormality based on the comparison result.
  • the method of calculating the deviation includes a step of selecting an appropriate vector from the first vector, and a step of selecting a predetermined number of vectors from the first vector in an ascending order of distance to the selected vector.
  • a step of extracting the vector from the vector a step of generating a translation vector obtained by translating the extraction vector and the selection vector by an appropriate time, and a step of extracting the generated translation vector, the selection vector, and the extraction.
  • the method includes a step of obtaining a difference vector from the vector and a step of calculating a variance of the obtained difference vector.
  • N is the number of data.
  • An appropriate embedding dimension D and sampling interval ⁇ t are set, and a vector X (t) expressed by the following equation (2) is generated from the time-series data. It has been proven that generating vectors X (t) for different embedding dimensions D is equivalent to obtaining time-series data for multiple types of physical quantities, respectively. Thus, abnormality detection can be performed with high accuracy.
  • the vector X (t) is translated by an appropriate time interval ⁇ to generate a vector ⁇ (t) expressed by the following equation (3).
  • Y (t) ⁇ X (t + ⁇ ⁇ t), x (t + TA t-At), ⁇ , ⁇ (t + ⁇ ⁇ t-(D- 1) ⁇ t) ⁇ ... (3 )
  • the N vectors X (t) force and the appropriate vector X (t (0)) are selected, and the selected vector X (t (0)) and the other vectors (N-1) are selected.
  • V (t (k)) X (t (k)) ⁇ (t (k))... (5)
  • FIG. 1 is a block diagram showing a configuration of an abnormality detection system according to the present invention
  • FIG. 2 is a flowchart showing an abnormality detection procedure using the computer of FIG. 1
  • FIG. 3 is a flowchart of FIG. Fig. 4 is a flowchart showing the abnormality detection procedure by the computer shown in Fig. 4
  • Fig. 4 is a flowchart showing the abnormality detection procedure by the computer shown in Fig. 1
  • Fig. 5 is the flowchart shown in Fig. 1.
  • 6 (A) and 6 (B) are graphs showing the time series data used in the comparative test
  • FIG. 7 is a flowchart showing the abnormality detection procedure according to the method of the present invention.
  • FIGS. 9 (A) and 9 (B) are graphs showing the time series data used in the comparative test
  • FIG. 10 is shown in FIGS. 9 (A) and 9 (B) by the method of the present invention.
  • Fig. 11 shows the results from the two time series data shown in Figs. 6 is a graph showing a result of creating a power spectrum.
  • FIG. 1 is a block diagram showing the configuration of the abnormality detection system according to the present invention.
  • reference numeral 1 denotes a sound sensor, an ultrasonic sensor, a vibration sensor, an optical sensor, a voltmeter, an ammeter, and the like.
  • Measuring instrument The type of measuring instrument 1 is selected according to the physical quantity to be measured, and measuring instrument 1 measures the physical quantity of the object to be detected for abnormality at a predetermined cycle.
  • the output of the measuring instrument 1 is supplied to an analog / digital (AZD) converter 2, where it is converted into a digital signal, which is then supplied to a memory 32 provided in a computer 3, where the time is output. Stored as sequence data.
  • APD analog / digital
  • the time series data stored in the memory 32 within a predetermined time is obtained by the central processing unit 31 by the vector deviation calculating unit 33, and the vector deviation calculating unit 33 outputs E (trans 'Calculate and give it to the failure diagnosis unit 34 that determines whether it is abnormal.
  • the failure condition input unit 35 has a first threshold and a second threshold set, and the failure diagnosis unit 34 has a vector deviation.
  • the E (trans) provided by the calculation unit 33 is compared with the first threshold or the second threshold set in the failure condition input unit 35, respectively, and E (trans) is the first ⁇ replacement sheet (Rule 26) If the value is equal to or more than the value, or if E (trans) is equal to or less than the second threshold value, it is determined that an abnormality has occurred, and the output unit 36 outputs an alarm to the external device 4 such as a system control device or an alarm device. .
  • FIGS. 2 to 5 are flowcharts showing the abnormality detection procedure using the computer 3 shown in FIG.
  • the vector deviation calculation part 33 of the computer 3 has an embedding dimension D, a sampling interval A t, a translation time TA t, a vector selection number M, and a repetition number K, Q.
  • a first threshold value and a second threshold value are set in the input unit 35 (step S 1).
  • the memory 32 of the computer 3 stores the digital signal given from the AZD converter 2 as time-series data (step S2), and stores the time-series data stored in the memory 32 within a predetermined time. Is supplied to the vector deviation calculation unit 33 by the central processing unit 31.
  • the vector deviation calculation unit 33 outputs the given time series data ⁇ X (t) ⁇ ,! From (N is the number of data), a vector X (t) represented by the following equation (2) is generated (step S3).
  • the vector deviation calculator 33 translates the vector X (t) by an appropriate time interval ⁇ , and 3) Generate a vector ⁇ (t) represented by the equation (step S4).
  • the vector deviation calculating unit 33 generates an arbitrary uniform random number (step S5), and uses the generated uniform random number from the N vectors X (t) generated in step S3.
  • M (M ⁇ N) vectors X (t (0)) are selected (step S6), and one vector X (t (0)) replacement sheet (rule 26) is appropriately selected from them. ) Is identified (step S7).
  • h x ⁇ [x (t (0)-d ⁇ t) — x (t (k)-d ⁇ t)] 2
  • the vector deviation calculator 33 calculates (K + 1) number of difference vectors V (t (k)) based on the following equation (5) (step S11), and calculates each difference vector V (t (k)). Based on t (k)) and the following equation (7), a first-order translation error E (trans)! is obtained (step S12).
  • V (t (k)) X (t (k)) -Y (t (k))-(5)
  • the vector deviation calculator 33 repeats the operations from step S7 to step S12 until it is determined in step S13 that the M primary translation errors E (trans), have been obtained. For the M vectors X (t (0)) selected in, obtain the first-order translation error E (trans) ,. Beq replacement paper (Rule 26) The torque deviation calculator 33 determines that M] translation errors E (trans) 1 have been obtained, calculates an intermediate value (median) thereof, and obtains a first-order average translation error E (trans). Is obtained (step S14).
  • the vector deviation calculation unit 33 performs the operations from step S5 to step S14 until it is determined in step S15 that Q primary average translation errors E (trans),, are obtained. Repetition, Q first-order average translation errors
  • the vector deviation calculation unit 33 gives it to the failure diagnosis unit 34.
  • the failure diagnosis unit 34 compares E (trans) given from the vector deviation calculation unit 33 with the first threshold and Z or the second threshold set by the failure condition input unit 35 (step S17). If E (trans) is equal to or greater than the first threshold, or if E (trans) is equal to or less than the second threshold, it is determined that there is an abnormality (step S18), and the external device 4 To output an alarm (Step S19) o
  • FIGS. 6 and 9 are graphs showing the time series data used in the comparative test, in which the vertical axis represents the signal intensity and the horizontal axis represents the time from the start of the measurement.
  • Fig. 6 is a graph showing the vibration of the motor rotating part of the fan in a time series.
  • Fig. 6 (A) shows the time series data when no abnormality occurs.
  • Figure (B) shows time-series data when an error occurs.
  • the signal strength on the vertical axis is converted into an electric signal by the vibration sensor, but the reflected wave of the vibration part may be detected by the optical sensor.
  • the voltage characteristics of the motor measured by a voltmeter or motor switching ( ⁇ ⁇ Article 26) 1 O
  • a similar oscillating component can be obtained from the current value of the evening current measured by an ammeter.
  • Fig. 9 is a graph showing the sound signal as time-series data.
  • Fig. 9 (A) shows the time-series data when no abnormality occurs
  • Fig. 9 (B) shows the time when abnormality occurs. This is the time series data at the time of performing.
  • the signal strength on the vertical axis is converted into an electric signal by a sound wave sensor.
  • FIG. 7 is a graph showing the result of calculating E (trans) from the two time-series data shown in FIGS. 6 (A) and (B) by the method of the present invention.
  • the symbols ⁇ indicate the abnormal cases. Note that the sampling interval ⁇ t is 1, the translation time ⁇ t T is 5, the vector selection number M is 300, the vector selection number K is 4, and the number of repetitions Q is The embedding dimension D is set to 5, 6, 6.7, 8, 9, and 10 in 20.
  • Fig. 8 is a graph showing the result of creating a power spectrum from the two time-series data shown in Figs. 6 (A) and (B) by the conventional method. Indicates a case where there is an abnormality.
  • the abnormal case and the non-abnormal case are clearly separated, and the abnormal condition is set by setting the threshold value to 0.05. Can be determined with high accuracy.
  • the power spectrum in the case of an abnormality and the power spectrum in the case of a non-abnormal case often overlapped, and a judgment error occurred between the two.
  • FIG. 10 is a graph showing the result of calculating E (trans) from the two time series data shown in FIGS. 9 (A) and (B) according to the present invention.
  • the sampling interval ⁇ t is set to 1 and the translation time T ⁇ t is replaced (Rule 26).
  • T to 5 vector selection ⁇ to 300, vector selection ⁇ to 4, repetition Q to 20 and embedding dimension D to 1 1, 12, 1 3, 1 4, 1 It is set to 5.
  • Fig. 11 is a graph showing the result of creating a power spectrum from the two time-series data shown in Figs. 9 ( ⁇ ) and ( ⁇ ) by the conventional method. The broken line shows the case where there is an abnormality.
  • the abnormal case and the non-abnormal case are clearly separated, and by setting 0.90 as the threshold value, the abnormal condition is determined. It can be determined with high accuracy whether or not there is.
  • FIG. 11 in the conventional method, there are many places where the power spectrums of the abnormal case and the non-abnormal case overlap, and a judgment error occurs between the two.

Landscapes

  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Testing Or Calibration Of Command Recording Devices (AREA)
  • Testing And Monitoring For Control Systems (AREA)
  • Testing Of Devices, Machine Parts, Or Other Structures Thereof (AREA)

Description

明 糸田
異常検出方法及び異常検出システム 技術分野
本発明は製造装置， 製造プラン ト又は製品等の異常を検出する方 法及びその実施に使用するシステムに関する。 背景技術
機械装置が発生する音響又は振動、 化学ブラン トにおける反応生 成の濃度、 電子デバィスに印加した電流又は電圧等を経時的に測定 して、 機械装置の異常、 反応制御の異常、 又は製品品質の異常等を 検出するこ とが行われている。
そのような異常検出方法として、 丹慶， 奥村， 佐藤， 小林 訳，
「ニュ ー メ リ カル レ シ ピ ' イ ン ' シー J , 第 1 2章一第 1 4章 （技 術評論社 1993年） には、 所定の周期で測定した時系列信号のスぺ ク トル解析によって異常を検出する方法が提案されている。
これは、 所定の周期で測定した時系列信号を、 適宜の周波数領域 において適当な周波数間隔で高速フーリエ変換し、 次の （ 1 ) 式を 満足する係数 a k , b k を求めて a k 2 + b k 2 を算出する。 この とき、 フーリエ変換の開始周波数及び終了周波数において、 （ 1 ) 式で示した f ( X ) の急激な立ち上がり又は立ち下がりを防止する ために、 適切なデータウイ ン ドウ， 例えば H ann i ngウイ ン ドウを設 定しておき、 それを f ( X ) に乗じて a k 2 + b k 2 の算出誤差を 低減している。 f ( t ) = ( a k sin wk t + b k cosoj k t ) ··· f l ) 横軸が 1 o g ωであり、 縦軸が l o g (スぺク トル密度） である 差替え用紙 （規則 26) 座摞領域に、 算出した a k 2 + b k 2 をプロ ッ 卜 してパワースぺク トルを作成し、 周期的な変動成分を除いた後、 最小自乗法によ り適 合する直線を求める。 そ して、 求めた直線のレベル又は傾き （パヮ 一スぺク トル指数） と予め定めた閾値との比較結果よ り、 異常の有 無を判断する。
しかしながら、 前述した如きスぺク トル解析によって異常を検出 する方法にあっては、 パワースぺク トルに適合する直線を高精度に 求めるこ とが困難であり、 時系列信号に含まれる測定ノ イズによつ て誤判定が発生する という問題があった。 また、 適切なデータウ イ ン ドウを設定するこ とが困難であり、 パワースぺク トルの誤差が大 きいという問題もあった。
本発明の主な目的は、 時系列デ一夕 とその予測モデルとを用いる こ とによって、 種々 の異常を高精度に検出し得る異常検出方法及び その実施に使用するシステムを提供するこ とにある。 発明の開示
本発明者は、 過去の複数の時系列信号に基づいて予測モデルを構 築しておき、 得られた時系列信号を予測モデルに与えて、 時系列信 号の動向を予測し、 その予測結果と実際の時系列デー夕 とを比較す るこ とによって、 何らかの異常検出が可能ではないかという着想に 至った。
そこで、 本発明者は、 （ 1 ) 線形自己回帰 （例えば、 G. E. P. Βο.χ and G. M. Jenkins, Time Series Analysis ： Forecasting and Control Holden-Day Inc. , 1976)， （ 2 ) ニューラルネ ッ ト ワーク （例えば、 A. S.し apedes and R. Farber, Los Alamos National Laboratory Report No. LA-UR 87-2662 ( 1987). No. LA-U 88-418 (1988) ) ， （ 3 ) 動 差替え用紙 （規則 26) 径基底関数ネ ッ トワーク （例えば、 M. Casdagli, Physica D, Vol.35, pp.335-356 (1989) ) , ( 4 ) シンプレ ッ クス投影法 （例えば、 Sugihara and R.M.May, Nature, Vol.344, pp.734-741 (1990)) 等 の数多く の予測モデルの構築を試みたが、 何れも不完全なものであ つ Ί乙。
試行錯誤の結果、 Wayland らによって提案された時系列解析アル ゴリ ズム (R. Wayland, D. Bromley, D. Pickett and A. Passamante. Physical Review Letters, Vol.70, pp.580 - 582 (1993)) に着目 した。
本発明に係る異常検出方法は、 被検出物が発する時間と共に変化 する物理量を所定の時間間隔で測定し、 得られた時系列データに基 づいて前記被検出物の異常を検出する方法において、 時系列データ から適当な次元の複数の第 1 ベク トルを生成し、 各第 1 ベク トルを 適宜時間だけ並進させた第 2べク トルを生成し、 生成した第 2べク トルと第 1 ベク トルとの偏差を演算し、 その算出値と設定した閾値 とを比較し、 その比較結果から異常であるか否かを判断する。
このような方法において、 偏差を演算する手法は、 第 1 ベク トル から適宜のべク トルを選択するステップと、 この選択べク トルに距 離が近い順に所定数のべク トルを第 1 べク トルから抽出するステツ プと、 この抽出べク トル及び選択べク トルを適宜時間だけ並進させ た並進べク トルをそれぞれ生成するステップと、 生成した並進べク トルと選択べク トル及び抽出べク トルとの差分べク トルをそれぞれ 求めるステップと、 求めた差分べク トルの分散を算出するステッ プ とを含んでいる。
以下、 本発明の異常検出方法の概要について説明する。 まず、 モ 一夕の回転音， 電子デバイスに印加した電流又は電圧等， 被検出物 差替え用紙 （規則 26) の物理量を適宜の周期で測定して、 時系列デ一夕 { X ( t ) } N
(Nはデータ数） を得る。 適宜の埋め込み次元 D及びサンプリ ング 間隔 Δ tを設定して、 次の （ 2 ) 式で表されるべク トル X ( t ) を 前記時系列データから生成する。 異なる埋め込み次元 Dに対してそ れぞれベク トル X ( t ) を生成する と、 複数種類の物理量に対する 時系列データをそれぞれ得たこ と と等価になるこ とが証明されてお り、 これによつて異常の検出を高精度に行う こ とができる。 このべ ク トル X ( t ) を適宜の時間間隔 ΤΔ ίだけ並進させて、 次の （ 3 ) 式で表されるべク トル Υ ( t ) を生成する。
X ( t ) = {x ( t ) , x ( t - A t ) , ·'·,
X ( t - （D - 1 ) 厶 t ) } … （ 2 )
Y ( t ) = { X ( t + Τ Δ t ) , x ( t + TA t - A t ) , ···, χ ( t + Τ Δ t - (D- 1 ) Δ t ) } … （ 3 )
N個のベク トル X ( t ) 力、ら適宜のベク トル X ( t ( 0 ) ) を選 択し、 選択したベク トル X ( t ( 0 ) ) とそれを除く （N— 1 ) 個 の他の全てのベク トル X ( t ) との間の距離 h x を次の （ 4 ) 式に 基づいて算出し、 距離が短いものから順に K個 （K < N) 抽出する。 h = Λ/Σ [ χ ( t (0) 一 d厶 t ) - x ( t (k)— d厶 t ) ] 2
… ( ) こ う して得られた （K+ l ) 個のベク トル X ( t ( k ) ) C k = 0 , 1 , ···， K) に対応するベク トル Y ( t ( k) ) を算出 し、 次 の （ 5 ) 式に基づいて差分べク トル V ( t ( k) ) を算出する。 そ して、 得られた （K + 1 ) 個の差分ベク トル V ( t ( k ) ) の分散 である並進誤差 E (trans)を次の （ 6 ) 式に基づいて求める こ とに よって、 時系列データの乱雑さの程度を E (trans)と して定量化す 差替え用紙 （規則 26) る
V ( t (k) ) =X ( t (k) ) Υ ( t ( k ) ) … （ 5 )
1 κ V ( t (k) ) - < V > I 2
E (trans) = ∑
K + 1 κ-° I < V > i 2
… ( 6 )
1 K
但し、 < V > ∑ V ( t (k) )
K + 1 k-°
本発明者が鋭意研究した結果、 時系列データに含まれるノ イズの 比率が増加する異常が発生した場合に、 E ( trans)の値か大き く な り、 時系列データに含まれるノイズのパワースペク トル指数が増加 する異常が発生した場合に、 E (trans)の値が小さ く なると、 パヮ 一スぺク トル指数より E ( trans)の方が鋭敏に変化するので認識し やすいという知見を得た。 そこで、 前者の場合の異常発生を判断す るための閾値と、 後者の場合の異常発生を判断するための閾値とを 設定し、 E (trans)が両閾値の間の領域を外れた場合、 異常である と判断する方法とそのシステムの発明に至った。 図面の簡単な説明
第 1図は本発明に係る異常検出システムの構成を示すプロ ッ ク図、 第 2図は第 1図のコ ン ピュータによる異常検出手順を示すフローチ ヤ ー ト、 第 3図は第 1図に示したコンピュータによる異常検出手順 を示すフローチヤ一ト、 第 4図は第 1 図に示したコンピュータによ る異常検出手順を示すフ ローチャー ト、 第 5図は第 1 図に示したコ ンビュー夕による異常検出手順を示すフローチャー ト、 第 6図 （A) (B) は比較試験に用いた時系列デ一夕を示すグラフ、 第 7図は本 発明方法によって第 6図 （A) ， （B) に示した 2つの時系列デー 差替え用紙 （規則 26) 夕から E (trans)を算出した結果を示すグラ フ、 第 8図は従来の方 法によって第 6図 （A) , ( B ) に示した 2つの時系列データから パワースぺク トルを作成した結果を示すグラフ、 第 9図 （A) ， ( B ) は比較試験に用いた時系列データを示すグラ フ、 第 1 0図は本 発明方法によって第 9図 （A) ， ( B ) に示した 2つの時系列デー 夕から E (trans)を算出した結果を示すグラ フ、 第 1 1 図は従来の 方法によって第 9図 （A) , ( B ) に示した 2つの時系列データか らパワースぺク トルを作成した結果を示すグラフである。 発明を実施するための最良の形態
以下、 本発明の実施の形態を図面に基づいて具体的に説明する。 第 1 図は本発明に係る異常検出システムの構成を示すブロ ッ ク図 であり、 図中 1 は音波セ ンサ， 超音波セ ンサ， 振動セ ンサ， 光セ ン サ， 電圧計又は電流計等の測定器である。 測定器 1 の種類は測定す べき物理量に応じて選択してあり、 測定器 1 は異常検出の対象であ る被検出物の物理量を所定の周期で測定する。 測定器 1 の出力はァ ナログ//ディ ジタル （ AZD ) 変換器 2に与えられ、 そこでデイ ジ 夕ル信号に変換され、 コ ンピュータ 3に備えられたメ モ リ 32に与え られ、 そこに時系列データと して記憶される。
所定の時間内にメ モ リ 32に記憶された時系列データは中央処理装 置 31によってべク トル偏差計算部 33> えられ、 べク トル偏差計算 部 33は後述するようにして E ( trans ' 算出し、 それを異常である か否かを判断する故障診断部 34に与える。 故障条件入力部 35には第 1 閾値及び第 2閾値が設定してあり、 故障診断部 34はベク トル偏差 計算部 33から与えられた E ( trans)と故障条件入力部 35で設定され た第 1 閾値又は第 2閾値とをそれぞれ比較し、 E (trans)が第 1 閎 差替え用紙 （規則 26) 値以上であった場合、 又は E (trans)が第 2閾値以下であった場合 は異常である と判断し、 出力部 36からシステム制御装置又は警報装 置等の外部装置 4へ警報を出力させる。
第 2図〜第 5図は第 1 図に示したコ ン ピュータ 3による異常検出 手順を示すフローチャー トである。 コ ン ピュータ 3のべク トル偏差 計算部 33には埋め込み次元 D、 サンプリ ング間隔 A t 、 並進時間 T A t の丁、 ベク トル選択数 M、 繰り返し回数 K, Qが設定してあり、 故障条件入力部 35には第 1 閾値及び第 2閾値が設定してある （ステ ップ S 1 ) 。 コ ンピュータ 3のメ モ リ 32は A ZD変換器 2から与え られたディ ジタル信号を時系列データ と して記憶し （ステップ S 2 ) 所定の時間内にメモ リ 32に記憶された時系列データは中央処理装置 31によってべク トル偏差計算部 33に与えられる。
べク トル偏差計算部 33は与えられた時系列デ一夕 { X ( t ) } ,一！ (Nはデータ数） から次の （ 2 ) 式で表されられるべク トル X ( t ) を生成する （ステップ S 3 ) 。
X ( t ) = { X ( t ) , x ( t - Δ t ) , ···,
χ ( t - （D - 1 ) Δ t ) } - ( 2 ) 更に、 べク トル偏差計算部 33はべク トル X ( t ) を適宜の時間間 隔 Τ Δ ΐだけ並進させ、 次の （ 3 ) 式で表されるベク トル Υ ( t ) を生成する （ステップ S 4 ) 。
Y ( t ) 二 { X ( t + Τ Δ t ) , χ ( t + Τ Δ t - Δ t ) , ·■·, χ ( t + ΤΔ ί - (D - l ) A t ) } ." ( 3 ) べク トル偏差計算部 33は任意の一様乱数を発生し （ステップ S 5 ) 発生した一様乱数を用いて、 ステップ S 3で生成した N個のべク ト ル X ( t ) から M個 （Mく N) のべク トル X ( t ( 0 ) ) を選択し (ステ ッ プ S 6 ) 、 その中から適宜に 1 つのベク トル X ( t ( 0 ) ) 差替え用紙 （規則 26) を特定する （ステップ S 7 ) 。 べク トル偏差計算部 33は、 特定した ベク トル X ( t ( 0 ) ) とそれを除く （N— 1 ) 個の他の全てのベ ク トル X ( t ) との間の距離 h x を次の （ 4 ) 式に基づいて算出し (ステップ S 8 ) 、 算出した距離が短いものから順に並べ換えて、 その先頭から K個 （Kく M) のべク トルを抽出 して、 べク トル X ( t (k) ) ( k = 1 , ·.·, Κ) を得る （ステップ S 9 ) 。 h x = ∑ [x ( t (0) - d Δ t ) — x ( t (k) - d Δ t ) ] 2
'V d - o
… ( 4 ) べク トル偏差計算部 33は、 ステップ S 7 , S 9で得たベク トル X ( t (k) ) ( k = 0 , 1 , ···, Κ) を時間間隔 ΤΔ tだけ並進さ せたべク トル Y ( t (k) ) を生成する （ステップ S 10) 。
べク トル偏差計算部 33は、 次の （ 5 ) 式に基づいて （K + 1 ) 個 の差分ベク トル V ( t (k) ) を算出し （ステップ S 11) 、 各差分 ベク トル V ( t (k) ) 及び次の （ 7 ) 式に基づいて、 1次並進誤 差 E (trans)! を求める （ステッ プ S 12) 。
V ( t ( k) ) =X ( t (k) ) -Y ( t ( k) ) - ( 5 )
1 κ I V ( t (k) )- < V > | 2
E (trans) i = ；
K + 1 k"° I < V > I 2
… ( 7 )
1 K
但し、 く V > = ∑ V ( t (k) )
K+ 1 k"°
べク トル偏差計算部 33は、 ステップ S 13で M個の 1 次並進誤差 E (trans), が得られたと判断されるまで、 ステップ S 7〜ステッ プ S 12までの操作を繰り返し、 ステップ 6で選択した M個のべク トル X ( t ( 0 ) ) について、 1次並進誤差 E ( trans) , を得る。 べク 差替え用紙 （規則 26) トル偏差計算部 33は、 M個の 】 次並進誤差 E ( trans) 1 が得られた と判断すると、 それらの中間値 （メディ アン） を算出して、 1 次平 均並進誤差 E ( trans) を得る （ステップ S 14) 。
そして、 べク トル偏差計算部 33は、 ステップ S 15で Q個の 1 次平 均並進誤差 E (trans), ,が得られたと判断されるまで、 ステップ S 5〜ステップ S 14までの操作を繰り返し、 Q個の 1 次平均並進誤差
E (trans)! ,が得られたと判断する と、 それらの平均値を算出 して E (trans)を得る （ステップ S 16) 。 これによつて、 時系列データ の E ( trans)が高精度に求められる。
E (trans)が得られると、 ベク トル偏差計算部 33はそれを故障診 断部 34に与える。 故障診断部 34は、 べク トル偏差計算部 33から与え られた E (trans)と故障条件入力部 35で設定された第 1 閾値及び Z 又は第 2閾値とをそれぞれ比較し （ステップ S 17) 、 E (trans)が 第 1 閾値以上であった場合、 又は E (trans)が第 2閾値以下であつ た場合は異常であると判断し （ステッ プ S 18) 、 出力部 36から外部 装置 4へ警報を出力させる （ステップ S 19) o
次に、 本発明方法による異常検出と従来の方法による異常検出と を比較した結果について説明する。
第 6図及び第 9図は比較試験に用いた時系列データを示すグラフ であり、 縦軸は信号強度を、 横軸は測定開始からの時間をそれぞれ 示している。 第 6図は、 フ ァ ンのモータ回転部の振動を時系列デ一 夕にして示したグラフであり、 第 6図 （A) は異常が発生していな いときの時系列データ、 第 6図 （ B ) は異常が発生したときの時系 列データである。 縦軸の信号強度は、 振動センサによって電気信号 に変換したものであるが、 光センサで振動部の反射波を検出 しても 良い。 また、 モータの電圧特性を電圧計で測定した電圧値又はモー 替换页 (细则第 26条） 1 O
夕の電流特性を電流計で測定した電流値からも同様な振動成分が得 られる。
また、 第 9図は音響信号を時系列データにして示したグラフであ り、 第 9図 （A) は異常が発生していないときの時系列データ、 第 9図 （ B ) は異常が発生したときの時系列データである。 縦軸の信 号強度は、 音波センサによって電気信号に変換したものである。
第 7図は本発明方法によって第 6図 （A) ， （ B) に示した 2つ の時系列データから E (trans)を算出した結果を示すグラフであり、 暴印は異常でない場合を、 〇印は異常である場合をそれぞれ示して いる。 なお、 サンプリ ング間隔 Δ t を 1 に、 並進時間 Τ Δ t の Tを 5に、 べク トル選択数 Mを 3 0 0に、 べク トル選択数 Kを 4 に、 繰 り返し回数 Qを 2 0に、 埋め込み次元 Dを 5 , 6. 7 , 8， 9 , 10 にしてある。 また、 第 8図は従来の方法によって第 6図 （A) , ( B ) に示した 2つの時系列データからパワースペク トルを作成した結果 を示すグラフであり、 実線は異常でない場合を、 破線は異常である 場合をそれぞれ示している。
第 7図から明らかな如く、 本発明方法にあっては、 異常である場 合と異常でない場合とが明瞭に分かれており、 閾値として 0. 0 5 を設定しておく ことによって、 異常であるか否かを高精度に判断す ることができる。 一方、 第 8図から明らかな如く、 従来の方法にあ つては、 異常である場合と異常でない場合とのパワースぺク トルか 重なっている箇所が多く、 両者の判断誤差が生じた。
また、 第 1 0図は本発明によって第 9図 （ A ) , ( B ) に示した 2つの時系列デ一夕から E (trans)を算出した結果を示すグラフで あり、 ·印は異常でない場合を、 〇印は異常である場合をそれぞれ 示している。 なお、 サンプリ ング間隔 Δ t を 1 に、 並進時間 T△ t 差替え用紙 （規則 26) の Tを 5に、 ベク トル選択数 Μを 3 0 0 に、 ベク トル選択数 Κを 4 に、 繰り返し回数 Qを 2 0 に、 埋め込み次元 Dを 1 1 , 12, 1 3 , 1 4, 1 5にしてある。 また、 第 1 1 図は従来の方法によって第 9図 （Α ) ( Β ) に示した 2つの時系列データからパワースぺク トルを作成し た結果を示すグラフであり、 実線は異常でない場合を、 破線は異常 である場合をそれぞれ示している。
第 1 0図から明らかな如く、 本発明方法にあっては、 異常である 場合と異常でない場合とが明瞭に分かれており、 閾値として 0 . 9 0を設定しておく こ とによって、 異常であるか否かを高精度に判断 することができる。 一方、 第 1 1 図から明らかな如く、 従来の方法 にあっては、 異常である場合と異常でない場合とのパワースぺク ト ルが重なっている箇所が多く、 両者の判断誤差が生じた。 産業上の利用可能性
以上詳述した如く、 本発明にあっては、 測定対象の物理量の変動 が非線形ダイナミ クスによって生じる場合においても、 種々 の異常 を高精度に検出し得るため、 製造プラ ン トの運転管理、 製造装置の 部品交換、 品質管理等を的確に行う ことができる等、 優れた効果を 奏する。
差替え用紙 （規貝 IJ26)

Claims

請 求 の 範 囲
1 . 被検出物が発する時間と共に変化する物理量を所定の時間間 隔で測定し、 得られた時系列データに基づいて前記被検出物の異常 を検出する方法において、 前記時系列データから適当な次元の複数 の第 1 べク トルを生成するステッ プと、 各第 1 べク トルを適宜時間 だけ並進させた第 2べク ト ルを生成するステップと、 生成した第 2 べク トルと第 1 べク トルとの偏差を演算するステッ プと、 その演算 値と設定した閾値とを比較するステップと、 その比較結果から前記 被検出物が異常であるか否かを判断するステップとを有する こ とを 特徴とする異常検出方法。
2 . 前記偏差を演算するステップは、 第 1 ベク トルから適宜のベ ク ト ルを選択する第 1 ステップと、 この選択べク ト ルに距離が近い 順に所定数のべク トルを第 1 べク トルから抽出する第 2 ステップと 、 この抽出べク トル及び前記選択べク トルを適宜時間だけ並進させ た並進ベク トルをそれぞれ生成する第 3 ステップと、 生成した並進 べク トルと前記選択べク トル及び抽出べク トルとの差分べク トルを それぞれ求める第 4 ステップと、 求めた差分べク トルの分散を算出 する第 5 ステップとを含むこ とを特徴とする請求項 1 記載の異常検 出方法。
3 . 前記選択ベク トルを変えて前記第 1 〜第 5 ステッ プを繰り返 すこ とにより複数の分散を算出し、 算出 した複数の分散の中間値又 は平均値を求めるこ とを複数回繰り返して、 複数の中間値又は平均 値を求め、 求めた複数の中間値又は平均値の平均値を前記閾値と比 較する演算値と設定する こ とを特徴とする請求項 2記載の異常検出 方法。
4 . 前記選択べク トルを変えて前記第 1 〜第 5 ステ ッ プを繰り返 差替え用紙 （規則 26) すこ とによ り複数の分散を算出し、 算出した複数の分散の中間値又 は平均値を求め、 求めた中間値又は平均値を前記閾値と比較する演 算値と設定するこ とを特徴とする請求項 2記載の異常検出方法。
5 . 被検出物が発する時間と共に変化する物理量を所定の時間間 隔で測定し、 得られた時系列データに基づいて前記被検出物の異常 を検出するシステムにおいて、 前記時系列データから適当な次元の 複数の第 1 ベク トルを生成する手段と、 各第 1 ベク トルを適宜時間 だけ並進させた第 2ベク トルを生成する手段と、 生成した第 2 べク トル と第 1 ベク トルとの偏差を演算する演算手段と、 その演算値と 設定した閾値とを比較する手段と、 その比較結果から前記被検出物 が異常であるか否かを判断する手段とを備えるこ とを特徴とする異 常検出システム。
6 . 前記演算手段は、 第 1 ベク トルから適宜のベク トルを選択す る手段と、 この選択べク トルに距離が近い順に所定数のべク トルを 第 1 ベク トルから抽出する手段と、 この抽出ベク トル及び前記選択 べク ト ルを適宜時間だけ並進させた並進べク ト ルをそれぞれ生成す る手段と、 生成した並進べク トルと前記選択べク トル及び抽出べク トルとの差分べク トルをそれぞれ求める手段と、 求めた差分べク ト ルの分散を算出する手段とを含むこ とを特徴とする請求項 5 記載の 異常検出システム。
7 . 被検出物の異常を検出するシステムにおいて、 前記被検出物 が発する時間と共に変化する物理量を所定の時間間隔で測定する測 定手段と、 該则定手段にて得られた測定値の時系列データを記億す る手段と、 記憶した時系列データから適当な次元の複数の第 1 べク トルを生成する手段と、 各第 1 ベク トルを適宜時間だけ並進させた 第 2べク トルを生成する手段と、 生成した第 2べク トル と第 1 べク 差替え用紙 （規則 26) トルとの偏差を演算する手段と、 その演算値と設定した閾値とを比 較する手段と、 その比較結果に応じて前記被検出物が異常であるこ とを検出する手段とを備えるこ とを特徴とする異常検出システム。
8 . 前記測定手段は、 音波センサ， 超音波センサ， 振動センサ， 光センサ， 電圧計及び電流計からなる群から選ばれた機器であるこ とを特徴とする請求項 7記載の異常検出システム。
差替え用紙 （規則 26)
PCT/JP1997/000367 1996-02-16 1997-02-13 Procede de detection d'anomalies et systeme de detection d'anomalies WO1997030329A1 (fr)

Priority Applications (2)

Application Number Priority Date Filing Date Title
AU16716/97A AU1671697A (en) 1996-02-16 1997-02-13 Method of detecting abnormality and abnormality detection system
EP97902662A EP0881467A4 (en) 1996-02-16 1997-02-13 METHOD AND SYSTEM FOR DETECTING ANOMALIES

Applications Claiming Priority (2)

Application Number Priority Date Filing Date Title
JP2973396 1996-02-16
JP8/29733 1996-02-16

Publications (1)

Publication Number Publication Date
WO1997030329A1 true WO1997030329A1 (fr) 1997-08-21

Family

ID=12284315

Family Applications (1)

Application Number Title Priority Date Filing Date
PCT/JP1997/000367 WO1997030329A1 (fr) 1996-02-16 1997-02-13 Procede de detection d'anomalies et systeme de detection d'anomalies

Country Status (6)

Country Link
EP (1) EP0881467A4 (ja)
KR (1) KR19990082532A (ja)
CN (1) CN1213429A (ja)
AU (1) AU1671697A (ja)
TW (1) TW327221B (ja)
WO (1) WO1997030329A1 (ja)

Cited By (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN114184274A (zh) * 2021-11-19 2022-03-15 哈尔滨工程大学 一种受高频噪声干扰振动信号的非线性特性诊断方法

Families Citing this family (15)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US6526405B1 (en) * 1999-12-17 2003-02-25 Microsoft Corporation Determining similarity between event types in sequences
US7152068B2 (en) * 2001-12-21 2006-12-19 Honeywell International Inc. Method and apparatus for retrieving time series data related to an activity
JP4200332B2 (ja) * 2006-08-29 2008-12-24 パナソニック電工株式会社 異常監視装置、異常監視方法
KR101029394B1 (ko) 2006-09-28 2011-04-15 미쓰비시덴키 가부시키가이샤 고장 검출 장치, 고장 검출 방법 및 고장 검출 프로그램
TWI474023B (zh) * 2008-12-10 2015-02-21 Ind Tech Res Inst 馬達故障診斷方法及其診斷裝置
DE102009028325A1 (de) * 2009-08-07 2011-02-10 Robert Bosch Gmbh Verfahren und Vorrichtung zur Detektion von Schwingungen
KR101659989B1 (ko) * 2014-12-10 2016-09-27 부산대학교 산학협력단 다차원 특징을 이용한 이상신호 분석장치 및 방법
CN105938011B (zh) * 2016-07-15 2017-08-11 武汉理工大学 一种水力除焦过程钻杆工作状态实时判定方法
CN106872168A (zh) * 2017-03-09 2017-06-20 武汉科技大学 一种用于滚动轴承的故障诊断方法
CN107703916B (zh) * 2017-09-30 2019-06-07 英特尔产品（成都）有限公司 电子压力调节器诊断系统及方法
CN109861857A (zh) * 2019-01-28 2019-06-07 网联清算有限公司 故障检测方法及装置
CN109919495A (zh) * 2019-03-12 2019-06-21 中国计量大学 基于阈值计算的光伏电站监测数据异常与程度判断方法
CN111984827A (zh) * 2019-05-24 2020-11-24 上海东方富联科技有限公司 门磁传感器数据异常检测方法及系统、存储介质及终端
KR102230463B1 (ko) * 2019-12-12 2021-03-19 울산대학교 산학협력단 설비 부품의 결함 진단 시스템 및 이의 방법
CN112710446B (zh) * 2020-12-21 2022-06-07 北京和中普方新能源科技有限公司 电动车电池系统振动试验的判断方法、系统及存储介质

Citations (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JPH02226024A (ja) * 1989-02-27 1990-09-07 Toshiba Corp 回転機械のラビング検出装置

Family Cites Families (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US4453407A (en) * 1980-04-17 1984-06-12 Hitachi, Ltd. Vibration diagnosis method and apparatus for rotary machines
JP2909134B2 (ja) * 1990-04-02 1999-06-23 株式会社東芝 回転機械の回転部異常検知装置

Patent Citations (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JPH02226024A (ja) * 1989-02-27 1990-09-07 Toshiba Corp 回転機械のラビング検出装置

Non-Patent Citations (1)

* Cited by examiner, † Cited by third party
Title
See also references of EP0881467A4 *

Cited By (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN114184274A (zh) * 2021-11-19 2022-03-15 哈尔滨工程大学 一种受高频噪声干扰振动信号的非线性特性诊断方法

Also Published As

Publication number Publication date
KR19990082532A (ko) 1999-11-25
AU1671697A (en) 1997-09-02
CN1213429A (zh) 1999-04-07
EP0881467A4 (en) 1999-05-12
EP0881467A1 (en) 1998-12-02
TW327221B (en) 1998-02-21

Similar Documents

Publication Publication Date Title
WO1997030329A1 (fr) Procede de detection d&#39;anomalies et systeme de detection d&#39;anomalies
Tanner Jr et al. Definitions of d′ and η as psychophysical measures
JP3040471B2 (ja) 力学系アナライザ
US11934184B2 (en) Automated analysis of non-stationary machine performance
CN102016607B (zh) 用于分析电力系统的波形信号的方法和装置
KR101025163B1 (ko) 진동 및 소음 전달경로 해석 시스템과 진동 및 소음 전달경로 해석 방법
CN106199267B (zh) 一种电力设备故障特征分析方法
JP6316510B1 (ja) 電動機の診断装置
Balouji et al. Partial discharge classification in power electronics applications using machine learning
JP2022094037A (ja) 異常温度検知装置、異常温度検知方法、及び異常温度検知プログラム
Aggoune et al. Fault detection in the distillation column process using Kullback Leibler divergence
Ferreira et al. An application of machine learning approach to fault detection of a synchronous machine
US8140282B2 (en) Determining a total length for conductive whiskers in computer systems
Roy et al. Software reliability allocation of digital relay for transmission line protection using a combined system hierarchy and fault tree approach
Aggoune et al. Modeling of a distillation column based on NARMAX and Hammerstein models
Wilson et al. Uncertainty Quantification of Capacitor Switching Transient Location Using Machine Learning
Reger et al. Pattern recognition in load profiles of electric drives in manufacturing plants
JP7230371B2 (ja) 異常検出装置、異常検出方法、異常検出プログラム及び異常検出システム
CN104267259B (zh) 电力谐波测量方法和系统
RU2242164C2 (ru) Способ выявления информативных параметров st-сегмента электрокардиосигнала и устройство для его осуществления
CN111679651A (zh) 用于故障引起的变结构变参数系统的辨识方法及系统
JP2002196783A (ja) 時系列信号の識別方法及び装置
EP0386820A1 (en) Method and arrangement for probabilistic pitch measuring
JP2002196784A (ja) 時系列信号の識別方法及び装置
Abid et al. Incipient Broken Bar Detection Using Information Entropy of Second Generation Wavelet Transform

Legal Events

Date Code Title Description
WWE Wipo information: entry into national phase

Ref document number: 97192893.2

Country of ref document: CN

AK Designated states

Kind code of ref document: A1

Designated state(s): AL AM AT AU AZ BA BB BG BR BY CA CH CN CU CZ DE DK EE ES FI GB GE HU IL IS JP KE KG KR KZ LC LK LR LS LT LU LV MD MG MK MN MW MX NO NZ PL PT RO RU SD SE SG SI SK TJ TM TR TT UA UG US UZ VN YU

AL Designated countries for regional patents

Kind code of ref document: A1

Designated state(s): KE LS MW SD SZ UG AT BE CH DE DK ES FI FR GB GR IE IT LU MC NL PT SE BF BJ CF CG CI CM GA GN ML MR NE SN TD

DFPE Request for preliminary examination filed prior to expiration of 19th month from priority date (pct application filed before 20040101)
121 Ep: the epo has been informed by wipo that ep was designated in this application
WWE Wipo information: entry into national phase

Ref document number: 1997902662

Country of ref document: EP

WWE Wipo information: entry into national phase

Ref document number: 1019980706268

Country of ref document: KR

WWP Wipo information: published in national office

Ref document number: 1997902662

Country of ref document: EP

REG Reference to national code

Ref country code: DE

Ref legal event code: 8642

WWP Wipo information: published in national office

Ref document number: 1019980706268

Country of ref document: KR

WWG Wipo information: grant in national office

Ref document number: 1019980706268

Country of ref document: KR

WWW Wipo information: withdrawn in national office

Ref document number: 1997902662

Country of ref document: EP