JP6076751B2 - 異常診断方法およびその装置 - Google Patents

異常診断方法およびその装置 Download PDF

Info

Publication number
JP6076751B2
JP6076751B2 JP2013009316A JP2013009316A JP6076751B2 JP 6076751 B2 JP6076751 B2 JP 6076751B2 JP 2013009316 A JP2013009316 A JP 2013009316A JP 2013009316 A JP2013009316 A JP 2013009316A JP 6076751 B2 JP6076751 B2 JP 6076751B2
Authority
JP
Japan
Prior art keywords
vector
feature
abnormality
sensor
calculating
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
JP2013009316A
Other languages
English (en)
Other versions
JP2014142697A (ja
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.)
Hitachi Ltd
Original Assignee
Hitachi 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 Hitachi Ltd filed Critical Hitachi Ltd
Priority to JP2013009316A priority Critical patent/JP6076751B2/ja
Priority to US14/761,748 priority patent/US9779495B2/en
Priority to PCT/JP2014/050547 priority patent/WO2014115615A1/ja
Publication of JP2014142697A publication Critical patent/JP2014142697A/ja
Application granted granted Critical
Publication of JP6076751B2 publication Critical patent/JP6076751B2/ja
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • G—PHYSICS
    • G05—CONTROLLING; REGULATING
    • G05B—CONTROL OR REGULATING SYSTEMS IN GENERAL; FUNCTIONAL ELEMENTS OF SUCH SYSTEMS; MONITORING OR TESTING ARRANGEMENTS FOR SUCH SYSTEMS OR ELEMENTS
    • G05B23/00—Testing or monitoring of control systems or parts thereof
    • G05B23/02—Electric testing or monitoring
    • G05B23/0205—Electric testing or monitoring by means of a monitoring system capable of detecting and responding to faults
    • G05B23/0218—Electric testing or monitoring by means of a monitoring system capable of detecting and responding to faults characterised by the fault detection method dealing with either existing or incipient faults
    • G05B23/0224—Process history based detection method, e.g. whereby history implies the availability of large amounts of data
    • G05B23/0227—Qualitative history assessment, whereby the type of data acted upon, e.g. waveforms, images or patterns, is not relevant, e.g. rule based assessment; if-then decisions
    • G05B23/0235—Qualitative history assessment, whereby the type of data acted upon, e.g. waveforms, images or patterns, is not relevant, e.g. rule based assessment; if-then decisions based on a comparison with predetermined threshold or range, e.g. "classical methods", carried out during normal operation; threshold adaptation or choice; when or how to compare with the threshold
    • F—MECHANICAL ENGINEERING; LIGHTING; HEATING; WEAPONS; BLASTING
    • F01—MACHINES OR ENGINES IN GENERAL; ENGINE PLANTS IN GENERAL; STEAM ENGINES
    • F01D—NON-POSITIVE DISPLACEMENT MACHINES OR ENGINES, e.g. STEAM TURBINES
    • F01D21/00—Shutting-down of machines or engines, e.g. in emergency; Regulating, controlling, or safety means not otherwise provided for
    • F01D21/003—Arrangements for testing or measuring
    • F—MECHANICAL ENGINEERING; LIGHTING; HEATING; WEAPONS; BLASTING
    • F01—MACHINES OR ENGINES IN GENERAL; ENGINE PLANTS IN GENERAL; STEAM ENGINES
    • F01D—NON-POSITIVE DISPLACEMENT MACHINES OR ENGINES, e.g. STEAM TURBINES
    • F01D21/00—Shutting-down of machines or engines, e.g. in emergency; Regulating, controlling, or safety means not otherwise provided for
    • F01D21/14—Shutting-down of machines or engines, e.g. in emergency; Regulating, controlling, or safety means not otherwise provided for responsive to other specific conditions
    • F—MECHANICAL ENGINEERING; LIGHTING; HEATING; WEAPONS; BLASTING
    • F01—MACHINES OR ENGINES IN GENERAL; ENGINE PLANTS IN GENERAL; STEAM ENGINES
    • F01D—NON-POSITIVE DISPLACEMENT MACHINES OR ENGINES, e.g. STEAM TURBINES
    • F01D25/00—Component parts, details, or accessories, not provided for in, or of interest apart from, other 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
    • G01M15/00—Testing of engines
    • G01M15/14—Testing gas-turbine engines or jet-propulsion engines
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F18/00—Pattern recognition
    • G06F18/20—Analysing
    • G06F18/22—Matching criteria, e.g. proximity measures
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F18/00—Pattern recognition
    • G06F18/20—Analysing
    • G06F18/24—Classification techniques
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F18/00—Pattern recognition
    • G06F18/20—Analysing
    • G06F18/24—Classification techniques
    • G06F18/241—Classification techniques relating to the classification model, e.g. parametric or non-parametric approaches
    • G06F18/2413—Classification techniques relating to the classification model, e.g. parametric or non-parametric approaches based on distances to training or reference patterns
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F18/00—Pattern recognition
    • G06F18/20—Analysing
    • G06F18/24—Classification techniques
    • G06F18/241—Classification techniques relating to the classification model, e.g. parametric or non-parametric approaches
    • G06F18/2413—Classification techniques relating to the classification model, e.g. parametric or non-parametric approaches based on distances to training or reference patterns
    • G06F18/24133—Distances to prototypes
    • G06F18/24137—Distances to cluster centroïds
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T11/00—2D [Two Dimensional] image generation
    • G06T11/20—Drawing from basic elements, e.g. lines or circles
    • G06T11/206—Drawing of charts or graphs
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T7/00—Image analysis
    • G06T7/0002—Inspection of images, e.g. flaw detection
    • G06T7/0004—Industrial image inspection
    • G06T7/001—Industrial image inspection using an image reference approach
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06V—IMAGE OR VIDEO RECOGNITION OR UNDERSTANDING
    • G06V10/00—Arrangements for image or video recognition or understanding
    • G06V10/70—Arrangements for image or video recognition or understanding using pattern recognition or machine learning
    • G06V10/74—Image or video pattern matching; Proximity measures in feature spaces
    • G06V10/761—Proximity, similarity or dissimilarity measures
    • F—MECHANICAL ENGINEERING; LIGHTING; HEATING; WEAPONS; BLASTING
    • F05—INDEXING SCHEMES RELATING TO ENGINES OR PUMPS IN VARIOUS SUBCLASSES OF CLASSES F01-F04
    • F05D—INDEXING SCHEME FOR ASPECTS RELATING TO NON-POSITIVE-DISPLACEMENT MACHINES OR ENGINES, GAS-TURBINES OR JET-PROPULSION PLANTS
    • F05D2260/00—Function
    • F05D2260/80—Diagnostics
    • F—MECHANICAL ENGINEERING; LIGHTING; HEATING; WEAPONS; BLASTING
    • F05—INDEXING SCHEMES RELATING TO ENGINES OR PUMPS IN VARIOUS SUBCLASSES OF CLASSES F01-F04
    • F05D—INDEXING SCHEME FOR ASPECTS RELATING TO NON-POSITIVE-DISPLACEMENT MACHINES OR ENGINES, GAS-TURBINES OR JET-PROPULSION PLANTS
    • F05D2270/00—Control
    • F05D2270/70—Type of control algorithm
    • F05D2270/708—Type of control algorithm with comparison tables
    • 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
    • Y02T—CLIMATE CHANGE MITIGATION TECHNOLOGIES RELATED TO TRANSPORTATION
    • Y02T50/00—Aeronautics or air transport
    • Y02T50/60—Efficient propulsion technologies, e.g. for aircraft

Landscapes

  • Engineering & Computer Science (AREA)
  • Theoretical Computer Science (AREA)
  • General Physics & Mathematics (AREA)
  • Physics & Mathematics (AREA)
  • Data Mining & Analysis (AREA)
  • General Engineering & Computer Science (AREA)
  • Computer Vision & Pattern Recognition (AREA)
  • Evolutionary Computation (AREA)
  • Artificial Intelligence (AREA)
  • Bioinformatics & Cheminformatics (AREA)
  • Life Sciences & Earth Sciences (AREA)
  • Evolutionary Biology (AREA)
  • Bioinformatics & Computational Biology (AREA)
  • Mechanical Engineering (AREA)
  • Automation & Control Theory (AREA)
  • Multimedia (AREA)
  • Computing Systems (AREA)
  • Health & Medical Sciences (AREA)
  • General Health & Medical Sciences (AREA)
  • Databases & Information Systems (AREA)
  • Software Systems (AREA)
  • Medical Informatics (AREA)
  • Quality & Reliability (AREA)
  • Chemical & Material Sciences (AREA)
  • Combustion & Propulsion (AREA)
  • Testing And Monitoring For Control Systems (AREA)

Description

本発明は、プラントや設備などの出力する多次元時系列データをもとに異常を早期に検知し、検知された異常に関連するセンサを特定する異常診断方法およびその装置に関する。
電力会社では、ガスタービンの廃熱などを利用して地域暖房用温水を供給したり、工場
向けに高圧蒸気や低圧蒸気を供給したりしている。石油化学会社では、ガスタービンなどを電源設備として運転している。このようにガスタービンなどを用いた各種プラントや設備において、設備の不具合あるいはその兆候を検知する異常検知は、社会へのダメージを最小限に抑えるためにも極めて重要である。
ガスタービンや蒸気タービンのみならず、水力発電所での水車、原子力発電所の原子炉、風力発電所の風車、航空機や重機のエンジン、鉄道車両や軌道、エスカレータ、エレベータ、機器・部品レベルでも、搭載電池の劣化・寿命など、上記のような予防保全を必要とする設備は枚挙に暇がない。
このため、対象設備やプラントに複数のセンサを取り付け、センサ毎の監視基準に従って正常か異常かを判断することが行われている。米国特許第６，９５２，６６２号（特許文献１）には、観測データと過去の学習データとの類似度を独自の方法で計算し、類似度に応じて重み付けされたデータの線形結合により推定値を算出して、推定値と観測データのはずれ度合をもとに異常を検出する方法が開示されている。また、特開２０１１−７０６３５号公報（特許文献２）には、過去の正常データから作成されたモデルとの比較によって算出される異常測度に基づいて異常の有無を検知する異常検知方法において、正常モデルを局所部分空間法によって作成することが開示されている。
しかし、異常を検知するのみでは対策や調査などの次のアクションが決められないため、どのセンサが異常に関連しているのかを診断したいというニーズがある。
このようなニーズに対応する技術として、特開２０１２−１３８０４４号公報（特許文献３）には、複数の測定変数から算出されるＴ２統計量およびＱ統計量に基づいて異常を検出し、異常検出後にこれらの統計量に対する各測定変数の寄与量に基づいて異常要因候補を列挙するプロセス状態監視装置が開示されている。
米国特許第６，９５２，６６２号明細書 特開２０１１−７０６３５号公報 特開２０１２−１３８０４４号公報
センサ毎の監視基準に従って正常か異常かを判断する方法では、監視基準を超えたセンサがすなわち異常に関連するセンサであるが、この方法では２種以上のセンサ信号の組合せとしての異常はそもそも検知できないという問題がある。
特許文献１や特許文献２に記載の方法によれば、学習データとして正常時のデータを与えることにより、学習にない観測データが観察されると、これらを異常として検出することができる。複数のセンサ信号をまとめて演算するため、センサ信号の組合せとしての異常も検出できる。しかし、異常に関連するセンサを正しく診断する方法については記載がない。
特許文献３に記載の方法によれば、複数の測定項目のＴ２統計量やＱ統計量を異常の指標とすることにより、２種以上のセンサ信号の組合せとしての異常も検知でき、統計量に対する各測定項目の寄与量に基づいて異常に関連するセンサを診断することが可能である。しかし、Ｔ２統計量やＱ統計量を指標とした異常検出は、運転と停止の切り替りや運転条件の切り替りが頻繁かつ不規則にある設備やプラントに適用すると、誤報が多くなってしまい、そのままでは適用困難である。特許文献２に記載の局所部分空間法は、学習データを網羅的に集めることにより、上記のようなケースでも誤報を抑制することが可能である。
しかし、特許文献１に記載の異常検知手法に、特許文献２に記載の寄与量に基づいて異常に関連するセンサを診断する方法を適用すると、誤る場合がある。この理由について簡単に説明する。局所部分空間法は、図２に示す通り、観測データ(q)に近い学習データ(x1，x2，x3)を用いて作られるアフィン部分空間への投影距離を異常の指標としている。つまり、観測データ(q)からアフィン部分空間への垂線の足(b)を基準として、bからqのベクトルの各要素の大きさが寄与量となる。あるセンサ信号に異常がある場合、これに基準(b)が追随して異常とは関係ない他のセンサの寄与量が大きくなることが起こりうる。
単純化した例を図３Ａ及び図３Ｂに示す。図３Ａにおいて、abcの３次元のセンサ信号において正常なデータが図示するように分布しており、観測データがq→q’→qと変化した場合、基準データはb→b’→bと変化する。時系列のグラフで表すと図３Ｂに示す通りとなるため、aとbの寄与量が大きくなるが、実際にはcが異常に関連していると考えられるため、誤って診断されることになる。このことは、局所部分空間法のみでなく、ＡＲＴ２のようなクラスタに基づく方法や特許文献１に記載の類似度ベースモデリングと呼ばれる手法など、観測データによって基準が変化する事例ベースの手法すべてにあてはまる。
そこで、本発明の目的は、上記課題を解決し、多次元時系列センサ信号に基づいて高感度に異常を検出し、検出した異常に関連するセンサを正しく診断する異常診断方法および装置を提供することにある。
上記課題を解決するために、本発明では、設備または装置に装着された複数のセンサか
ら出力される複数のセンサ信号に基づいて設備または装置の異常を診断する方法及びその
装置において、予め指定された学習期間において、複数のセンサ信号をそれぞれ正準化して得た特徴値から多次元の特徴ベクトルを求め、この多次元の特徴ベクトルから注目ベクトルを選び、この注目ベクトルの近傍ベクトルを用いて基準ベクトルを作成し、多次元の特徴ベクトルの基準ベクトルまでの距離に基づいて異常測度を算出し、算出した異常測度に基づいてしきい値を設定し、多次元の特徴ベクトルの特徴値から算出した値の頻度分布に基づいて分布密度画像を生成し、設備または装置の異常を診断する時において、複数のセンサ信号から時刻毎に多次元特徴ベクトルを抽出し、求めた予め指定された学習期間の多次元の特徴ベクトルの中から時刻毎の多次元の特徴ベクトルに近い特徴ベクトルを選択し、前記選択した特徴ベクトルを用いて基準特徴ベクトルを算出し、抽出した時刻毎の多次元特徴ベクトルの算出した基準特徴ベクトルまでの距離に基づき異常測度を算出し、算出した異常測度と設定したしきい値とを比較することにより異常を検出し、複数のセンサ信号の学習期間において求めた特徴ベクトルを用いて作成した２次元の分布密度に基づいて検出した異常に関連するセンサを複数のセンサの中から特定するようにした。
又、上記課題を解決するために、本発明では、設備または装置に装着された複数のセン
サから出力される複数のセンサ信号をもとに学習データを作成して蓄積する工程と、設備
または装置に装着された複数のセンサから新たに出力された複数のセンサ信号の異常診断
を行う工程とを含む設備の異常を診断する方法及びその装置において、学習データを作成
して蓄積する工程は、複数のセンサ信号から特徴ベクトルを抽出する工程と、抽出した特
徴ベクトルを学習データとして蓄積する工程と、学習データとして蓄積した特徴ベクトルの各々について交差検証用のグループ番号付けをし、蓄積した学習データの中から注目ベクトルを選択する工程と、蓄積した学習データの中の注目ベクトルと異なるグループ番号付けをした特徴ベクトルの中から選択された注目ベクトルの近傍ベクトルを選択してこの選択した近傍ベクトルを用いて学習用の基準ベクトルを作成する工程と、作成した学習用の基準ベクトルまでの学習データとして蓄積した特徴ベクトルからの距離に基づいて異常測度を算出する工程と、算出した異常測度に基づきしきい値を算出する工程と、特徴ベクトルから全組合せの２次元の特徴分布密度を算出する工程とを有し、センサ信号の異常診断を行う工程は、設備または装置に装着された複数のセンサから出力される複数のセンサ信号から特徴ベクトルを抽出して観測ベクトルとする工程と、抽出した観測ベクトルに応じて蓄積した学習データの中から所定数の特徴ベクトルを選択する工程と、学習データの中から選択された所定数の特徴ベクトルを用いて異常診断用の基準ベクトルを作成する工程と、観測ベクトルから異常診断用の基準ベクトルまでの距離に基づいて観測ベクトルの異常測度を算出する工程と、算出した異常測度と学習データを作成して蓄積する工程において算出したしきい値に基づき観測ベクトルが異常か正常かを判定する工程と、判定する工程において異常と判定された観測ベクトルと２次元の特徴分布密度に基づき異常に関連するセンサを特定する工程とを有して構成した。
さらにまた、上記課題を解決するために、本発明では、設備または装置に装着された複
数のセンサから出力される複数のセンサ信号をもとに学習データを作成して蓄積する工程
と、前記設備または装置に装着された複数のセンサから新たに出力された複数のセンサ信
号の異常診断を行う工程とを含む設備の異常を診断する方法及びその装置において、学習
データを作成して蓄積する工程は、設備または装置から出力されるイベント信号に基づい
て稼動状態別のモード分割を行う工程と、複数のセンサ信号から特徴ベクトルを抽出する
工程と、抽出した特徴ベクトルを学習データとして蓄積する工程と、学習データとして蓄積した特徴ベクトルのおのおのについて交差検証用のグループ番号付けをして特徴ベクトルに応じて蓄積した学習データの中から特徴ベクトルを選択する工程と、蓄積した学習データの中の選択された特徴ベクトルと異なるグループ番号付けをした特徴ベクトルの中から選択された特徴べクトルの近傍ベクトルを選択してこの選択した近傍ベクトルを用いて学習用の基準ベクトルを作成する工程と、作成した学習用の基準ベクトルまでの学習データとして蓄積した特徴ベクトルからの距離に基づいて異常測度を算出する工程と、算出した異常測度に基づきモード毎にしきい値を算出する工程と、学習データとして蓄積した特徴ベクトルから全組合せの２次元の特徴分布密度をモード毎に算出する工程を有し、センサ信号の異常診断を行う工程は、設備または装置に装着された複数のセンサから出力される複数のセンサ信号から特徴ベクトルを抽出して観測ベクトルとする工程と、抽出した観測ベクトルに応じて蓄積した学習データの中から所定数の特徴ベクトルを選択する工程と、学習データの中から選択された所定数の特徴ベクトルの近傍ベクトルを用いて異常診断用の基準ベクトルを作成する工程と、観測ベクトルの異常診断用の基準ベクトルまでの距離に基づいて観測ベクトルの異常測度を算出する工程と、算出した異常測度と分割したモード別に算出したしきい値に基づき観測ベクトルが異常か正常かを判定する工程と、判定する工程において異常と判定された観測ベクトルに対応するセンサ信号とモード別に算出した２次元の特徴分布密度に基づき異常に関連するセンサを特定する工程とを有して構成した。
本発明によれば、各時刻のデータの異常検出を、学習データとその時刻のデータから算出される基準データからの距離に基づいて行うため、高感度な異常検出が可能である。異常が検出された時刻の関連センサ特定において、寄与量が大きいセンサをそのまま関連センサとするのではなく、センサ信号の２次元の分布密度に基づいて異常に関連するセンサを特定するため、正しく診断することができる。
以上の手法を適用したシステムにより、ガスタービンや蒸気タービンなどの設備のみならず、水力発電所での水車、原子力発電所の原子炉、風力発電所の風車、航空機や重機のエンジン、鉄道車両や軌道、エスカレータ、エレベータ、そして機器・部品レベルでは、搭載電池の劣化・寿命など、種々の設備・部品において、対象の異常を早期に検出するのみならず、関連センサに応じた調査あるいは対策を迅速に行うことが可能となる。
本発明の実施例１に関する設備状態監視システムの概略の構成を示すブロック図である。 局所部分空間法を説明する図である。 寄与量に基づくセンサ信号特定における誤りの例を説明する図で、正常データの分布を示す３次元のグラフである。 寄与量に基づくセンサ信号特定における誤りの例を説明する図で、３つのセンサ信号の基準値に対する観測値の変化を示すセンサ信号のグラフである。 センサ信号の例を示す信号リストの表である。 センサ信号入力からしきい値設定までの処理の流れを示すフロー図である。 ２次元の特徴の分布密度を算出する処理の流れを示すフロー図である。 センサ信号入力から異常判定までの処理の流れを示すフロー図である。 異常に関連するセンサを特定する処理の流れを示すフロー図の例である。 異常に関連するセンサを特定する処理の流れを示すフロー図の別の例である。 異常に関連する３個のセンサを特定する処理の流れを示すフロー図である。 異常が連続する区間毎に関連するセンサを特定する処理の流れを示すフロー図である。 分布密度画像に異常区間の対応点をプロットした例を表す図である。 分布密度画像に異常区間の対応点をプロットした例を表す図である。 分布密度画像に異常区間の対応点をプロットした例を表す図である。 分布密度画像に異常区間の対応点をプロットした例を表す図である。 レシピ設定のためのＧＵＩの１例を表す表示画面の正面図である。 結果表示画面の例で複数の時系列データを表示した画面の正面図である。 結果表示画面の例で複数の時系列データを拡大表示した画面の正面図である。 結果表示画面の例で異常関連センサの分布密度と複数の時系列データを表示した画面の正面図である。 レシピ設定のためのテスト結果の一覧を表示した画面の正面図である。 結果表示対象指定のためのＧＵＩの１例を表す表示画面の正面図である。 結果表示画面に含まれる期間表示ウィンドウの表示例である。 本発明の実施例２に関する設備状態監視システムの概略の構成を示すブロック図である。 本発明の実施例２に関するイベント信号の例を示す信号リストの表である。 本発明の実施例２に関するイベント信号に基づくモード分割処理の流れを示すフロー図である。 本発明の実施例２に関する設備の可動状態を分割して４種のモードの何れかに分類した状態を示すイベント信号の模式図である。
本発明は、設備または装置に装着された複数のセンサから出力される複数のセンサ信号に基づいて設備または装置の異常を診断する方法及びその装置において、センサ信号から時刻毎に多次元の特徴ベクトルを抽出し、予め指定された学習期間の特徴ベクトルの集合と各時刻の特徴ベクトルに基づいて各時刻の基準ベクトルを算出し、各時刻の特徴ベクトルと基準特徴ベクトルとの差に基づいて異常測度を算出し、異常測度と所定のしきい値との比較により異常を検出し、異常が検出された時刻における特徴値の２次元の分布密度に基づいて異常に関連するセンサを特定できるようにしたものである。
以下に、本発明の実施例を図を用いて説明する。
以下、図面を用いて本発明の内容を詳細に説明する。
図１に、本発明の異常診断方法を実現するシステムの一構成例を示す。
本システムは、設備１０１から出力されるセンサ信号１０２を蓄積するセンサ信号蓄積部１０３、センサ信号１０２をもとに特徴ベクトルを抽出する特徴ベクトル抽出部１０４、予め指定された学習期間の特徴ベクトルの集合と各時刻の特徴ベクトルに基づいて各時刻の基準特徴ベクトルを算出する基準ベクトル算出部１０５、各時刻の特徴ベクトルと基準特徴ベクトルの差に基づき異常測度を算出する異常測度算出部１０６、予め指定された学習期間の異常測度に基づきしきい値を算出するしきい値算出部１０７、異常測度と算出されたしきい値との比較により異常を検出する異常検出部１０８、学習期間のセンサ信号の２次元の分布密度を算出する分布密度算出部１０９、異常が検出された時刻について、異常に関連するセンサを特定する関連センサ特定部１１０を備えて構成される。
本システムの動作には、蓄積されたデータを用いて学習データの生成、保存を行う「学習」と入力信号に基づき異常を検出し、関連センサの特定を行う「異常診断」の二つのフェーズがある。基本的に前者はオフラインの処理、後者はオンラインの処理である。ただし、後者をオフラインの処理とすることも可能である。以下の説明では、それらを学習時、異常診断時という言葉で区別する。
状態監視の対象とする設備１０１は、ガスタービンや蒸気タービンなどの設備やプラントである。設備１０１は、その状態を表すセンサ信号１０２を出力する。センサ信号１０２はセンサ信号蓄積部１０３に蓄積されている。センサ信号１０２をリスト化して表形式に表した例を図４に示す。センサ信号１０２は一定間隔毎に取得される多次元時系列信号であり、それをリスト化した表は、図４に示すように、日時の欄４０１と設備１０１に設けられた複数のセンサ値のデータの欄４０２からなる。センサの種類は、数百から数千といった数になる場合もあり、例えば、シリンダ、オイル、冷却水などの温度、オイルや冷却水の圧力、軸の回転速度、室温、運転時間などである。出力や状態を表すのみならず、何かをある値に制御するための制御信号の場合もある。
学習時の処理の流れを図５を用いて説明する。学習時はセンサ信号蓄積部１０３に蓄積されたデータのうち指定された期間のデータを用いて、特徴ベクトルを抽出する。ここで抽出されたデータを学習データと呼ぶ。交差検証によって学習データの異常測度を算出し、全学習データの異常測度に基づき異常判定のしきい値を算出する。また、分布密度算出部１０９において特徴値の２次元の分布密度を算出する。
まず、特徴ベクトル抽出部１０４、基準ベクトル算出部１０５、異常測度算出部１０６、しきい値算出部１０７における処理の流れを図５を用いて説明する。はじめに、特徴ベクトル抽出部１０４において、センサ信号蓄積部１０３から学習期間として指定された期間のセンサ信号１０２を入力し（Ｓ５０１）、センサ信号毎に正準化した後（Ｓ５０２）、特徴ベクトルの抽出を行う（Ｓ５０３）。次に、基準ベクトル算出部１０５において、抽出された特徴ベクトルに交差検証用のグループ番号を付与する（Ｓ５０４）。グループは、１日分を１グループとするなど、期間を指定して決めるか、予め指定されたグループ数に等分に分けて決める。
抽出された特徴ベクトルから、１個目の特徴ベクトルに注目し（Ｓ５０５）、注目ベクトルと異なるグループの特徴ベクトルの中から、注目ベクトルに近い順に予め指定された数の特徴ベクトルを選択し（Ｓ５０６）、選択した特徴ベクトルを用いて基準ベクトルを作成する（Ｓ５０７）。異常測度算出部１０６において、注目した特徴ベクトルの基準ベクトルまでの距離に基づいて異常測度を算出する（Ｓ５０８）。全ベクトルの異常測度算出が終了していれば（Ｓ５０９）、しきい値算出部１０７において、全ベクトルの異常測度に基づいてしきい値を設定する（Ｓ５１１）。ステップＳ５０９において全ての異常測度算出が終了していなければ、次の特徴ベクトルに注目し（Ｓ５１０）、ステップＳ５０６からＳ５０９の処理を繰り返す。
次に、各ステップについて詳細に説明する。
ステップＳ５０２においては、各センサ信号の正準化を行う。例えば、指定された期間の平均と標準偏差を用いて、平均を０、分散を１となるように変換する。異常診断時に同じ変換ができるよう、各センサ信号の平均と標準偏差を記憶しておく。あるいは、各センサ信号の指定された期間の最大値と最小値を用いて最大が１、最小が０となるように変換する。あるいは、最大値と最小値の代わりに予め設定した上限値と下限値を用いてもよい。異常診断時に同じ変換ができるよう、各センサ信号の最大値と最小値または上限値と下限値を記憶しておく。センサ信号の正準化は、単位およびスケールの異なるセンサ信号を同時に扱うためのものである。
ステップＳ５０３においては、時刻毎に特徴ベクトル抽出を行う。センサ信号を正準化したものをそのまま並べることが考えられるが、ある時刻に対して±1，±2，…のウィンドウを設け，[ウィンドウ幅（3，5，…）×センサ数]の特徴ベクトルにより、データの時間変化を表す特徴を抽出することもできる。また、離散ウェーブレット変換(DWT: Discrete Wavelet Transform)を施して、周波数成分に分解してもよい。さらに、ステップＳ５０３において、特徴選択を行う。最低限の処理として、分散が非常に小さいセンサ信号および単調増加するセンサ信号を除く必要がある。
また、相関解析による無効信号を削除することも考えられる。これは、多次元時系列信号に対して相関解析を行い、相関値が１に近い複数の信号があるなど、極めて類似性が高い場合に、これらは冗長だとして、この複数の信号から重複する信号を削除し、重複しないものを残す方法である。このほか、ユーザが指定するようにしてもよい。また、長期変動が大きい特徴を除くことも考えられる。長期変動が大きい特徴を用いることは正常状態の状態数を多くすることにつながり、学習データの不足を引き起こすためである。例えば、１周期期間毎の平均と分散を算出し、それらのばらつきによって長期変動の大きさを推定できる。
基準ベクトル作成手法としては、局所部分空間法(LSC: Local Sub-space Classifier)や投影距離法(PDM: Projection Distance Method)が考えられる。
局所部分空間法は、注目ベクトルqのk−近傍ベクトルを用いてｋ−１次元のアフィン部分空間を作成する方法である。図２にｋ＝３の場合の例を示す。k-近傍ベクトルはすなわちステップＳ５０６において選択された特徴ベクトルである。指定する数は、ｋの値である。
図２に示すように、異常測度は図に示す投影距離で表されるため、注目ベクトルqに最も近いアフィン部分空間上の点bを求めればよい。評価データqとそのk−近傍ベクトルxi( i＝1，…，k )からbを算出するには、qをk個並べた行列Qとxiを並べた行列Xから
Figure 0006076751
 により相関行列Cを求め、
Figure 0006076751
 によりbを計算する。ここまでが、ステップＳ５０７における基準ベクトル作成にあたる。
異常測度ｄはｑとｂの間の距離であるから次式で表される。
Figure 0006076751
なお、図２ではｋ＝３の場合を説明したが、特徴ベクトルの次元数より十分小さければいくつでもよい。ｋ＝１の場合は、最近傍法と等価の処理になる。
投影距離法は、選択された特徴ベクトルに対し独自の原点をもつ部分空間すなわちアフィン部分空間（分散最大の空間）を作成する方法である。ステップＳ５０６において指定する数はいくつでも良いが、大きすぎるとベクトルの選択、部分空間の算出ともに時間がかかってしまうので数十から数百とするとよい。
アフィン部分空間の算出方法について説明する。まず、選択された特徴ベクトルの平均μと共分散行列Σ を求め、次にΣの固有値問題を解いて値の大きい方から予め指定したr個の固有値に対応する固有ベクトルを並べた行列Ｕをアフィン部分空間の正規直交基底とする。ｒは特徴ベクトルの次元より小さくかつ選択データ数より小さい数とする。あるいはｒを固定した数とせず、固有値の大きい方から累積した寄与率が予め指定した割合を超えたときの値としてもよい。異常測度は、注目ベクトルのアフィン部分空間への投影距離とする。
この他、注目ベクトルqのk-近傍ベクトルの平均ベクトルまでの距離を異常測度とする局所平均距離法や、ガウシアンプロセスなどを用いてもよい。
ステップＳ５１１におけるしきい値を設定する方法について説明する。学習期間の全特徴ベクトルの異常測度を昇順にソートし、予め指定した１に近い比率に到達する値を求める。この値を基準としてオフセットを加える、定数倍するなどの処理によりしきい値を算出する。オフセット０、倍率を１とすれば、この値がしきい値となる。算出したしきい値は、図示していないが、学習データと対応付けて記録しておく。
上記に示したように、注目ベクトルの近傍ベクトルを用いて基準ベクトルを作成することにより、状態が複雑に変化する設備に対しても適切な基準により精度の高い異常測度を算出することが可能である。学習データの交差検証により算出された異常測度に基づいてしきい値を算出するため、誤報を抑制することができる。
次に、分布密度算出部１０９における処理の流れを図６を用いて説明する。
始めに、学習期間の特徴ベクトルを入力する（Ｓ６０１）。最初の特徴に注目し（Ｓ６０２）、学習期間のデータの最大値（ＭＡＸ）と最小値（ＭＩＮ）を求める（Ｓ６０３）。次に、最小値から最大値を指定された数Ｎで分割する際の刻み幅Ｓを算出する（Ｓ６０４）。Ｓ＝（ＭＡＸ−ＭＩＮ）／Ｎ で計算できる。次に、最小値と最大値から外側に範囲を広げて分布密度算出の処理範囲を算出する（Ｓ６０５）。広げる範囲は例えばＭＩＮをＭＩＮ−Ｓ×Ｍ、ＭＡＸをＭＡＸ＋Ｓ×Ｍに変更する。ここでＭは予め決められた１以上の整数とする。
次に学習期間の全データについて、特徴値（Ｆ）からビン番号（ＢＮＯ）を次式で算出する（Ｓ６０６）。
ＢＮＯ＝ＩＮＴ（（Ｆ−ＭＩＮ）・Ｎ／（ＭＡＸ−ＭＩＮ））
ただしＩＮＴ（Ｘ）はＸの整数部を表す。
以上の処理について未処理の特徴が残っている場合は（Ｓ６０７）、次の特徴に注目し（Ｓ６０８）、ステップＳ６０３からＳ６０６の処理を行う。ステップＳ６０７において全特徴について処理していれば、ステップＳ６０９に進む。
ステップＳ６０９においては、特徴２個の最初の組合せに注目する。２個の特徴は同じものであってもよい。次に分布密度算出用の２次元配列を確保し、すべての要素に０をセットする（Ｓ６１０）。配列のサイズはＮ＋２Ｍである。学習期間の全データについて、２個の特徴のビン番号に対応する配列の要素に１を加算する（Ｓ６１１）。この処理により、特徴２個による２次元の頻度分布（ヒストグラム）が算出される。この頻度分布を画像に変換して保存する（Ｓ６１２）。変換方法については後述する。上記処理が特徴２個の全ての組合せについて実施済みであれば（Ｓ６１３）、処理を終了し（Ｓ６１５）、そうでなければ次の組合せに注目し（Ｓ６１４）、ステップＳ６１０からＳ６１２までの処理を行う。図示はしていないが、２次元配列のサイズおよびステップＳ６０５で算出した各特徴の最小値と最大値を記録しておく。
ステップＳ６１２における、画像変換方法の例を説明する。始めに配列要素の最大値すなわち最大頻度を求める。画像サイズは配列サイズと同じとし、各要素の値から対応する座標の画素値を例えば、[２５５×配列の要素値／最大頻度]とする。２５５は画素値を８ビットで表す場合の最大値であり、この値を用いれば、そのままビットマップ形式で保存できる。あるいは、画素値を[２５５×ＬＯＧ（配列の要素値＋１）／ＬＯＧ（最大頻度＋１）]とする。ただしＬＯＧ（Ｘ）はＸの対数を表す。このような変換式を用いれば、最大頻度が大きい場合も非ゼロの頻度に非ゼロの画素値を対応させることが可能になる。
上記処理により得られた画像は、２次元の特徴空間上で密度が高いところが高い画素値で表されているため、分布密度画像と呼ぶこととする。画像の作り方は、上記方法に限定されない。例えば単純な頻度分布ではなく、１個のデータにガウス分布や他の重みつきフィルタを割り当て、それを重畳するようにしてもよい。あるいは、上記方法で得られた画像に所定サイズの最大値フィルタをかけたり、平均フィルタ、その他の重みつきフィルタをかけたりしてもよい。また必ずしも画像形式で保存する必要はなく２次元配列をテキスト形式で保存してもよい。
異常診断時の処理の流れを図７ないし図９を用いて説明する。異常診断時はセンサ信号蓄積部１０３に蓄積されたデータのうち指定された期間のデータあるいは新たに観測されたデータの異常測度を算出し、正常か異常かの判定を行い、異常と判定された時刻について関連するセンサを特定する。
図７は特徴ベクトル抽出部１０４、基準ベクトル算出部１０５、異常測度算出部１０６、異常検出部１０８における異常診断時の処理の流れを説明する図である。始めに、特徴ベクトル抽出部１０４において、センサ信号蓄積部１０３または設備１０１からセンサ信号を入力し（Ｓ７０１）、センサ信号毎に正準化した後（Ｓ７０２）、特徴ベクトルの抽出を行う（Ｓ７０３）。センサ信号の正準化には、図５のステップＳ５０２に示す処理において、学習データの正準化に用いたものと同じパラメータを用いる。特徴ベクトルの抽出は、ステップＳ５０３の処理と同じ方法で行う。したがって、ステップＳ５０３において特徴選択を実行した場合は同じ特徴を選択する。ここで抽出された特徴ベクトルを、学習データと区別するために観測ベクトルと呼ぶこととする。
次に、基準ベクトル作成部１０５において、学習データの特徴ベクトルの中から、観測ベクトルに近い順に予め指定された数の特徴ベクトルを選択し（Ｓ７０４）、それらの特徴ベクトルを用いて基準ベクトルを作成する（Ｓ７０５）。異常測度算出部１０６において、観測ベクトルの基準ベクトルまでの距離に基づいて異常測度を算出する（Ｓ７０６）。異常検出部１０８において、学習時に算出したしきい値と異常測度を比較して、異常測度がしきい値より大きければ異常、そうでなければ正常と判定する（Ｓ７０７）。
図８は、関連センサ特定部１１０における処理の流れの一実施例を説明する図である。この処理は、ステップＳ７０７において、異常と判定されたデータについて実施する。始めに、観測ベクトルと基準ベクトルから、各特徴の観測値と基準値の差を算出し（Ｓ８０１）、これが最大となる特徴を探す（Ｓ８０２）。これを特徴Ａとしておく。次に観測ベクトルから、図６で示した分布密度算出処理で使用した２次元配列のサイズおよびステップＳ６０５で算出した各特徴の最小値と最大値を用いて、各特徴のビン番号を算出する（Ｓ８０３）。一方の特徴が特徴Ａである分布密度画像を読み込み、各特徴と特徴Ａのビン番号に対応する座標の画素値を読み（Ｓ８０４）、これが最小となる特徴を探す（Ｓ８０５）。これを特徴Ｂとし、特徴Ａと特徴Ｂの分布密度画像に観測ベクトルの対応点をプロットする（Ｓ８０６）。例えば、カラー画像に変換して、特徴Ａと特徴Ｂのビン番号に対応する座標に使用していない色をつける。Ｓ８０４における画素値が０であれば、あるいは予め決められた値より小さければ、特徴Ａと特徴Ｂを関連センサとする（Ｓ８０７）。
図９は、関連センサ特定部１１０における処理の流れの別の一実施例を説明する図である。図８の処理と同様、ステップＳ７０７において、異常と判定されたデータについて実施する。始めに、観測ベクトルから、図６で示した分布密度算出処理で使用した２次元配列のサイズおよびステップＳ６０５で算出した各特徴の最小値と最大値を用いて、各特徴のビン番号を算出する（Ｓ９０１）。
次に全ての分布密度画像を読み込み、それぞれの分布密度画像について２個のビン番号に対応する座標の画素値を読み（Ｓ９０２）、特徴毎にその値の合計siを算出する（Ｓ９０３）。つまり、特徴iと特徴jの分布密度画像から読み取った値をI(i，j)とすると、siは次式で計算される。
Figure 0006076751
このsiが最小となる特徴を探し（Ｓ９０４）、これを特徴Ｃとしておく。次に一方の特徴が特徴Ｃである分布密度画像について、各特徴と特徴Ｃのビン番号に対応する座標の画素値を読み（Ｓ９０５）、これが最小となる特徴を探す（Ｓ９０６）。これを特徴Ｄとし、ステップＳ８０６と同様の方法で特徴Ｃと特徴Ｄの分布密度画像に観測ベクトルの対応点をプロットする（Ｓ９０７）。Ｓ９０５における画素値が０であれば、あるいは予め決められた値より小さければ、特徴Ｃと特徴Ｄを関連センサとする（Ｓ９０８）。
また別の実施例では、ステップＳ９０１で各特徴のビン番号を算出し、ステップＳ９０２で全ての分布密度画像について観測ベクトルに対応する座標の画素値を読んだ後、これが最小となる特徴の組合せを探し、これらを特徴Ｅおよび特徴Ｆとする。ステップＳ８０６と同様の方法で特徴Ｅと特徴Ｆの分布密度画像に観測ベクトルの対応点をプロットし、特徴Ｅと特徴Ｆの分布密度画像から読んだ画素値が０であれば、あるいは予め決められた値より小さければ、特徴Ｅと特徴Ｆを関連センサとする。
以上、説明した関連センサ特定の方法はいずれかひとつを実施してもよいし、関連センサが特定できるまでシーケンシャルに実施してもよいし、全て実施して観測ベクトルをプロットした分布密度画像を表示して、ユーザに決めさせてもよい。
さらに、上記方法で関連センサが特定できない場合に詳細に解析する方法について、図１０を用いて説明する。始めに、学習期間の特徴ベクトルを入力する（Ｓ１００１）。次に学習期間の全データについて、図６で示した分布密度算出処理で使用した２次元配列のサイズおよびステップＳ６０５で算出した各特徴の最小値と最大値を用いて、各特徴のビン番号を算出する（Ｓ１００２）。同様に、観測ベクトルについて各特徴のビン番号を算出する（Ｓ１００３）。
次に、学習期間の特徴ベクトルから前述の特徴Ａのビン番号が観測ベクトルと同じデータのみ抽出し（Ｓ１００４）、これらを用いてステップＳ６１０からＳ６１２に従い、特徴Ａを除く組合せの部分的分布密度画像を作成する（Ｓ１００５）。全ての部分的分布密度画像について２個のビン番号に対応する座標の画素値を読み（Ｓ１００６）、特徴毎にその値の合計を算出する（Ｓ１００７）。これが最小となる特徴を探し（Ｓ１００８）、これを特徴Ｇとしておく。
一方の特徴が特徴Ｇである部分的分布密度画像について、各特徴と特徴Ｇのビン番号に対応する座標の画素値を読み（Ｓ１００９）、これが最小となる特徴を探し特徴Ｈとする（Ｓ１０１０）。特徴Ｇと特徴Ｈの部分的分布密度画像に観測ベクトルの対応点をプロットする（Ｓ１０１１）。Ｓ１００９における画素値が０であれば、あるいは予め決められた値より小さければ、特徴Ａと特徴Ｇと特徴Ｈを関連センサとする（Ｓ１０１２）。なお、これらの処理では特徴Ａのかわりに前述の特徴Ｃを用いてもよい。
上記方法は、異常検出した時刻毎に関連センサを特定する方法であるが、連続して異常検出された場合に適用すると冗長な処理になってしまう。そこで、異常が連続して検出されている間はバッファに情報をためておき、まとめて処理を行うことが考えられる。この場合、異常の継続時間や累積異常度でセンサ特定を行う対象を絞り込んでもよい。
連続して検出された一連の異常区間について関連するセンサを特定する方法を、図１１を用いて説明する。まず、異常検出部１０８においてステップＳ７０７で異常判定された結果をもとに、連続して異常検出された区間を特定する（Ｓ１１０１）。その間、観測ベクトルと基準ベクトルはバッファに蓄積しておく。なお「連続」というのは必ずしも厳密な意味ではなく、予め決められた長さの中断を含んでいてもよい。次に、上記区間に含まれる全ての観測ベクトルについて、特徴毎に観測値と基準値の差の絶対値の和を算出し（Ｓ１１０２）、これが最大となる特徴を探す（Ｓ１１０３）。これを特徴Ａとしておく。
次に、上記区間に含まれる全ての観測ベクトルについて、図６で示した分布密度算出処理で使用した２次元配列のサイズおよびステップＳ６０５で算出した各特徴の最小値と最大値を用いて、各特徴のビン番号を算出する（Ｓ１１０４）。一方の特徴が特徴Ａである分布密度画像を読み込み、上記区間に含まれる全ての観測ベクトルについて、各特徴と特徴Ａのビン番号に対応する座標の画素値を読む（Ｓ１１０５）。上記区間に含まれる全ての観測ベクトルについて特徴毎に非０の画素をカウントする（Ｓ１１０６）。特徴iについてカウントした値をcnt(i)とする。また、上記区間に含まれる全ての観測ベクトルについて特徴毎に画素値の合計を算出する（Ｓ１１０７）。特徴iについて算出した値をsum(i)とする。cnt(i)が最小となる特徴を探し、複数ある場合はそのうちsum(i)が最小となる特徴を探す（Ｓ１１０８）。これを特徴Ｂとし、特徴Ａと特徴Ｂの分布密度画像に異常区間に含まれる観測ベクトルの対応点をプロットする（Ｓ１１０９）。Ｓ１１０６における非０の画素数が区間に含まれる観測ベクトル数より小さければ、つまり０の画素数が１以上あれば、特徴Ａと特徴Ｂを関連センサとする（Ｓ１１１０）。
あるいは、cnt(i)が最小となる特徴が複数ある場合は、学習データの相関の高い特徴を探して特徴Ｂとしてもよい。相関の高さは例えば特徴Ａともう一方の特徴の相関値に基づいて決める。または、分布密度画像の非０の画素の数が少ないものを相関が高いと考える。
なお、上記に説明した図１１の処理フローは、図８に示す処理を複数の観測ベクトルに拡張した処理であるが、図９または図１０に示す処理を複数の観測ベクトルに拡張して用いてもよいし、それらを組み合わせて用いてもよい。
ステップＳ１１０９にて作成される分布密度画像の例を図１２Ａ〜Ｃに示す。横軸と縦軸はそれぞれの特徴の相対的な値を表し、色の濃い部分が学習データの密度が高い部分である。星印は異常区間に含まれる観測データを表す。実際には色をつけて表示すればよい。
図１２Ａは、特徴ａ、ｂの分布密度を表す。学習データの分布からは、特徴ａ、ｂとも低い状態と特徴ａ、ｂとも高い状態の定常状態があり状態遷移時は特徴ａに遅れて特徴ｂが上昇していることがわかる。時間ずれの関係は安定している。特徴ａ上昇のタイミングに対して特徴ｂの上昇が早く起こっていることが異常として検出されている。あるいは下降タイミングのずれを検出しているかもしれないが、どちらであるかはセンサ信号の時系列グラフから確認できる。
図１２Ｂは、特徴ｃ、ｄの分布密度を表す。特徴ｃ、ｄとも低い状態、特徴ｃが高く特徴ｄが低い状態、特徴ｃ、ｄとも高い状態とそれらの中間の状態はあるが特徴ｄが高く特徴ｃが低い状態にはならない。そのような状態が異常として検出されている。
図１２Ｃは、横軸も縦軸も特徴eの値を表したものである。実質的には２次元ではないが、同じ形式の画像として表すことにより、表示やセンサ特定の処理を簡便にすることができる。学習データは直線上に分布する。特徴eの値が学習データのどの値よりも高いことが異常として検出されている。
上記の通り、予め学習データを用いて全ての組合せの２次元特徴分布密度を算出しておき、異常が検出されたときの観測ベクトルの対応点の分布密度に基づいて異常に関連するセンサを特定するため、図３Ａ及び図３Ｂに示すような状況である場合でも正しいセンサを特定することが可能である。また、２次元特徴分布密度を画像で表し、これをユーザに提示することにより、それらのセンサが異常に関連することが理解されやすくなる。また、学習データの分布の傾向を知ることができる。例えば、図１２Ｄに示すように、学習データの密度の低い部分が広い範囲に存在する場合に、星印１２０１で示す異常を検出したとしてもそれが本当に異常かどうかは疑わなくてはならない。その場合、学習期間を増やして密度が高くなるようにするか、この組合せでは一緒に処理しないことにする、具体的には特徴ｇを処理対象からはずすようにするとよい。
以上の方法を実現するシステムのＧＵＩの実施例を説明する。
学習期間および処理パラメータ設定のためのＧＵＩの例を、図１３に示す。以下の説明ではこの設定のことを単にレシピ設定と呼ぶことにする。また、過去のセンサ信号１０２は設備ＩＤおよび時刻と対応付けられてデータベースに保存されているものとする。レシピ設定画面１３０１では、対象装置、学習期間、使用センサ、基準算出パラメータ、しきい値設定パラメータを入力する。設備ＩＤ入力ウィンドウ１３０２には、対象とする設備のＩＤを入力する。設備リスト表示ボタン１３０３押下により図示はしていないがデータベースに保存されているデータの装置ＩＤのリストが表示されるので、リストから選択入力する。
学習期間入力ウィンドウ１３０４には、学習データを抽出したい期間の開始日と終了日を入力する。センサ選択ウィンドウ１３０５には、使用するセンサを入力する。リスト表示ボタン１３０６のクリックによりセンサリスト１３０７が表示されるので、リストから選択入力する。リストから複数選択することも可能である。基準算出パラメータ入力ウィンドウ１３０８には、正常モデル作成において指定するパラメータを入力する。図は正常モデルとして局所部分空間を採用した場合の例であり、モデル作成に使う近傍ベクトル数と正則化パラメータを入力する。正則化パラメータは、数２において相関行列Cの逆行列が求められないことを防ぐため、対角成分に加算する小さい数である。
しきい値設定パラメータ入力ウィンドウ１３０９には、図５に示す処理における交差検証のグループをどのように決めるかをラジオボタンで選択する。１日を１グループとするか、指定された数のグループに分割するかを選び、後者の場合はグループ数を入力する。また、ステップＳ５１１のしきい値設定処理において累積ヒストグラムに適用する比率を値入力する。レシピ名入力ウィンドウ１３１０には、入力された情報に対応付けるユニークな名前を入力する。全ての情報を入力したらテスト期間入力ウィンドウ１３１１にテスト対象期間を入力し、テストボタン１３１２の押下により、レシピのテストを行う。
この操作により、同じレシピ名で実行したテストの通し番号が採番される。次に、特徴ベクトル抽出部１０４において、指定した学習期間のセンサ信号１０２から特徴ベクトルの抽出を行う。図５で説明したステップＳ５０２の正準化においては、指定した学習期間の全データを用いて平均と標準偏差を求める。この平均と標準偏差の値もセンサ毎にレシピ名およびテスト番号に対応付けて保存しておく。
図５で説明したステップＳ５０７に対応して、交差検証により基準ベクトル算出部１０５において各時刻の基準ベクトルを作成し、ステップＳ５０８に対応して、異常測度算出部１０６において異常測度を算出する。更に、ステップＳ５１１に対応して、しきい値算出部１０７において異常判定しきい値を算出し、レシピ名およびテスト番号と対応付けて保存する。分布密度算出部１０９において、図６に示した処理により２次元の特徴の分布密度画像を作成し、保存する。併せて、ステップＳ６０５で算出した処理範囲を全ての特徴について保存しておく。その他、装置ＩＤ情報、使用センサ情報、学習期間、特徴ベクトル抽出に用いるパラメータ、基準作成パラメータをレシピ名およびテスト番号と対応付けて保存する。
次に、指定したテスト期間のセンサ信号１０２を入力として、図７に示した異常検知の処理を行う。検出された異常について、図８乃至図１１の何れかに示した手順で関連センサ特定の処理を行う。検出された異常は区間毎に通し番号が付加される。区間毎に、通し番号、区間開始時刻、区間終了時刻、複数の関連センサ名と観測ベクトルの対応点がプロットされた分布密度画像名が記録される。
テストの結果をユーザに示すためのＧＵＩの例を図１４Ａおよび図１４Ｂおよび図１４Ｃに示す。各画面の上部に表示されたタブを選択することにより、結果全体表示画面１４０１と結果拡大表示画面１４０２と異常診断結果表示画面１４０３を切り換えることができる。
図１４Ａには、結果全体表示画面１４０１を示す。結果全体表示画面１４０１には、指定された全期間の異常測度、しきい値、判定結果とセンサ信号の時系列グラフを表示する。期間表示ウィンドウ１４０４には、指定された学習期間およびテスト期間が表示される。異常測度表示ウィンドウ１４０５には、指定された学習期間およびテスト期間の異常測度としきい値と判定結果が表示される。センサ信号表示ウィンドウ１４０６には、指定された期間の指定されたセンサの出力値が表示される。センサの指定は、センサ名選択ウィンドウ１４０７への入力によって行う。ただし、ユーザが指定する前は、先頭のセンサが選択されている。カーソル１４０８は、拡大表示の時の起点を表し、マウス操作により移動できる。
表示日数指定ウィンドウ１４０９には、この画面では使用しないが、結果拡大表示画面１４０２での、拡大表示の起点から終点までの日数が表示される。この画面で入力することもできる。日付表示ウィンドウ１４１０には、カーソル位置の日付が、表示される。終了ボタン１４１１押下により結果全体表示画面１４０１、結果拡大表示画面１４０２、異常診断結果表示画面１４０３とも消去し終了する。
図１４Ｂには、結果拡大表示画面１４０２を示す。結果拡大表示画面１４０２には、結果全体表示画面１４０１において、カーソル１４０８で示された日付を起点として、指定された日数の異常測度、しきい値、判定結果とセンサ信号の時系列グラフを表示する。期間表示ウィンドウ１４０４には、結果全体表示画面１４０１と同じ情報が表示される。異常測度表示ウィンドウ１４０５およびセンサ信号表示ウィンドウ１４０６には、結果全体表示画面１４０１と同様の情報が、拡大表示される。表示日数指定ウィンドウ１４０９で、拡大表示の起点から終点までの日数を指定する。
日付表示ウィンドウ１４１０には、拡大表示の起点の日付が表示されている。スクロールバー１４１２で表示の起点を変更することも可能であり、この変更はカーソル１４０８の位置と日付表示ウィンドウ１４１０の表示に反映される。スクロールバー表示領域１４１３の全体の長さは結果全体表示画面１４０１に表示されている全期間に相当する。また、スクロールバー１４１２の長さは表示日数指定ウィンドウ１４０９で指定された日数に相当し、スクロールバー１４１２の左端部が拡大表示の起点に対応する。異常が検出された時刻には、異常測度グラフの対応する位置に異常区間番号１４１４がバルーン表示される。終了ボタン１４１１押下により終了する。
図１４Ｃには、異常診断結果表示画面１４０３の例を示す。異常に関連するセンサの２次元の特徴分布密度と時系列グラフを併せて表示する。異常診断結果表示画面１４０３に切り替ったときには、結果拡大表示画面１４０２に表示されている異常のうち最も番号の若い異常を表示対象とする。日付表示ウィンドウ１４１５には、表示対象の異常の日付が表示される。異常区間番号表示ウィンドウ１４１６には、表示対象の異常の区間番号が表示される。ウィンドウ１４１６では数値入力により、異常区間番号を指定することができる。時刻表示ウィンドウ１４１７には、表示対象の異常の区間開始時刻と区間終了時刻が表示される。
分布密度表示ウィンドウ１４１８には、例えば図１１のステップＳ１１０９において保存された画像を表示する。この画像は、ステップＳ１１１０において特定された２個のセンサの２次元頻度分布を表したものであり、横軸は”SensorX”の規格化された値、縦軸は”SensorY” の規格化された値を表す。１４１９は、学習期間の分布を表しており、濃く表示されているほど高頻度であることを意味する。１４２０は、表示対象の異常区間の観測ベクトルに対応するデータである。異常表示指示ウィンドウ１４２２で、チェックボックスのチェックをはずすことにより、表示しないようにすることもでき、その特徴値の学習期間のデータの分布密度を確認できる。処理範囲表示ウィンドウ１４２１には２個のセンサのステップＳ６０５において算出された処理範囲が表示される。これらの値は、それぞれ分布密度画像の左端に対応する”SensorX”の値、右端に対応する”SensorX”の値、下端に対応する”SensorY”の値、上端に対応する”SensorY”の値である。
異常測度表示ウィンドウ１４２３には、表示対象の異常区間を含む期間の異常測度としきい値と判定結果が表示される。表示期間は信号の変化が観測可能なよう十分に拡大される長さを予め決めておく。第一関連センサ信号表示ウィンドウ１４２４には、ステップＳ１１１０において特定された２個のセンサのうち一方の、異常測度表示ウィンドウ１４２３と同じ期間の出力値がグラフ表示される。第二関連センサ信号表示ウィンドウ１４２５には、もう一方のセンサ信号出力値がグラフ表示される。基準表示指示ウィンドウ１４２６のチェックボックスにチェックを入れることにより、図示していないが第一関連センサ信号表示ウィンドウ１４２４および第二関連センサ信号表示ウィンドウ１４２５に基準値を異なる色で重ねて表示することができる。
ボタン１４２７及び１４２８は、表示する異常診断結果を切り替えるためのボタンで、ボタン１４２７「次」をクリックすると、現在表示されている異常のデータよりも異常診断結果表示画面１４０３において次に若い番号の異常のデータが表示される。また、ぼたん１４２８「前」をクリックすると、現在表示されている異常のデータよりも異常診断結果表示画面１４０３において１つ若い番号の異常のデータが表示される。
又、図示はしていないが、分布密度表示ウィンドウ１４１８には、ステップＳ１１０９において保存された表数の画像を並べて表示するようにしてもよい。
本図では２個のセンサの時系列グラフを別々に表示したが、１個のウィンドウに異なる色で重ねて表示してもよい。終了ボタン１４１１押下により終了する。
図１４Ａ〜Ｃに示すいずれかの画面で、終了ボタン１４１１押下により異常検出結果および診断結果の確認が終了したら、図１３に示すレシピ設定画面１３０１の表示に戻る。テスト番号表示ウィンドウ１３１３には、上記のテストで採番された番号が表示されている。確認した内容に問題があれば、学習期間や選択センサ、パラメータなどを変更し、テストボタン１３１２の押下により、再度テストを行う。あるいは、一度行ったテストの結果を再度確認することもできる。テスト番号表示ウィンドウ１３１３からテスト番号を選択入力し、表示ボタン１３１４を押下する。この操作により、レシピ名とテスト番号に対応付けて保存された情報をロードし結果全体表示画面１４０１を表示する。タブの切り替えにより結果拡大表示画面１４０２または異常診断結果表示画面１４０３を表示させることもできる。確認が済んだら終了ボタン１４１１押下により、レシピ設定画面１３０１の表示に戻る。
登録ボタン１３１５の押下により、上記レシピ名とテスト番号表示ウィンドウ１３１３に表示中のテスト番号に対応付けて保存されている情報をレシピ名と対応付けて登録し、終了する。キャンセルボタン１３１６が押下された場合は、何も保存しないで終了する。
また、テスト結果一覧ボタン１３１７が押下された場合は、図１５に示す、テスト結果一覧表示画面１５０１を表示する。テスト結果リスト１５０２には、全てのテストの学習期間、テスト期間、選択センサ番号、基準作成パラメータ、しきい値設定パラメータなどのレシピ情報と、しきい値、異常区間数などのテスト結果情報を表示する。リストの左端に選択チェックボタンがあり、いずれか一つのみ選択することができる。詳細表示ボタン１５０３押下により、レシピ名とテスト番号に対応付けて保存された情報をロードし、結果全体表示画面１４０１を表示する。タブの切り替えにより結果拡大表示画面１４０２または異常診断結果表示画面１４０３を表示させることもできる。確認が済んだら終了ボタン１４１１押下により、テスト結果一覧表示画面１５０１の表示に戻る。登録ボタン１５０４の押下により、選択中のテスト番号に対応付けて保存されている情報をレシピ名と対応付けて登録し、テスト結果一覧表示画面１５０１の表示およびレシピ設定画面１３０１の表示を終了する。戻るボタン１５０５が押下された場合は、レシピの登録は行わずにレシピ設定画面１３０１の表示に戻る。
登録されたレシピは、活性か不活性かのラベルをつけて管理され、新しく観測されたデータに対しては、装置ＩＤが一致する活性なレシピの情報を用いて図７ないし図１１を用いて説明した何れかの特徴ベクトル抽出から異常検出、関連センサ特定までの処理を行い、結果をレシピ名と対応付けて保存しておく。
以上の異常検診断処理の結果をユーザに示すためのＧＵＩの例を、図１６に示す。
図１６は、表示対象を指定するＧＵＩの例である。表示対象指定画面１６０１から表示対象の設備、レシピおよび期間を指定する。初めに、装置ＩＤ選択ウィンドウ１６０２により装置ＩＤを選択する。次に、レシピ名選択ウィンドウ１６０３により、装置ＩＤを対象としたレシピのリストから表示対象のレシピを選択する。データ記録期間表示部１６０４には、入力されたレシピを用いて処理され、記録が残されている期間の開始日と終了日が表示される。結果表示期間指定ウィンドウ１６０５には、結果を表示したい期間の開始日と終了日を入力する。表示センサ指定ウィンドウ１６０６には、表示したいセンサの名を入力する。表示ボタン１６０７押下により図１４Ａに示す結果全体表示画面１４０１を表示する。終了ボタン１６０８押下により終了する。
結果表示にかかわるＧＵＩの画面および操作は、図１４Ａ〜Ｃに示すテスト結果表示にかかわるＧＵＩとほぼ同じであるため、異なる部分のみ説明する。結果全体表示画面１４０１、結果拡大表示画面１４０２における期間表示ウィンドウ１４０４に対応する表示ウィンドウとしては、図１７に示すように結果表示期間指定ウィンドウ１６０５で指定された表示期間が表示される。結果全体表示画面１４０１に対応する画面において、異常測度表示ウィンドウ１４０５に対応するウィンドウには、指定された表示期間の異常測度としきい値と判定結果が表示される。センサ信号表示ウィンドウ１４０６に対応するウィンドウには、指定された期間の、表示センサ指定ウィンドウ１６０６により指定されたセンサの出力値が表示される。表示対象センサは、センサ名選択ウィンドウ１４０８に対応するウィンドウへの入力によって変更することも可能である。
上記実施例は学習データ設定をオフライン、異常診断処理をリアルタイム、結果表示をオフラインでそれぞれ処理するものであるが、結果表示もリアルタイムに行うことが可能である。その場合、表示期間の長さ、表示対象とするレシピ、表示対象とする情報を予め定めておき、一定時間毎に最新の情報を表示するよう構成すればよい。
逆に、任意の期間を設定し、レシピを選択して、オフラインで異常診断処理を行う機能を付加したものも本発明の範囲に含まれる。
本実施例によれば、状態が複雑に偏けする設備に対しても精度の高い異常測度を算出することができる。又、誤報の発生を抑制することができる。
更に、本実施例によれば、正常なセンサの信号から得られるデータが異常なセンサの信号のデータの影響を受けるような場合であっても、検出した異常に関連するセンサを正しく診断することが可能になった。
また、検査の結果得られる２次元特徴分布密度を画像で表示するので、異常に関連するセンサの見分けが容易になる。さらに、学習データの分布の傾向を画面上で確認することができるようになった。
以上、設備から出力されるセンサ信号に基づき異常診断する方法の実施例を説明したが、別の実施例として、さらに、設備から出力されるイベント信号も利用して異常診断する方法を説明する。
図１８に、本発明の異常診断方法を実現するシステムの第二の実施例を示す。本実施例における構成は、実施例１で説明した図１に示した構成に、モード分割部１８０２を加えたものとなっている。実施例１と同じ構成については同じ番号を付してある。モード分割部１８０２は、設備からイベント信号１８０１を入力し、これに基づき設備の稼動状態を表すモードに分割する。モード分割の結果は、しきい値算出部１０７´に入力され、図５で説明した処理フローが分割されたモードごとに実行されてステップＳ５１１のしきい値設定がモード別に実施される。
また、図７で説明した処理フローも分割されたモードごとに実行されて、異常検出部１０８´におけるステップＳ７０７の異常判定を、対応するモードのしきい値を用いて行う。また、分布密度算出部１０９´でもモード分割の結果を入力し、分布密度画像作成をモード毎に行い、関連センサ特定部１１０´におけるセンサ特定の処理をモード毎に行う。
上記以外は、全て前述の実施例１で説明した方法と同様であるため、実施例１と異なるイベント信号に基づくモード分割方法の実施例を図１９Ａ〜Ｃを用いて説明する。イベント信号の例を図１９Ａに示す。図１９Ａに示したリスト１９２０は、不定期に出力される設備の操作・故障・警告を表す信号であり、時刻と操作・故障・警告を表す文字列またはコードからなる。
図１９Ｂに示すように、このイベント信号をモード分割部１８０２に入力し（Ｓ１９０１）、所定の文字列またはコードの検索により起動シーケンスと停止シーケンスの切り出しを行う（Ｓ１９０２）。その結果をもとに、停止シーケンスの終了時刻から起動シーケンスの開始時刻までの「定常OFF」モード１９１１、起動シーケンス中の「起動」モード１９１２、起動シーケンスの終了時刻から停止シーケンスの開始時刻までの「定常ON」モード１９１３、停止シーケンス中の「停止」モード１９１４の４つの稼動状態に分割する（Ｓ１９０３）。図１９Ｃに例を示す。
シーケンス切り出しのためには、予めシーケンスの開始イベントおよび終了イベントを指定しておき、イベント信号の先頭から最後まで以下の要領でスキャンしながら切り出していく。
(1)シーケンスの途中でない場合は、開始イベントを探索する。見つかったらシーケンスの開始とする。
(2)シーケンスの途中の場合は、終了イベントを探索する。見つかったらシーケンスの終了とする。ここで終了イベントとは、指定の終了イベントのほか、故障、警告、指定の開始イベントとする。
以上のように、イベント信号を利用することにより、多様な稼動状態を正確に分けることができ、モード別にしきい値を設定することにより、「起動」モード１９１２および「停止」モード１９１４の過渡期において学習データ不足により感度を落とす必要がある場合でも、「定常OFF」モード１９１１および「定常ON」モード１９１３では高感度な異常検知が可能になる。
即ち、本実施例によれば、設備の多様な稼働状態に応じてしきい値を設定することができ、高感度な異常検出を可能にし、更に検出した異常に関連するセンサを正しく診断することが可能になった。また、学習データの分布密度をモード毎に算出することにより、異なる状態のデータが含まれなくなるため、診断が容易になり、分布密度画像の理解も容易になる。
１０１・・・設備 １０２・・・センサ信号 １０３・・・センサ信号蓄積部 １０４・・・特徴ベクトル抽出部 １０５・・・基準ベクトル作成部 １０６・・・異常測度算出部 １０７・・・しきい値算出部 １０８・・・異常検出部 １０９・・・分布密度算出部 １１０・・・関連センサ特定部 １３０１・・・レシピ設定画面 １４０１・・・結果全体表示画面 １４０２・・・結果拡大表示画面 １４０３・・・異常診断結果表示画面 １５０１・・・テスト結果一覧表示画面 １６０１・・・表示対象指定画面 １８０１・・・イベント信号 １８０２・・・モード分割部。

Claims (11)

  1. 設備または装置に装着された複数のセンサから出力される複数のセンサ信号に基づいて
    前記設備または装置の異常を診断する方法であって、予め指定された学習期間において、前記複数のセンサ信号をそれぞれ正準化して得た特徴値から多次元の特徴ベクトルを求め、前記多次元の特徴ベクトルから注目ベクトルを選び、前記注目ベクトルの近傍ベクトルを用いて基準ベクトルを作成し、前記多次元の特徴ベクトルの前記基準ベクトルまでの距離に基づいて異常測度を算出し、前記算出した異常測度に基づいてしきい値を設定し、前記多次元の特徴ベクトルの前記特徴値から算出した値の頻度分布に基づいて分布密度画像を生成し、前記設備または装置の異常を診断する時において、前記複数のセンサ信号から時刻毎に多次元特徴ベクトルを抽出し、前記求めた予め指定された学習期間の前記多次元の特徴ベクトルの中から前記時刻毎に抽出した多次元の特徴ベクトルに近い特徴ベクトルを選択し、前記選択した特徴ベクトルを用いて基準特徴ベクトルを算出し、前記抽出した時刻毎の多次元特徴ベクトルの前記基準特徴ベクトルまでの距離に基づき異常測度を算出し、前記算出した異常測度と前記設定したしきい値とを比較することにより異常を検出し、前記複数のセンサ信号と前記学習期間において求めた特徴ベクトルを用いて生成した２次元の分布密度に基づいて前記検出した異常に関連するセンサを前記複数のセンサの中から特定することを特徴とする異常診断方法。
  2. 設備または装置に装着された複数のセンサから出力される複数のセンサ信号をもとに学
    習データを作成して蓄積する工程と、前記設備または装置に装着された複数のセンサから
    新たに出力された複数のセンサ信号の異常診断を行う工程とを含む設備の異常を診断する
    方法であって、
    前記学習データを作成して蓄積する工程は、前記複数のセンサ信号から特徴ベクトルを
    抽出する工程と、前記抽出した特徴ベクトルを学習データとして蓄積する工程と、前記学
    習データとして蓄積した特徴ベクトルの各々について交差検証用のグループ番号付けをし、前記蓄積した学習データの中から注目ベクトルを選択する工程と、前記蓄積した学習データの中の前記注目ベクトルと異なる前記グループ番号付けをした特徴ベクトルの中から前記選択された注目ベクトルの近傍ベクトルを選択してこの選択した近傍ベクトルを用いて学習用の基準ベクトルを作成する工程と、前記作成した学習用の基準ベクトルまでの前記学習データとして蓄積した特徴ベクトルの距離に基づいて異常測度を算出する工程と、前記算出した前記異常測度に基づきしきい値を算出する工程と、前記特徴ベクトルから全組合せの２次元の特徴分布密度を算出する工程とを有し、
    前記センサ信号の異常診断を行う工程は、前記設備または装置に装着された前記複数の
    センサから出力される複数のセンサ信号から特徴ベクトルを抽出して観測ベクトルとする
    工程と、前記抽出した観測ベクトルに応じて前記蓄積した学習データの中から所定数の特
    徴ベクトルを選択する工程と、前記学習データの中から選択された所定数の特徴ベクトル
    を用いて異常診断用の基準ベクトルを作成する工程と、前記観測ベクトルの前記作成した異常診断用の基準ベクトルまでの距離に基づいて前記観測ベクトルの異常測度を算出する工程と、前記算出した異常測度と前記学習データを作成して蓄積する工程において算出したしきい値に基づき前記観測ベクトルが異常か正常かを判定する工程と、前記判定する工程において異常と判定された観測ベクトルと前記２次元の特徴分布密度に基づき異常に関連するセンサを特定する工程とを有する
    ことを特徴とする異常診断方法。
  3. 前記学習用の基準ベクトルを作成する工程と、前記異常診断用の基準ベクトルを作成す
    る工程とは、何れも局所部分空間法を利用すること特徴とする請求項２に記載の異常診断方法。
  4. 前記センサ信号の異常診断を行う工程における前記異常に関連するセンサを特定する工
    程は、前記観測ベクトルと前記基準ベクトルに基づき観測値と基準値の差を算出する工程
    と、前記算出された差が最大となる特徴を構成するセンサを第一の関連センサとする工程
    と、前記第一の関連センサにかかわる前記２次元の特徴分布密度に基づき、前記観測ベク
    トルに対応する値の密度が最小となる前記第一の関連センサのペアとなる特徴を構成する
    センサを第二の関連センサとする工程からなる請求項２に記載の異常診断方法。
  5. 前記センサ信号の異常診断を行う工程における前記異常に関連するセンサを特定する工
    程は、全ての前記２次元の特徴分布密度に基づき、前記観測ベクトルに対応する値の密度
    を算出する工程と、前記特徴毎に前記算出した密度を合計した値が最小となる特徴を構成
    するセンサを第一の関連センサとする工程と、前記第一の関連センサにかかわる前記２次
    元の特徴分布密度に基づき、前記観測ベクトルに対応する値の密度が最小となる前記第一
    の関連センサのペアとなる特徴を構成するセンサを第二の関連センサとする工程からなる
    請求項２に記載の異常診断方法。
  6. 前記センサ信号の異常診断を行う工程における前記異常に関連するセンサを特定する工
    程は、前記学習データのうち前記第一の関連センサのセンサ信号が異常と判定された観測
    ベクトルと近いデータを用いて前記第一の関連センサを除く全組合せの２次元の部分的特
    徴分布密度を算出する工程と、全ての前記２次元の部分的特徴分布密度に基づき、前記観
    測ベクトルに対応する値の密度を算出する工程と、特徴毎に前記密度を合計した値が最小
    となる特徴を構成するセンサを第二の関連センサとする工程と、前記第二の関連センサに
    かかわる前記２次元の部分的特徴分布密度に基づき、前記観測ベクトルに対応する値の密
    度が最小となる前記第二の関連センサのペアとなる特徴を構成するセンサを第三の関連セ
    ンサとする工程とをさらに含む請求項４または請求項５に記載の異常診断方法。
  7. 前記センサ信号の異常診断を行う工程における前記異常に関連するセンサを特定する工
    程は、連続して異常検出された区間を特定する工程と、前記区間に含まれる全ての観測ベ
    クトルについて、前記観測ベクトルと前記基準ベクトルに基づき観測値と基準値の差の絶
    対値の和を算出する工程と、前記算出された和が最大となる特徴を構成するセンサを第一
    の関連センサとする工程と、前記区間に含まれる全ての観測ベクトルについて、前記第一
    の関連センサにかかわる前記２次元の特徴分布密度に基づき、前記観測ベクトルに対応す
    る値の密度の和および密度が非０となる回数を算出する工程と、前記算出された回数が最
    小かつそのうち前記算出された和が最小となる前記第一の関連センサのペアとなる特徴を
    構成するセンサを第二の関連センサとする工程からなる請求項２に記載の異常診断方法。
  8. 前記異常測度および前記しきい値および前記判定結果の時系列グラフと、前記特定され
    た異常に関連するセンサから出力されたセンサ信号の時系列グラフと、前記特定された２
    次元の特徴の学習データの分布密度を表す画像に、前記異常と判定された観測ベクトルに
    対応する点を重ねてプロットした画像を表示する工程をさらに含むことを特徴とする請求
    項２に記載の異常診断方法。
  9. 設備または装置に装着された複数のセンサから出力される複数のセンサ信号をもとに学
    習データを作成して蓄積する工程と、前記設備または装置に装着された複数のセンサから
    新たに出力された複数のセンサ信号の異常診断を行う工程とを含む設備の異常を診断する
    方法であって、
    前記学習データを作成して蓄積する工程は、設備または装置から出力されるイベント信
    号に基づいて稼動状態別のモード分割を行う工程と、前記複数のセンサ信号から特徴ベク
    トルを抽出する工程と、前記抽出した特徴ベクトルを学習データとして蓄積する工程と、
    前記学習データとして蓄積した特徴ベクトルのおのおのについて交差検証用のグループ番号付けをして前記特徴ベクトルに応じて前記蓄積した学習データの中から特徴ベクトルを選択する工程と、前記蓄積した学習データの中の前記選択された特徴ベクトルと異なる前記グループ番号付けをした特徴ベクトルの中から前記選択された特徴ベクトルの近傍ベクトルを選択してこの選択した近傍ベクトルを用いて学習用の基準ベクトルを作成する工程と、前記作成した学習用の基準ベクトルまでの前記学習データとして蓄積した特徴ベクトルの距離に基づいて異常測度を算出する工程と、前記算出した異常測度に基づき前記モード毎にしきい値を算出する工程と、前記学習データとして蓄積した特徴ベクトルから全組合せの２次元の特徴分布密度を前記モード毎に算出する工程を有し、
    前記センサ信号の異常診断を行う工程は、前記設備または装置に装着された前記複数の
    センサから出力される複数のセンサ信号から特徴ベクトルを抽出して観測ベクトルとする
    工程と、
    前記抽出した観測ベクトルに応じて前記蓄積した学習データの中から注目ベクトルを選択する工程と、前記学習データの中から選択された特徴ベクトルの近傍ベクトルを用いて異常診断用の基準ベクトルを作成する工程と、前記観測ベクトルの前記作成した異常診断用の基準ベクトルまでの距離に基づいて前記観測ベクトルの異常測度を算出する工程と、前記算出した異常測度と前記分割したモード別に算出した前記しきい値に基づき前記観測ベクトルが異常か正常かを判定する工程と、前記判定する工程において異常と判定された観測ベクトルに対応するセンサ信号と前記モード別に算出した前記２次元の特徴分布密度に基づき異常に関連するセンサを特定する工程とを有する
    ことを特徴とする異常診断方法。
  10. 設備または装置に装着されたセンサから出力されるセンサ信号に基づいて前記設備また
    は装置の異常を診断する装置であって、前記センサ信号を蓄積する生データ蓄積手段と、
    前記センサ信号から特徴ベクトルを抽出する特徴ベクトル抽出手段と、予め指定された学習期間において前記センサ信号から抽出した前記特徴ベクトルから注目ベクトルを選んで前記選んだ注目ベクトルの近傍ベクトルを用いて基準ベクトルを算出して前記設備または装置の異常を診断する期間においては前記診断している各時刻において前記センサ信号から抽出した前記特徴ベクトルに 近い特徴ベクトルを選択し、前記選択した特徴ベクトルを用いて前記各時刻における基準特徴ベクトルを算出する基準ベクトル算出手段と、前記予め指定された学習期間においては前記多次元の特徴ベクトルの前記基準ベクトルまでの距離に基づいて異常測度を算出し前記設備または装置の異常を診断する期間においては前記各時刻の前記特徴ベクトルの前記基準特徴ベクトルまでの距離に基づき異常測度を算出する異常測度算出手段と、前記予め指定された学習期間における前記異常測度に基づいてしきい値を算出するしきい値算出手段と、前記設備または装置の異常を診断する期間において前記異常測度算出手段で算出した前記異常測度と前記しきい値算出手段で算出した前記しきい値との比較により異常を検出する異常検出手段と、前記予め指定された学習期間における前記センサ信号の２次元の分布密度を算出する分布密度算出手段と、前記異常検出手段により検出された異常に対応するセンサ信号と前記分布密度算出手段で算出した２次元の分布密度に基づいて異常に関連するセンサを特定する関連センサ特定手段とを備えたことを特徴とする異常診断装置。
  11. 設備または装置に装着されたセンサから出力されるセンサ信号およびイベント信号に基
    づいて前記設備または装置の異常を診断する装置であって、前記センサ信号を蓄積する生
    データ蓄積手段と、前記イベント信号に基づいて稼動状態別のモード分割を行うモード分
    割手段と、前記センサ信号から特徴ベクトルを抽出する特徴ベクトル抽出手段と、予め指定された学習期間において前記センサ信号から抽出した前記特徴ベクトル から注目ベクトルを選んで前記選んだ注目ベクトルの近傍ベクトルを用いて基準ベクトルを算出して前記設備または装置の異常を診断する期間においては前記診断している各時刻において前記センサ信号から抽出した前記特徴ベクトルに 近い特徴ベクトルを選択し、前記選択した特徴ベクトルを用いて前記各時刻における基準特徴ベクトルを算出する基準ベクトル算出手段と、前記予め指定された学習期間においては前記多次元の特徴ベクトルの前記基準ベクトルまでの距離に基づいて異常測度を算出し前記設備または装置の異常を診断する期間においては前記各時刻の前記特徴ベクトルの前記基準特徴ベクトルまでの距離に基づき異常測度を算出する異常測度算出手段と、前記予め指定された学習期間における前記異常測度に基づいて前記モード分割手段で分割したモード毎にしきい値を算出するしきい値算出手段と、前記設備または装置の異常を診断する期間において前記異常測度算出手段で算出した前記異常測度と前記しきい値算出手段で算出したモード毎のしきい値との比較により異常を検出する異常検出手段と、前記予め指定された学習期間における前記センサ信号の２次元の分布密度を前記モード毎に算出する分布密度算出手段と、前記異常検出手段により検出された異常に対応するセンサ信号のうち前記モードに対応する前記センサ信号と前記分布密度算出手段で算出した２次元の分布密度に基づいて異常に関連するセンサを特定する関連センサ特定手段とを備えたことを特徴とする異常診断装置。
JP2013009316A 2013-01-22 2013-01-22 異常診断方法およびその装置 Active JP6076751B2 (ja)

Priority Applications (3)

Application Number Priority Date Filing Date Title
JP2013009316A JP6076751B2 (ja) 2013-01-22 2013-01-22 異常診断方法およびその装置
US14/761,748 US9779495B2 (en) 2013-01-22 2014-01-15 Anomaly diagnosis method and apparatus
PCT/JP2014/050547 WO2014115615A1 (ja) 2013-01-22 2014-01-15 異常診断方法およびその装置

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
JP2013009316A JP6076751B2 (ja) 2013-01-22 2013-01-22 異常診断方法およびその装置

Publications (2)

Publication Number Publication Date
JP2014142697A JP2014142697A (ja) 2014-08-07
JP6076751B2 true JP6076751B2 (ja) 2017-02-08

Family

ID=51227407

Family Applications (1)

Application Number Title Priority Date Filing Date
JP2013009316A Active JP6076751B2 (ja) 2013-01-22 2013-01-22 異常診断方法およびその装置

Country Status (3)

Country Link
US (1) US9779495B2 (ja)
JP (1) JP6076751B2 (ja)
WO (1) WO2014115615A1 (ja)

Cited By (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
EP3462267A1 (en) 2017-09-27 2019-04-03 Panasonic Intellectual Property Corporation of America Anomaly diagnosis method and anomaly diagnosis apparatus

Families Citing this family (62)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
FR3015757B1 (fr) * 2013-12-23 2019-05-31 Electricite De France Procede d'estimation quantitative du colmatage des plaques d'un generateur de vapeur
JP6223936B2 (ja) * 2014-09-12 2017-11-01 株式会社日立ハイテクノロジーズ 異常傾向検出方法およびシステム
FR3032273B1 (fr) * 2015-01-30 2019-06-21 Safran Aircraft Engines Procede, systeme et programme d'ordinateur pour phase d'apprentissage d'une analyse acoustique ou vibratoire d'une machine
JP5849167B1 (ja) * 2015-04-09 2016-01-27 株式会社日立パワーソリューションズ 異常検知方法およびその装置
KR101647826B1 (ko) * 2015-05-27 2016-08-11 한국에너지기술연구원 전력시스템의 안정도를 판별하는 장치 및 그 전력시스템의 안정도를 관리하는 방법
WO2017058688A1 (en) * 2015-09-30 2017-04-06 Schlumberger Technology Corporation Downhole tool analysis using anomaly detection of measurement data
KR102384742B1 (ko) * 2015-10-27 2022-04-07 삼성에스디에스 주식회사 센서 데이터를 이용한 이상 감지 장치 및 방법
US10073421B2 (en) * 2015-11-17 2018-09-11 Rockwell Automation Technologies, Inc. Predictive monitoring and diagnostics systems and methods
WO2017150286A1 (ja) * 2016-02-29 2017-09-08 日本電気株式会社 システム分析装置、システム分析方法、及び、コンピュータ読み取り可能な記録媒体
GB2554038B (en) * 2016-05-04 2019-05-22 Interactive Coventry Ltd A method for monitoring the operational state of a system
EP3258333A1 (en) * 2016-06-17 2017-12-20 Siemens Aktiengesellschaft Method and system for monitoring sensor data of rotating equipment
US10565513B2 (en) * 2016-09-19 2020-02-18 Applied Materials, Inc. Time-series fault detection, fault classification, and transition analysis using a K-nearest-neighbor and logistic regression approach
US10677943B2 (en) * 2016-12-16 2020-06-09 Smiths Detection, Llc System and method for monitoring a computed tomography imaging system
KR101971553B1 (ko) * 2017-03-21 2019-04-23 (주)심플랫폼 사물인터넷 기반 기기 관리 시스템 및 방법
CN110678820B (zh) * 2017-05-26 2022-10-18 三菱电机楼宇解决方案株式会社 异常重要度计算***及异常重要度计算装置
JP6772963B2 (ja) * 2017-06-05 2020-10-21 トヨタ自動車株式会社 異常診断装置及び異常診断方法
WO2019012653A1 (ja) 2017-07-13 2019-01-17 日本電気株式会社 学習システム、分析システム、学習方法及び記憶媒体
JP7017363B2 (ja) * 2017-10-06 2022-02-08 株式会社日立パワーソリューションズ 異常検知装置および異常検知方法
JP6693938B2 (ja) * 2017-11-17 2020-05-13 ファナック株式会社 外観検査装置
WO2019160068A1 (ja) * 2018-02-16 2019-08-22 日本電気株式会社 異常音判定基準作成装置および異常音検知装置
US11972178B2 (en) * 2018-02-27 2024-04-30 Falkonry Inc. System and method for explanation of condition predictions in complex systems
KR101967339B1 (ko) * 2018-03-06 2019-04-09 단국대학교 산학협력단 심층학습 기반의 ａｄａｓ 센서 고장 진단 및 백업을 위한 장치 및 방법
JP7013993B2 (ja) * 2018-03-26 2022-02-01 トヨタ自動車株式会社 診断装置及び診断方法
JP7043320B2 (ja) * 2018-03-30 2022-03-29 株式会社小松製作所 状態分析装置および状態分析方法
US10553046B2 (en) 2018-04-05 2020-02-04 GM Global Technology Operations LLC Vehicle prognostics and remedial response
CN109002833B (zh) * 2018-06-12 2019-08-27 国家卫生健康委科学技术研究所 一种微液滴数据分析方法及***
WO2019244203A1 (ja) 2018-06-18 2019-12-26 三菱電機株式会社 診断装置、診断方法及びプログラム
CN109005513B (zh) * 2018-06-26 2021-03-19 北京酷云互动科技有限公司 手机终端关联方法和手机终端关联***
JP7032261B2 (ja) * 2018-07-24 2022-03-08 日立グローバルライフソリューションズ株式会社 異常検知システムおよび異常検知方法
US11122065B2 (en) * 2018-08-14 2021-09-14 Vmware, Inc. Adaptive anomaly detection for computer systems
CN112567306A (zh) * 2018-08-31 2021-03-26 东芝三菱电机产业***株式会社 制造过程监视装置
KR101971278B1 (ko) * 2018-12-13 2019-04-26 주식회사 알고리고 인공신경망을 이용한 비정상 데이터 구분 장치
EP3871120A1 (en) * 2018-12-17 2021-09-01 Huawei Technologies Co., Ltd. Apparatus and method for detecting an anomaly among successive events and computer program product therefor
CN109507990B (zh) * 2018-12-25 2021-06-15 中南大学 一种故障溯源方法及***
JP7463055B2 (ja) * 2018-12-27 2024-04-08 株式会社Ihi 異常診断装置、異常診断方法、異常診断プログラム、及び記録媒体
US10922906B2 (en) 2019-03-28 2021-02-16 GM Global Technology Operations LLC Monitoring and diagnosing vehicle system problems using machine learning classifiers
WO2020202567A1 (ja) * 2019-04-05 2020-10-08 株式会社Ｉｈｉ原動機 振動音響解析方法及び装置と機器異常部位推定方法及び装置
KR102119891B1 (ko) * 2019-04-11 2020-06-05 주식회사 알고리고 미세 변화 데이터 및 공간 데이터를 기초로 한 인공신경망을 이용한 비정상 데이터 구분 장치
KR102198224B1 (ko) * 2019-04-11 2021-01-05 주식회사 알고리고 공간 데이터를 기초로 한 인공신경망을 이용한 비정상 데이터 구분 장치
JP7218867B2 (ja) * 2019-05-24 2023-02-07 Ｋｙｂ株式会社 異常検知装置及び異常検知方法
KR20200141812A (ko) 2019-06-11 2020-12-21 삼성전자주식회사 뉴럴 네트워크를 이용하여 이상 신호를 감지하는 방법 및 장치
JP7344015B2 (ja) * 2019-06-13 2023-09-13 株式会社日立ハイテクソリューションズ 異常検知装置及び異常検知方法
DE102019209797A1 (de) * 2019-07-03 2021-01-07 Thyssenkrupp Ag Verfahren und Einrichtung zur Ermittlung des fahrbetriebsbedingten Zustandes von Karosseriekomponenten eines Fahrzeugs sowie eines entsprechenden Fahrerverhaltens
EP4035008A4 (en) * 2019-09-27 2023-10-11 Tata Consultancy Services Limited METHOD AND SYSTEM FOR DIAGNOSING ANOMALIES IN A MANUFACTURING PLANT
CN110866908B (zh) * 2019-11-12 2021-03-26 腾讯科技（深圳）有限公司 图像处理方法、装置、服务器及存储介质
CN110831050B (zh) * 2019-11-21 2022-09-30 武汉宝久智控科技有限公司 一种传感器节点控制方法及***
DE102019219858A1 (de) * 2019-12-17 2021-06-17 Zf Friedrichshafen Ag Verfahren zur Bestimmung eines Betriebszustands von Systemkomponenten eines Fahrmischers
CN111061235B (zh) * 2019-12-20 2020-12-15 中核控制***工程有限公司 一种具有故障预警功能的dcs设备诊断方法
US20230091197A1 (en) * 2020-02-27 2023-03-23 Siemens Healthcare Diagnostics Inc. Automatic sensor trace validation using machine learning
CN111539482B (zh) * 2020-04-28 2021-03-09 三峡大学 基于rbf核函数的空间多维风电功率数据降维及重构方法
US11954309B2 (en) * 2020-05-04 2024-04-09 Adobe Inc. Systems for predicting a terminal event
CN111912638B (zh) * 2020-06-13 2021-12-21 宁波大学 一种在线故障根源识别的精馏塔故障诊断方法
DE102020123721A1 (de) 2020-09-11 2022-03-17 Volkswagen Aktiengesellschaft Zustandsbestimmungsverfahren für eine Fahrzeugkomponente, Schaltung, Kraftfahrzeug
CN112257754B (zh) * 2020-09-24 2023-07-28 北京航天测控技术有限公司 航天器运行状态的分析方法和装置
CN113027658B (zh) * 2021-03-24 2022-08-02 华中科技大学 一种水轮机转轮实时状态评估方法及其应用
CN113442691A (zh) * 2021-06-18 2021-09-28 科大讯飞股份有限公司 智能车膜的控制方法、装置和存储介质及电子设备
EP4375826A1 (en) * 2021-07-19 2024-05-29 Kabushiki Kaisha Toshiba Information processing device, information processing method, and program
JP2023042945A (ja) * 2021-09-15 2023-03-28 株式会社東芝 監視装置、方法およびプログラム
CN114228637B (zh) * 2021-12-02 2024-02-20 科大讯飞股份有限公司 一种车辆断电保护方法、装置、存储介质及设备
CN116311374B (zh) * 2023-03-27 2023-10-20 淮阴工学院 一种化工厂工人异常行为识别与预警方法及***
CN116957543B (zh) * 2023-09-19 2023-12-22 成都秦川物联网科技股份有限公司 基于大数据的智慧燃气设备管理方法和物联网***
CN118051796B (zh) * 2024-04-16 2024-06-18 自贡市第一人民医院 一种消毒供应中心监测数据智能分析方法

Family Cites Families (7)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US5680409A (en) * 1995-08-11 1997-10-21 Fisher-Rosemount Systems, Inc. Method and apparatus for detecting and identifying faulty sensors in a process
US6952662B2 (en) 2000-03-30 2005-10-04 Smartsignal Corporation Signal differentiation system using improved non-linear operator
JP4041395B2 (ja) * 2002-12-27 2008-01-30 東京瓦斯株式会社 故障検出装置及び故障検出方法
JP5431235B2 (ja) 2009-08-28 2014-03-05 株式会社日立製作所 設備状態監視方法およびその装置
JP5813317B2 (ja) 2010-12-28 2015-11-17 株式会社東芝 プロセス状態監視装置
JP5259797B2 (ja) * 2011-09-05 2013-08-07 株式会社東芝 学習型プロセス異常診断装置、およびオペレータ判断推測結果収集装置
JP2013143009A (ja) * 2012-01-11 2013-07-22 Hitachi Ltd 設備状態監視方法およびその装置

Cited By (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
EP3462267A1 (en) 2017-09-27 2019-04-03 Panasonic Intellectual Property Corporation of America Anomaly diagnosis method and anomaly diagnosis apparatus

Also Published As

Publication number Publication date
US9779495B2 (en) 2017-10-03
US20150363925A1 (en) 2015-12-17
WO2014115615A1 (ja) 2014-07-31
JP2014142697A (ja) 2014-08-07

Similar Documents

Publication Publication Date Title
JP6076751B2 (ja) 異常診断方法およびその装置
JP6216242B2 (ja) 異常検知方法およびその装置
JP5301717B1 (ja) 設備状態監視方法およびその装置
EP2905665B1 (en) Information processing apparatus, diagnosis method, and program
JP5342708B1 (ja) 異常検知方法及びその装置
JP7017363B2 (ja) 異常検知装置および異常検知方法
JP5945350B2 (ja) 設備状態監視方法およびその装置
JP5431235B2 (ja) 設備状態監視方法およびその装置
JP5331774B2 (ja) 設備状態監視方法およびその装置並びに設備状態監視用プログラム
JP5538597B2 (ja) 異常検知方法及び異常検知システム
JP6223936B2 (ja) 異常傾向検出方法およびシステム
JP2013143009A (ja) 設備状態監視方法およびその装置
JP6076421B2 (ja) 設備状態監視方法およびその装置
WO2013011745A1 (ja) 設備状態監視方法およびその装置
JP5849167B1 (ja) 異常検知方法およびその装置
JP2010191556A (ja) 異常検知方法及び異常検知システム
JP2023106472A (ja) 異常検知装置及び異常検知方法
JP6213655B2 (ja) 時系列データの異常監視装置、異常監視方法及びプログラム

Legal Events

Date Code Title Description
RD02 Notification of acceptance of power of attorney

Free format text: JAPANESE INTERMEDIATE CODE: A7422

Effective date: 20140908

A621 Written request for application examination

Free format text: JAPANESE INTERMEDIATE CODE: A621

Effective date: 20150430

A131 Notification of reasons for refusal

Free format text: JAPANESE INTERMEDIATE CODE: A131

Effective date: 20160223

A521 Request for written amendment filed

Free format text: JAPANESE INTERMEDIATE CODE: A523

Effective date: 20160418

A131 Notification of reasons for refusal

Free format text: JAPANESE INTERMEDIATE CODE: A131

Effective date: 20160712

A521 Request for written amendment filed

Free format text: JAPANESE INTERMEDIATE CODE: A523

Effective date: 20160831

TRDD Decision of grant or rejection written
A01 Written decision to grant a patent or to grant a registration (utility model)

Free format text: JAPANESE INTERMEDIATE CODE: A01

Effective date: 20161227

A61 First payment of annual fees (during grant procedure)

Free format text: JAPANESE INTERMEDIATE CODE: A61

Effective date: 20170111

R150 Certificate of patent or registration of utility model

Ref document number: 6076751

Country of ref document: JP

Free format text: JAPANESE INTERMEDIATE CODE: R150