JP6690670B2 - 判定装置、判定方法および判定プログラム - Google Patents

Description

本件は、判定装置、判定方法および判定プログラムに関する。

異常を判定する技術が開発されている。例えば、モデルとの乖離度を判断基準に用いて異常を判定する技術（例えば、特許文献１参照）や、モデルとの差分の積分を判断基準に用いて異常を判定する技術（例えば、特許文献２参照）が開示されている。

国際公開第２０１０／０８２３２２号 特公平４−２５５６５号公報

しかしながら、上記技術では、異常の予兆を判定することは困難である。

本件は上記課題に鑑みなされたものであり、異常の予兆を判定することができる判定装置、判定方法および判定プログラムを提供することを目的とする。

一つの態様では、判定装置は、複数の第１センサ検出値の基準モデルを作成するモデル作成部と、所定の時点から、第２センサ検出値と作成された前記基準モデルとから求められた推定誤差が第１閾値を超過するか否かを判定し、前記第２センサ検出値と作成された前記基準モデルとから求められた前記推定誤差の積算値が第２閾値を超過するか否かを判定し、前記推定誤差が前記第１閾値を超えるまでの第１時間が第１の所定時間よりも短いか否かを判定し、前記推定誤差の積算値が前記第２閾値を超えるまでの第２時間が第２の所定時間より短いか否かを判定する判定部と、前記判定部が、前記推定誤差が前記第１閾値を超過したと判定し、且つ、前記第１時間が前記第１の所定時間より短いと判定した場合、及び前記積算値が前記第２閾値を超過したと判定し、且つ、前記第２時間が前記第２の所定時間より短いと判定した場合、異常に係る信号を出力する出力部と、を備える。

異常の予兆を判定することができる。

（ａ）および（ｂ）は光ファイバによる温度測定方式の例である。 目的変数および説明変数群の一例を例示する図である。 （ａ）は燃料種Ａを用いた場合の推定値と実測値との差を温度差として例示したイメージ図であり、（ｂ）は燃料種Ｂを用いた場合の推定値と実測値との差を温度差として例示したイメージ図である。 閾値と異常判定との関係を例示する図である。 （ａ）は第１実施形態に係る判定装置の概略図であり、（ｂ）は判定部のハードウェア構成を説明するためのブロック図である。 説明変数群を例示する。 目的変数の推定式の設定、ならびに乖離度および推定有効時間の閾値の設定の処理を例示するフローチャートである。 異常判定を行う際に実行されるフローチャートの一例である。 （ａ）は第２実施形態に係る判定装置の概略図であり、（ｂ）は温度センサの例である。 異常判定を行う際に実行されるフローチャートの一例である。 （ａ）は目的変数の実測値であり、（ｂ）は説明変数の実測値である。 （ａ）は推定誤差の瞬時値であり、（ｂ）は推定誤差の積算値であり、（ｃ）は推定有効時間である。 温度センサの例である。 （ａ）はマハラノビス距離であり、（ｂ）は推定有効時間である。 変形例２の判定システムを例示する図である。 （ａ）および（ｂ）はセンサ部および測定器を例示する図である。 （ａ）〜（ｃ）はセンサ部を例示する図である。 無次元化手続きを表すフローチャートを例示する図である。 （ａ）および（ｂ）は異常予兆を例示する図である。 異常判定を行う際に実行されるフローチャートの一例である。 比較例を例示するフローチャートである。 比較例の結果を例示する図である。 図２０に記載の方法の結果を例示する図である。 比較例の結果を例示する図である。 図２０に記載の方法の結果を例示する図である。 異常判定を行う際に実行されるフローチャートの一例である。 異常判定の結果を例示する図である。 センシングデータの規格化値を例示する図である。

まず、異常判定の概略について説明する。

化学プラント、製油工場、火力発電所等において、腐食等に起因してガス漏れ等が生じることがある。そこで、早期に異常を判定できることが望まれている。例えば、主配管に枝配管が溶接された配管系に対して、温度センサを配置し、ガスや液体の漏洩を温度変化として早期に検出することが考えられる。または、冷却水の配管温度を監視することで、万が一の冷却不具合等が発生しても、火災になる前に早期の温度異常を検知することができる。なお、以下における異常の「予兆の検知」とは、「ある対象を監視する場合に、何らかの手段を用いてその動作の一部を可視化し、その可視化された状態が正常な状態と異なるのではないか、と判断される場合に異常の予兆を検知したとする」というものである。

温度異常の検知に、例えば、ラマン散乱光を測定して温度情報とする光ファイバによる温度測定方式が挙げられる。例えば、図１（ａ）で例示するように、光ファイバを枝配管に敷設することで、漏えいを温度変化として早期に検出することができる。また、図１（ｂ）で例示するように、ボイラの冷却水配管に光ファイバを敷設することで、冷却水温度を監視することができる。それにより、冷却不具合等が発生しても、火災になるまえに早期の温度異常を検知することができる。

しかしながら、設備管理としてより好適なのは「事故を最小限に食い止める異常の検知」ではなく、「異常となる兆候を予測することで事故を未然に防止すること」にある。そこで、種々の運転制御用のセンサデータ等のパラメータを説明変数とし、後付のセンサ情報を目的変数とした回帰分析を行ったり、後付のセンサ情報同士の相関分析を行い、問題のない状態か否かを統計的に判断する手法が採られるようになってきた。これは、運転制御用のセンサデータ等のパラメータを一括管理できる程度にネットワーク技術が向上してきたこと、コンピュータパワーの増加により、種々の回帰分析手法や相関分析手法をリアルタイムに適用できるようになってきたことに起因する。

しかしながら、上記手法は、いかに「推定誤差」（推定値と実測値との差）を最小化するかにフォーカスされている。したがって、「推定誤差」がどの程度になると異常と判定するのかについては、実用的な議論がなされていない。例えば、火力発電所においては、油種、石炭種等はパターン分類されるものではなく、産出国やその混合量等により、毎回特性が変わる。したがって、特性が変わるたびに推定のための初期化作業を行う必要が生じる。しかしながら、その精度により「推定誤差」が影響を受けるため、たとえば、一定時間に「推定誤差」がどの程度変化したか等をもって異常と判断する、といった設定の簡単な方法のみが採用されやすくなる。その結果、あいまいな閾値が設定されて、非実用的なシステムになりやすい。

図２は、目的変数および説明変数群の一例を例示する図である。図２の例では、目的変数１〜３は、ボイラの外壁メタルの各箇所の温度である。説明変数群は、目的変数と相関を有するセンサの出力値である。目的変数１〜３の実測値は、外壁メタルに設置された温度センサなどを用いて取得することができる。図２で例示されるように、目的変数１〜３の推定式は、各説明変数の係数および定数を設定することによって得られる。これらの係数および定数については、過去の各説明変数および目的変数の実測値に基づいて、最小２乗回帰（Ordinary Least Mean Square)や主成分回帰（Principal Components Regression）、偏最小２乗回帰（Partial Least Squares)等により設定することができる。説明変数情報が収集されるたびに、推定式に対して説明変数が入力され、推定値が導出される。

同時刻に取得された目的変数の実測値と推定式の推定値とを比較することによって、システムが正常なのか異常なのか判定することができる。推定式の係数および定数を設定するには、一定期間の直近の過去データが必要である。この直近の期間を「モデリング期間」と称する。これに対して、実際に推定値と実測値との比較を行う期間を「スコアリング期間」と称する。スコアリング期間において、「推定誤差」が一定の値を超えれば推定から外れている状況が発生していると考えられる。

図３（ａ）は、燃料種Ａを用いた場合のスコアリング期間における目的変数１および目的変数２の推定値と実測値との差を温度差として例示したイメージ図である。図３（ｂ）は、燃料種Ｂを用いた場合のスコアリング期間における目的変数１および目的変数２の推定値と実測値との差を温度差として例示したイメージ図である。

燃料種Ａの場合は、しきい値が３σで±０．３℃であるとする。燃料種Ｂの場合は、しきい値が３σで±０．９℃であるとする。しかしながら、この３σを超えた状態が直ちに異常といえるか否かの客観的な判断は難しい。なぜならば０．３％の確率でこの値を超える可能性があるからである。逆に４σを設定した場合、異常を見逃す可能性も考えられる。

図４は、閾値と異常判定との関係を例示する図である。図４で例示するように、閾値を比較的大きい閾値１に設定した場合は、真に問題が発生した場合に対して遅れが生じてしまい、対策が遅れてしまう。閾値を閾値１よりも小さい閾値２に設定した場合には、真に問題が発生する状況とは異なる状況で異常と判定されてしまうため、こちらも予兆検知の機能を果たせない。つまり、「十分に高精度な推定を行い、かつ、目的変数ごとに有用な閾値を設定し、その閾値で異常判定を行う」ことが出来ないと、有用な予兆検知は困難である。定期点検や油種の割合変更などの閾値に影響を与える内容の更新ごとに適正な閾値を検討することは事業継続において現実的ではない。

以下の実施形態では、異常の予兆を判定することができる判定装置、判定方法および判定プログラムについて説明する。

（第１実施形態）
図５（ａ）は、第１実施形態に係る判定装置１００のブロック図である。本実施形態においては、判定装置１００は、一例として石炭の燃焼サイクルを利用する火力発電設備に設置される。判定装置１００は、説明変数取得部１０、複数の温度センサ２０ａ〜２０ｃ、判定部３０などを備える。判定部３０は、モデル作成部３１、閾値設定部３２、異常判定部３３および出力部３４を備える。

図５（ｂ）は、判定部３０のハードウェア構成を説明するためのブロック図である。図５（ｂ）で例示するように、判定部３０は、ＣＰＵ１０１、ＲＡＭ１０２、記憶装置１０３、インタフェース１０４などを備える。これらの各機器は、バスなどによって接続されている。ＣＰＵ（Ｃｅｎｔｒａｌ Ｐｒｏｃｅｓｓｉｎｇ Ｕｎｉｔ）１０１は、中央演算処理装置である。ＣＰＵ１０１は、１以上のコアを含む。ＲＡＭ（Ｒａｎｄｏｍ Ａｃｃｅｓｓ Ｍｅｍｏｒｙ）１０２は、ＣＰＵ１０１が実行するプログラム、ＣＰＵ１０１が処理するデータなどを一時的に記憶する揮発性メモリである。記憶装置１０３は、不揮発性記憶装置である。記憶装置１０３として、例えば、ＲＯＭ（Ｒｅａｄ Ｏｎｌｙ Ｍｅｍｏｒｙ）、フラッシュメモリなどのソリッド・ステート・ドライブ（ＳＳＤ）、ハードディスクドライブに駆動されるハードディスクなどを用いることができる。ＣＰＵ１０１が記憶装置１０３に記憶されている判定プログラムを実行することによって、判定部３０にモデル作成部３１、閾値設定部３２、異常判定部３３および出力部３４が実現される。なお、モデル作成部３１、閾値設定部３２、異常判定部３３および出力部３４は、専用の回路などのハードウェアであってもよい。

説明変数取得部１０は、各説明変数を取得する。図６は、説明変数群を例示する。図６で例示するように、説明変数は、電力量、給炭量、内部温度１、内部温度２、空気流量、圧力１、圧力２、圧力３、通風口１、通風口２、運転回転数、運転率、運転振動数などを含む。各説明変数は、目的変数１〜３（温度センサ２０ａ〜２０ｃの検出温度）と相関を有している。また、各説明変数は、互いに独立とみなせる（多重共線性が低い）ことが好ましい。電力量は、火力発電によって得られる発電電力である。給炭量は、炉に供給される石炭量である。内部温度１，２は、例えば炉の内部のいずれかの箇所の温度である。空気流量は、炉に供給される空気の流量である。圧力１〜３は、例えば炉に接続される配管内の圧力である。通風口１，２は、通風口の温度などである。運転回転数、運転率および運転振動数は、炉の運転回転数、運転率および運転振動数などである。これらの説明変数は、各センサの出力値である。

複数の温度センサ２０ａ〜２０ｃは、例えば炉の外部の壁面の互いに異なる箇所に設置される。本実施形態においては、３箇所に温度センサが設置される。温度センサ２０ａ〜２０ｃの温度測定方式は、例えば、光ファイバ内のラマン散乱光を用いた方式を用いることができる。例えば、正確な温度を測定するために、ほぼ同一の温度と見做すことができる小領域に２ｍ程度の長さの光ファイバを捲回する。このようにすることで、各捲回部がそれぞれ１つの温度センサとして機能する。本実施形態においては、温度センサ２０ａ〜２０ｃが検出する壁面温度１〜３を目的変数１〜３の実測値として用いる。

モデル作成部３１は、説明変数取得部１０が取得した説明変数および温度センサ２０ａ〜２０ｃの検出値を用いて、目的変数１〜３の推定式を作成する。この推定式は、目的変数１〜３の基準モデルである。推定式は、図２で例示した推定式と同様に、各説明変数の係数および定数を設定することによって得られる。これらの係数および定数については、過去の各説明変数および温度センサ２０ａ〜２０ｃの実測値に基づいて、最小２乗回帰、主成分回帰、偏最小２乗回帰等の回帰分析により設定することができる。

閾値設定部３２は、基準モデルに対する目的変数の実測値の乖離度と推定有効時間とに閾値を設定する。基準モデルに対する目的変数の実測値の乖離度は、例えば、推定誤差＝（目的変数の実測値）−（目的変数の推定値）、推定誤差のデータ更新ごとの積算値、目的変数の実測値と目的変数の推定値の比などである。本実施形態においては、上記乖離度として、推定誤差、および推定誤差の積算値を用いる。目的変数の推定値は、推定式に説明変数を入力することによって得られる数値である。目的変数１に係る推定誤差は、（温度センサ２０ａの実測値）−（目的変数１の推定値）である。目的変数２に係る推定誤差は、（温度センサ２０ｂの実測値）−（目的変数２の推定値）である。目的変数３に係る推定誤差は、（温度センサ２０ｃの実測値）−（目的変数３の推定値）である。推定有効時間とは、基準モデルを用いたスコアリング期間において推定誤差の測定を開始してから上記いずれかの乖離度が閾値を超えるまでの時間である。異常判定部３３は、この推定有効時間が閾値を下回るか否かを判定することで、異常判定を行う。出力部３４は、異常判定部３３が異常の判定をした場合に、異常に係る信号を出力する。

図７は、目的変数の推定式の設定、ならびに乖離度および推定有効時間の閾値の設定の処理を例示するフローチャートである。図７で例示するように、閾値設定部３２は、まず、初期条件更新フラグを検知する（ステップＳ１）。初期条件更新フラグは、目的変数の推定式、ならびに乖離度および推定有効時間の閾値の更新のトリガーとなるフラグである。次に、閾値設定部３２は、適当な許容値１，２を設定する（ステップＳ２）。許容値１は、推定誤差に対する閾値である。許容値２は、推定誤差の積算値に対する閾値である。

次に、モデル作成部３１は、モデリング期間のデータセットを収集する（ステップＳ３）。このデータセットには、モデリング期間における所定の時間ごとの説明変数および温度センサ２０ａ〜２０ｃの検出値（実測値）が含まれる。次に、モデル作成部３１は、ステップＳ３で収集したデータセットを用いて、目的変数１〜３の推定式の係数および定数を決定する（ステップＳ４）。ステップＳ４の実行によって、目的変数１〜３の推定式が設定される。

次に、閾値設定部３２は、推定誤差の測定を開始し（スコアリング期間）、例えばスコアリング期間およびモデリング期間の最初の６０回（３０秒周期の測定の場合には３０分）の推定誤差の平均値および標準偏差を求める（ステップＳ５）。次に、閾値設定部３２は、平均値＋１σ値〜３σ値を許容値１として再設定する。また、閾値設定部３２は、おおよそ推定有効時間が測定周期の６０〜２４０倍（３０秒間隔の測定において３０分〜２時間）程度になるように許容値２を再設定する（ステップＳ６）。ステップＳ６は、３０秒間隔の測定においてスコアリング期間が３０分になる前に再推定になるようならば、許容値１，２が小さいということなので、緩和することに相当する。

次に、閾値設定部３２は、ステップＳ６の再設定後、推定有効時間の仮の閾値を設定し、推定誤差の仮測定を開始する（ステップＳ７）。次に、モデル作成部３１は、推定誤差が許容値１を超えるか推定誤差の積算値が許容値２を超えた場合に基準モデルの再作成を行うことを繰り返す。閾値設定部３２は、当該再作成が３０回行われる程度のデータが蓄積されたか否かを判定する（ステップＳ８）。ステップＳ８で「Ｎｏ」と判定された場合、ステップＳ８が再度実行される。ステップＳ８で「Ｙｅｓ」と判定された場合、閾値設定部３２は、ステップＳ８で得られた推定有効時間の平均値と標準偏差を求め、例えば３σ値を用いて推定有効時間の閾値を再設定する（ステップＳ９）。なお、ステップＳ６〜ステップＳ９において、仮の有効時間でも、閾値を下回れば、出力部３４は異常に係る信号を出力する。

許容値１，２および推定有効時間の閾値については、いったんデータを蓄積してしまえばあらかじめ組み込んだプログラムで過去にさかのぼって何度でも設定できるものである。したがって、許容値１，２を再設定したらそこから再度データ蓄積を行って検算する、同じく推定有効時間の閾値を再設定したらそこから再度データ蓄積を行って検算する、ということは不要である。したがって、定期点検や油種の割合変更などの後にシステムに対して、「変更がなされた」という情報のみを入力を行えばよいシステムが構築可能である。

図８は、閾値設定部３２によって許容値１，２および推定有効時間の閾値が設定された後に異常判定部３３が異常判定を行う際に実行されるフローチャートの一例である。異常判定部３３は、閾値設定部３２によって許容値１，２および推定有効時間の閾値が設定された後のデータセットを収集する（ステップＳ１１）。次に、異常判定部３３は、推定誤差が許容値１を超えたか、または推定誤差の積算値が許容値２を超えたか否かを判定する（ステップＳ１２）。ステップＳ１２でいずれも超えていないと判定された場合、ステップＳ１２が再度実行される。ステップＳ１２でいずれかが超えたと判定された場合、異常判定部３３は、推定有効時間が所定の時間（例えば３０秒間隔での測定で１０分）より短いか否かを判定する（ステップＳ１３）。ステップＳ１３で「Ｙｅｓ」と判定された場合、出力部３４は、異常に係る信号を出力する（ステップＳ１４）。ステップＳ１３で「Ｎｏ」と判定された場合、異常判定部３３は、その時点からの過去データ（例えば３０秒間隔の測定で１時間）を用いて再度推定を行い、推定式の係数および定数を更新する（ステップＳ１５）。その後、ステップＳ１１から再度実行される。

図８の例では、許容値１は推定誤差の閾値であり、許容値２は推定誤差のデータ更新ごとの積算値である。一定の推定が出来ていれば長期間で平均すれば推定誤差の平均値はほぼ０になるが、推定と異なる状況が生じてくると正負どちらかの値が大きくなり始める。この変化に対して設定された閾値が許容値２である。たとえ推定誤差の平均値がゼロ近傍になっていたとしても、突発的に生じ始めた事象の場合、推定誤差が大きくなるため、異常とみなすことができる。この場合の値が許容値１である。

許容値１，２は、あいまいさを含んでいるため、測定される推定誤差が許容値1を超過し、または推定誤差の積算値が許容値２を超過することがある。しかしながら、モデル作成部３１は、その時点で再度推定を行い、推定式の係数および定数を更新し、その時点を推定開始時間としてカウントを再開する。次回に許容値１または許容値２を超過するまでの期間を「推定有効時間」とし、その「推定有効時間」に対して閾値を設定するのが図８の方式である。図８の例では許容値を２つ設定したが、１つでも構わないし、もっと多数の設定を行っても構わない。たとえば、１つないし複数個の許容値を超過したうえで、それとは別の１つないし複数個の許容値を超過した段階で「推定有効時間」を確定し、再度推定を行う、という方法でも構わない。

本実施形態においては、閾値を超過する現象が頻繁に発生することが許容される。その代り、その頻度がどの程度の間隔で発生するのかが重要視される。発電所やプラント、その他の想定適用先においても、現象は徐々に変化する。したがって、基準モデルの作成後において推定誤差が小さくなる傾向にある。その状態で推定誤差が大きければそれは定式で考慮されない事象が発生しているということである。したがって、「推定有効時間」が短くなっているということは、閾値の設定や推定のあいまいさを考慮しても異常な状況が発生していると判断できるものである。すなわち、本実施形態によれば、異常の予兆を判定することができる。

また、本実施形態によれば、基準モデルに対する温度センサ２０ａ〜２０ｃの乖離度の閾値は、基準モデルの作成後の一定期間の乖離度を基に決定されることになる。この場合、乖離度の閾値の設定精度が向上する。また、本実施形態によれば、推定有効時間の閾値は、乖離度が閾値を超えるまでの時間のばらつきを基に決定される。この場合、推定有効時間の閾値の設定精度が向上する。

（第２実施形態）
第１実施形態においては、センサ検出値と、センサ検出値と相関を有する複数の他のセンサの検出値とを用いて基準モデルを作成したが、これに限られない。第２実施形態においては、複数のセンサの検出値同士の相関関係を用いて、センサ検出値の基準モデルを作成する。

まず、複数の温度センサが検出する温度推移データの相関関係から、早期に異常の兆候を見出す予兆検知方法について検討してみる。この場合、対象とするデータ群の平均および分散共分散行列から算出するマハラノビス平方距離を求める方法、データの中心および拡がりをロバストに推定し、マハラノビス平方距離に相当する統計量を求めるＭＳＤ法などを用いることができる。

これらは、「外れ値検出」と称されるものである。このマハラノビス平方距離を求める手続きを用いて具体的に早期に異常を見出す予兆検知を行う手順を、時刻ｔにおけるｎ個の温度データＴ１（ｔ），Ｔ２（t），Ｔ３（ｔ），…，Ｔｎ（t）として以下の（１）〜（３）に示す。
（１）モデリング期間（現在時刻よりも古い、データが蓄積された一定期間）を設定し、その期間の各ｎ個の温度センサの温度データの平均および温度データ群の不偏分散共分散行列と逆行列を求める。
（２）モデリング期間内の各時刻Ｔｍ（ｍ＝０，１，２，…）の温度データＴ１（ｍ）〜Ｔｎ（Ｔｍ）に対して、（１）を用いて「モデリング期間のｎ個の温度センサの温度データ群」に対するマハラノビス平方距離を求める。それらの値の標準偏差（３σ等）から異常となる閾値を設定する。
（３）新たな時刻ＴのデータセットＴ１（Ｔ）〜Ｔｎ（Ｔ）が得られるたびに、この「モデリング期間のｎ個の温度センサの温度データ群」との間のマハラノビス平方距離を求め、それが閾値以下か否かを判定する。

ＭＳＤ法を用いる場合も、モデリング期間を設定し、その期間から閾値を決定して、新たなデータセットと閾値とを逐次的に比較するという手順は同じである。つまり、各測定データ同士の相関を用いた「外れ値検出」による予兆検知においても、閾値の設定により、予兆検知の精度が決まってしまう。すなわち、「十分に高精度な推定を行い、かつ、有用な閾値を設定し、その閾値で異常判定を行う」ことが出来ないと、有用な予兆検知は困難である。そこで、第２実施形態においても、異常の予兆を判定することができる判定装置、判定方法および判定プログラムについて説明する。

図９（ａ）は、第２実施形態に係る判定装置１００ａの概略図である。判定装置１００ａが第１実施形態の判定装置１００と異なる点は、説明変数取得部１０が設けられておらず、温度センサ２０ａ〜２０ｃの代わりに温度センサ２０が設けられている点である。温度センサ２０は、温度値が互いに相関を有する複数個所の温度を検出する。例えば、図９（ｂ）で例示するように、温度センサ２０は、同一光ファイバの異なる長さ位置の後方散乱光によって得られる結果を基に各箇所の温度を検出する。図９（ｂ）の例では、各捲回部がそれぞれ個別の温度センサとして機能する。判定部３０の構成は第１実施形態と同様である。

本実施形態においては、モデル作成部３１は、温度センサ２０が検出する各箇所の検出値の平均および分散共分散行列から算出するマハラノビス平方距離を求めることで、基準モデルを作成する。また、モデル作成部３１は、上記検出値の中心および拡がりをロバストに推定し、マハラノビス平方距離に相当する統計量を求めるＭＳＤ法などを用いることで、基準モデルを作成する。ここでの基準モデルは、複数の検出値同士の相関関係の大きさ（各２つのセンサ同士に注目すれば、方位）を反映したばらつきの中心である。

閾値設定部３２は、基準モデルに対するセンサ検出値の実測値の乖離度と推定有効時間とに閾値を設定する。基準モデルに対する目的変数の実測値の乖離度は、マハラノビス平方距離、マハラノビス平方距離のデータ更新ごとの積算値などである。本実施形態においては、上記乖離度として、マハラノビス平方距離およびマハラノビス平方距離の積算値を用いる。推定有効時間とは、基準モデルの作成後のスコアリング期間においてマハラノビス平方距離の測定を開始してから上記いずれかの乖離度が閾値を超えるまでの時間である。異常判定部３３は、この推定有効時間が閾値を下回るか否かを判定することで、異常判定を行う。出力部３４は、異常判定部３３が異常の判定をした場合に、異常に係る信号を出力する。

基準モデルの設定、ならびに乖離度および推定有効時間の閾値の設定については、図７と同様の処理によって行うことができる。なお、閾値設定部３２は、許容値１，２の他に、許容値３を予め設定しておく。図１０は、閾値設定部３２によって許容値１，２および推定有効時間の閾値が設定された後に異常判定部３３が異常判定を行う際に実行されるフローチャートの一例である。異常判定部３３は、閾値設定部３２によって許容値１，２および推定有効時間の閾値が設定された後のデータセットを収集する（ステップＳ２１）。このデータセットは、温度センサ２０の各箇所の検出値ごとに収集される。

次に、異常判定部３３は、いずれかの箇所のマハラノビス平方距離が許容値１を超えたか、または当該箇所のマハラノビス平方距離の積算値が許容値２を超えたか否かを判定する（ステップＳ２２）。ステップＳ２２でいずれも超えていないと判定された場合、ステップＳ２２が再度実行される。ステップＳ２２でいずれかが超えたと判定された場合、異常判定部３３は、マハラノビス平方距離が許容値３を超えたか否かを判定する（ステップＳ２３）。

ステップＳ２３で「Ｙｅｓ」と判定された場合、出力部３４は、異常に係る信号を出力する（ステップＳ２４）。ステップＳ２３で「Ｎｏ」と判定された場合、異常判定部３３は、推定有効時間が所定の時間（例えば３０秒間隔での測定で１０分）未満であるか否かを判定する（ステップＳ２５）。ステップＳ２５で「Ｙｅｓ」と判定された場合、出力部３４は、異常に係る信号を出力する（ステップＳ２６）。

ステップＳ２５で「Ｎｏ」と判定された場合、異常判定部３３は、その時点から過去データ（例えば３０秒間隔の測定で１時間）を用いて再度マハラノビス平方距離やＭＳＤ法におけるマハラノビス平方距離に相当する統計量を導出するためのパラメータを求める（ステップＳ２７）。この場合のパラメータとして、温度センサ２０の各箇所における新規モデリング期間における平均値や不偏分散共分散行列および逆行列などが含まれる。その後、ステップＳ２１から再度実行される。

図１０の例は、３つの許容値を設定し、マハラノビス平方距離が許容値１を超えるか積算値が許容値２を超えた場合は新しいデータ群を用いて再度推定を行うが、その前に、許容値３を超えていたら異常と判定する、という例である。たとえば、マハラノビス平方距離を導出する方法の場合、逐次求められるマハラノビス平方距離と許容値１とを比較し、そのマハラノビス平方距離の積算値と許容値２を比較する、といった方法が挙げられる。また、たとえば温度を複数点測定していて、その複数点の温度が揃って危険な温度に推移した場合などではマハラノビス平方距離は小さいが、異常状態と言える。このような状況を回避するために、許容値３を設定し、測定データそのものと比較する。許容値３は固定数値とせず、測定点数それぞれに対して固有の値を決めてもよい。許容値の設定数や設定方法は図１０の例に限られない。

本実施形態においても、推定有効時間が短くなった場合に異常と判定されるため、異常の予兆を判定することができる。また、基準モデルに対する温度センサ２０の各検出値の乖離度の閾値は、基準モデルの作成後の一定期間の乖離度を基に決定されることになる。この場合、乖離度の閾値の設定精度が向上する。また、推定有効時間の閾値は、乖離度が閾値を超えるまでの時間のばらつきを基に決定される。この場合、推定有効時間の閾値の設定精度が向上する。

上記実施形態に従い、具体的な実施例について説明する。実施例１においては、第１実施形態に従い、図６と同様の目的変数および説明変数を用いた。実施例１に係るシステムは、壁面温度１〜３が適正な範囲に保たれているのかを説明変数群から予測し、最も効率的な状態で運転サイクルを回すことを目的としている。運転サイクルが高温・高圧になりすぎればボイラ以前の燃焼が発生し、低温になりすぎればボイラでの燃焼効率が低下する。最適な温度・圧力を、燃焼を回避してコントロールすることが求められる。

これら説明変数に対して、図７の処理に従って、閾値を以下の通り設定した。本実施例においては、目的変数１〜３に対して、同じ数値を用いた。
許容値１ ±２℃
許容値２ ±１０℃
有効時間の閾値 ２０分未満

毎回の推定毎のモデリング時間は１時間とし、推定には主成分回帰を用いた。主成分回帰の具体的な方法を以下に示す。
（１）モデリング期間の説明変数群と目的変数１の値を用いて１４行１４列の分散共分散行列を生成する。具体的には、目的変数１の分散が１４行１４列目になるように最後段に配置する。
（２）目的変数の行／列を除く、１３行１３列に対して、逆行列を生成する。
（３）（２）で求めた逆行列と（１）で求めた１４列目の１４行目の要素を除く積を求めて、１３個の数値を得て、これらを説明変数の各係数とする。
（４）モデリング期間の各説明変数の平均値を算出し、この各平均値に（３）で求めた係数を掛けてそれらの和を取る。
（５）モデリング期間の目的変数の平均値を算出してから（４）で求めた値を差し引き、これを推定式の定数とする。

本実施例においては、システム異常をシステム運用停止で代用する。図１１（ａ）および図１１（ｂ）で例示するように、１７：２０前後からシステムの運転状態が変更されるが、特にそれに伴う異常は発生していないものの、１８：１８に各説明変数が急激に変動する。これはシステムの運用を停止したからである。この１８：１８のシステム運用停止を数分間でも早く推定できればシステムとして有効である。

図１２（ａ）は、推定誤差の瞬時値である。図１２（ｂ）は、推定誤差の積算値である。図１２（ｃ）は、推定有効時間である。図１２（ｂ）で積算値がゼロにリセットされている位置と図１２（ｃ）のデータ位置とが同じである。これは、その時点で図１２（ａ）の瞬時値が許容値１を超えるか図１２（ｂ）の積算値が許容値２を超過し、再推定計算を行ったからである。図１２（ａ）を見ると、前述したとおり、再推定開始直後の推定誤差が小さくなっていることがわかる。また、図１１（ａ）、図１１（ｂ）および図１２（ｃ）を比較すると、推定有効時間が１８：１８に向かって徐々に低下し、１７：５５から１６分後の１８：１２時点で異常と判断できる条件を満たしているが、この異常の自動判定が合理的であるのがわかる。１８：１８よりも６分前に異常確認が出来たことになり、システムの速やかな停止指令が可能になる。

実施例２は、第２実施形態に従った例である。図１３で例示するように、捲回部の集合を４つ作成し、炉の壁面にそれぞれ密着させた。これら４つの集合のそれぞれを温度センサとして用いた。図１３において、○で描かれた部分が捲回部であり、各捲回部は同一の光ファイバで結ばれている。また、網目模様を利用して温度分布が描かれている。粗い網目部分では温度が低く、細かい網目部分では温度が高くなっている。各機器の内部で部分的に過度な熱の蓄積がはじまれば、部分的に温度が上昇するので、それにより異常を検知することができる。

本実施例では、各４つの捲回部集合からそれぞれ３か所の捲回部を抽出し、計１２個の捲回部の相関分析を行うことで異常の予兆検知を行う。具体的には、各捲回部集合の左下をローカル座標の基準原点（Ｘ，Ｚ）＝（０，０）として、各捲回部が含まれる座標を左下と右上の２点で付与することで領域を設定する。
つまり、
捲回部集合１ 領域１ （Ｘ１ａ１，Ｚ１ａ１），（Ｘ１ａ２，Ｚ１ａ２）
領域２ （Ｘ１ｂ１，Ｚ１ｂ１），（Ｘ１ｂ２，Ｚ１ｂ２）
領域３ （Ｘ１ｃ１，Ｚ１ｃ１），（Ｘ１ｃ２，Ｚ１ｃ２）
捲回部集合２ 領域１ （Ｘ２ａ１，Ｚ２ａ１），（Ｘ２ｂ２，Ｚ２ａ２）
領域２ （Ｘ２ｂ１，Ｚ２ｂ１），（Ｘ２ｂ２，Ｚ２ｂ２）
といったように、領域を指定する。

さらに、各領域に含まれる光ファイバの各位置の温度からそれらの平均値や最高値、最低値などをもとめ、それを各領域の温度データとする。本システムでは、これら１２個の温度に対してそれぞれ閾値を設け、これを許容値３とする。本実施例において、許容値３に相当する温度を超過することはなかったので、許容値１，２と有効時間閾値のみ設定した。

許容値１は、新規データセットの、その時点で用いているモデリング期間のデータに対するマハラノビス平方距離の値に対して設定したものである。許容値２は、１サンプル前のマハラノビス平方距離と新規データセットのマハラノビス平方距離の平均値に対して設定したものである。

積分値を用いないのは、以下の理由による。すなわち、マハラノビス平方距離は、新規データセットがモデリング期間のデータセットの重心に対してどれだけ離れているかを示すものである。マハラノビス平方距離は、どの方向にというベクトル的な成分も考慮して加算すれば重心のまわりで一定の距離をもって回転するようにデータセットが更新される場合と、重心を横切って異なる象限に移行するような場合を異なる数値で示すことができる。しかしながら、スカラー的な加算ではこれらは同じに見えてしまうと考えられるためである。

そこで、少なくとも徐々にデータセットが外れてくれば近傍のデータは同じ方向に行くという仮定の下で、平均値を許容値２として採用することとした。
許容値１ ６０
許容値２ ５０
有効時間の閾値 ２０分未満

マハラノビス平方距離の求め方は前述のとおりである。具体的に、マハラノビス平方距離を求めて予兆検知した結果を図１４（ａ）および図１４（ｂ）で例示する。図１４（ａ）で例示するように、スコアリング期間においてはモデリング期間直後が最もマハラノビス平方距離は小さくなる。当該マハラノビス距離は、時間が経過していくにしたがって増加するが、その増加傾向はそれぞれの時間帯で異なる。１７：０４の再推定では５０分となった推定有効時間は、その後短くなる一方となり、１７：５４の再推定の推定有効時間は１４分、１８：０８の再推定の推定有効時間は８分である。閾値２０分と比較すると１７：５４の１４分後の１８：０８に異常確認が出来たことになる。これは、システムが停止した１８：１８よりも１０分早い段階であり、本例ではシステム停止について記載したが、実際に何らかの事故の兆候が見られる場合に置いても、初動対応が可能な時間が確保できる可能性が見いだせたことになる。

（変形例１）
上記第１実施形態および第２実施形態において、複数の温度センサ２０ａ〜２０ｃまたは温度センサ２０の温度測定方式は、光ファイバ内のラマン散乱光を用いた方式であるが、これに限られるものではない。例えば、温度センサ２０ａ〜２０ｃまたは温度センサ２０として、熱電対、測温抵抗体、カメラ撮像タイプの赤外線サーモグラフィ等が用いられる。

ただし、熱電対や測温抵抗体を用いる場合には、測定ポイント各１点につき、導線同士及び各導線と壁面との間の絶縁が保持された２本の導線が必要である。また、赤外線サーモグラフィを用いる場合には、測定する面の保温等がなされておらず、遮蔽物もない外側から撮像可能な状態となっていること、測定ポイントに該当する位置近傍の複数の温度データから平均化などで測定ポイント１点にまとめること、輝度から温度への正確な変換のために外壁面の放射率をあらかじめ得ておくこと、などが必要である。

また、上記第１実施形態および実施例１では温度を目的変数とし、その他のセンシングデータを説明変数としたが、その他のセンシングデータから目的変数とするものを決めて、測定する温度もその他のセンシングデータとともに説明変数とするなどでも構わない。どのような値に注目するか、ということだけであるため、たとえば、電力量を目的変数として、異常予兆が検知されれば発電システムの効率の劣化が発生している、という使用方法でも構わない。これは、第２実施形態および第２実施例においてもいえることである。例えば、温度センサ２０ではなく、その他のセンシングデータ同士に対して上記実施形態を適用しても構わない。これに関しては、後述の実施例４で述べる。

（変形例２）
図１５は、変形例２の判定システムを例示する図である。実施形態２においては、判定部３０は、温度センサ２０から直接データを取得している。これに対して、変形例２の判定システムは、判定部の機能を有するサーバが、電気通信回線を通じて温度センサからデータを取得するものである。

変形例２の判定システムは、温度センサ２０、サーバ２０２、監視サーバ２０３を含む。温度センサ２０は、測定対象物の温度データを取得するためのセンサ部２１と、センサ部２１から測定データを取得し温度データ生成する測定器２２とを有する。

温度センサ２０は、インターネットなどの電気通信回線２０１を通じてサーバ２０２と接続された構成を有する。また、センサ部２１が設置された測定対象物を監視する監視サーバ２０３は、電気通信回線２０１に接続されている。サーバ２０２は、図５（ｂ）のＣＰＵ１０１、ＲＡＭ１０２、記憶装置１０３、インタフェース１０４などを備え、判定部３０としての機能を実現する。

このような判定ステムにおいて、例えば、日本に設置されているサーバ２０２は、他国の発電所の石炭バンカで測定された測定データを受信し、石炭バンカにおける異常発熱の予兆を検知する。サーバ２０２から出力された結果が、監視サーバ２０３に送信される。

本変形例は、実施形態１に適用することもできる。例えば、温度センサ２０の代わりに複数の温度センサ２０ａ〜２０ｃを用いればよい。

（変形例３）
図１６（ａ）および図１６（ｂ）は、センサ部２１および測定器２２を例示する図である。図１６（ａ）で例示するように、センサ部２１は、例えば、石炭を砕いて微粉を製造するための微粉炭機のミル中間ハウジング部４０の外壁に取り付けられている。図１６（ａ）で例示するように、ミル中間ハウジング部４０は、石炭４１が落下して一時的に貯留される貯留部４２と、貯留部４２に貯留されている石炭４１を粉砕する粉砕リング４３およびローラ４４とを備える。粉砕によって得られた微粉炭４５は、一次空気室４６の空気によって上昇する。

図１６（ｂ）で例示するように、測定器２２は、レーザ１１、ビームスプリッタ１２、光スイッチ１３、フィルタ１４、複数の検出器１５ａ，１５ｂ、演算部１６などを備える。レーザ１１は、半導体レーザなどの光源であり、所定の波長範囲のレーザ光を出射する。例えば、レーザ１１は、所定の時間間隔で光パルス（レーザパルス）を出射する。ビームスプリッタ１２は、レーザ１１が出射した光パルスを光スイッチ１３に入射する。光スイッチ１３は、入射された光パルスの出射先（チャネル）を切り替えるスイッチである。ダブルエンド方式では、光スイッチ１３は、センサ部２１の光ファイバ２３の第１端および第２端に一定周期で交互に光パルスを入射する。シングルエンド方式では、光スイッチ１３は、光ファイバ２３の第１端または第２端のいずれか一方に光パルスを入射する。光ファイバ２３は、温度測定対象の所定の経路に沿って配置されている。

光ファイバ２３に入射した光パルスは、光ファイバ２３を伝搬する。光パルスは、伝搬方向に進行する前方散乱光および帰還方向に進行する後方散乱光（戻り光）を生成しながら徐々に減衰して光ファイバ２３内を伝搬する。後方散乱光は、光スイッチ１３を通過してビームスプリッタ１２に再度入射する。ビームスプリッタ１２に入射した後方散乱光は、フィルタ１４に対して出射される。フィルタ１４は、ＷＤＭカプラなどであり、後方散乱光を長波長成分（ストークス成分）と短波長成分（アンチストークス成分）とを抽出する。検出器１５ａ，１５ｂは、受光素子である。検出器１５ａは、後方散乱光の短波長成分の受光強度を電気信号に変換して演算部１６に送信する。検出器１５ｂは、後方散乱光の長波長成分の受光強度を電気信号に変換して演算部１６に送信する。演算部１６は、ストークス成分およびアンチストークス成分を用いて、光ファイバ２３の延伸方向の温度分布を測定する。

図１７（ａ）は、センサ部２１の透過図であり、図１７（ｂ）のシート２４ｂを透過した図である。図１７（ｂ）は、図１７（ａ）のＡ−Ａ線断面図である。センサ部２１は、光ファイバ２３が所定の位置に配置されたファイバシートである。図１７（ａ）および図１７（ｂ）で例示するように、センサ部２１は、光ファイバ２３を挟んで保持する一対のシート２４ａ，２４ｂと、シート２４ａ，２４ｂ間の隙間を保つガラスクロステープ２５と、光ファイバ２３の捲回部２６の位置を大まかに決めるスリット入り金属管２７とを有する。

光ファイバ２３は、複数回巻きまわされた状態でシート２４ａ，２４ｂに保持された捲回部２６ａ〜２６ｈ（以下、捲回部２６と総称する）を有する。図１７（ｃ）は、光ファイバ２３が１重に巻きまわされた例を表す。捲回部２６ａ〜２６ｈは、１本の光ファイバ２３により構成されるか、上段下段で各１本ずつ、計２本の光ファイバ２３で構成される。後者の場合、たとえば、図１７（ａ）における上段下段接続部で融着接続する。シート２４ａは測定対象に接している。シート２４ａには粘着テープ２８が設けられている。それにより、シート２４ａを温度測定対象物に貼り付けることができる。

捲回部２６ａ〜２６ｈは、それぞれ、例えば２〜８回、巻き回されている。例えば４回巻き回された捲回部である場合を考えると、光ファイバ２３の直径は耐熱温度により異なるが、０．１６〜０．４ｍｍであるので、金属管２７の内径は光ファイバ２３の直径の２倍以上の１〜２ｍｍ程度が必要となる。金属管２７の板厚は０．５ｍｍ程度なので、シート２４ａから２４ｂに向かって２〜３ｍｍ程度の厚みを持っている。

第２実施形態において、各領域に含まれる光ファイバの各位置の温度からそれらの平均値や最高値、最低値などを求め、それを各領域の温度データとしているが、変形例３の捲回部２６は図１７（ｂ）のごとく厚みを持つ構成となっており、捲回部２６のうち、測定対象から離れた部分の温度は測定対象の実際の温度から大きく外れる場合がある。

そこで、捲回部２６ａ〜２６ｈに含まれる複数の測定点が示す温度値の中から、あらかじめ定められた数点、例えば、最も温度の高い値から順に５点を選択して平均値を求め、各圏回部に各捲回部に対応する領域の温度とする。これにより、測定される温度の精度が高まり、測定対象の温度の異常状態の予兆を更に精度よく検出することができる。

（第３実施形態）
第２実施形態とは別の目的を保持する実施形態として、第３実施形態について説明する。図１０で例示した第２実施形態では、さまざまなセンシングデータを用いることを前提としていたものの、その具体的な方法は示していなかった。そこで、本実施形態ではまず、その方法の一例を示す。この方法自体は、特許第５３０８５０１号においても、確率変数への変換として開示されているものと同様で、通常の多変量解析では一般的な方法である。

図１８は、時刻ｔにおけるＮ個のセンシングデータＳ１（ｔ）〜ＳＮ（ｔ）を用いてモデリングをする場合の無次元化手続きを表すフローチャートを例示する図である。この無次元化手続きは、モデル作成部３１が図７のステップＳ３およびステップＳ４を実行する場合、異常判定部３３が図８のＳ１５または図１０のステップＳ２７を実行する場合に行われる。なお、Ｎ個のセンシングデータＳ１（ｔ）〜ＳＮ（ｔ）は、各目的変数および各説明変数である。ここでは、モデル作成部３１を主体とする場合について説明する。

モデル作成部３１は、現在の基準時刻ｔから所定の過去時間（モデリング期間）における各センシングデータの平均値と標準偏差を求める（ステップＳ３１）。時刻ｔにおけるＮ個の各センシングデータＳ１（ｔ）〜ＳＮ（ｔ）のその時点ｔ０からの、モデリングに用いる所定の過去時間ΔＴ間の各平均値Ｓ１＿ａｖｅ〜ＳＮ＿ａｖｅは、下記式で表される。
Ｓ１＿ａｖｅ＝Ａｖｅｒａｇｅ（Ｓ１（ｔ０），…，Ｓ１（ｔ０−ΔＴ）），
Ｓ２＿ａｖｅ＝Ａｖｅｒａｇｅ（Ｓ２（ｔ０），…，Ｓ２（ｔ０−ΔＴ）），
…
ＳＮ＿ａｖｅ＝Ａｖｅｒａｇｅ（ＳＮ（ｔ０），…，ＳＮ（ｔ０−ΔＴ））
標準偏差Ｓ１＿ｓｉｇｍａ〜ＳＮ＿ｓｉｇｍａは、下記式で表される。
Ｓ１＿ｓｉｇｍａ＝Ｓｔａｎｄａｒｄ Ｄｅｖｉａｔｉｏｎ（Ｓ１（ｔ０），…，Ｓ１（ｔ０−ΔＴ）），
Ｓ２＿ｓｉｇｍａ＝Ｓｔａｎｄａｒｄ Ｄｅｖｉａｔｉｏｎ（Ｓ２（ｔ０），…，Ｓ２（ｔ０−ΔＴ）），
…
ＳＮ＿ｓｉｇｍａ＝Ｓｔａｎｄａｒｄ Ｄｅｖｉａｔｉｏｎ（ＳＮ（ｔ０），…，ＳＮ（ｔ０−ΔＴ））

次に、モデル作成部３１は、パラメータ導出の際にはまず各センシングデータ値から各平均値を差し引いてゼロ基準の値とし、さらに標準偏差値もしくはその数倍の値で割り算を行う（ステップＳ３２）。それにより、各センシングデータが無次元化される。異常判定部３３が無次元化手続きを行った場合には、図８のステップＳ１１または図１０のステップＳ２１を実行する際に、新しいデータ群に対して上記の無次元化処理を行う。

上記無次元化処理によれば、各センシングデータはいずれもモデリング期間においては、平均値がゼロになり、また、値のばらつきも一様になる。値のばらつきは後のパラメータ導出の方法により是正可能なので、標準偏差値は用いずに、算出した各平均値で各センシングデータ値の割り算を行い、その結果得られた値でパラメータ導出してもよい。その場合は、モデリング期間における平均値は１になる。この無次元化手続きにより、図６で例示された壁面温度１〜３、電力量、給炭量、圧力１〜３等々の次元の異なるセンシングデータをたとえば図８や図１０の方法で扱うことが可能になる。

ところで、プラントにおける異常予兆とは、ほとんどの場合、何らかの事象が１発だけ起きてすぐ収束するものではなく、図１９（ａ）で例示する石炭の自然発火現象や、図１９（ｂ）で例示する軸受劣化故障のように、じわじわと悪化してあるとき加速度的に頻発しはじめて、それが連続現象となって事故に至る、というものであり、実際の発電所・プラントなどでみられる、突発的にある短時間のみ何らかの異常が発生するが、その後すぐに正常に収まる、という場合などとは切り分けて考える必要がある。なお、図１９（ａ）においては、石炭種Ａ〜Ｈのいずれにおいても、加速度的に自然発火現象が頻発することが表されている。

そこで、図２０の手順を本実施形態の１つの方法として開示する。装置構成は、実施形態２と同様である。異常判定部３３は、閾値設定部３２によって許容値１，２および推定有効時間の閾値が設定された後のデータセットを収集する（ステップＳ４１）。このデータセットは、温度センサ２０の各箇所の検出値ごとに収集される。

次に、異常判定部３３は、いずれかの箇所のマハラノビス平方距離が許容値１を超え、かつ当該箇所のマハラノビス平方距離の積算値が許容値２を超えたか否かを判定する（ステップＳ４２）。ステップＳ４２で「Ｎｏ」と判定された場合、ステップＳ４１から再度実行される。ステップＳ４２で「Ｙｅｓ」と判定された場合、異常判定部３３は、マハラノビス平方距離が許容値３を超えたか否かを判定する（ステップＳ４３）。

ステップＳ４３で「Ｙｅｓ」と判定された場合、出力部３４は、異常に係る信号を出力する（ステップＳ４４）。ステップＳ４３で「Ｎｏ」と判定された場合、異常判定部３３は、その時点から過去データ（例えば３０秒間隔の測定で１時間）を用いて再度マハラノビス平方距離やＭＳＤ法におけるマハラノビス平方距離に相当する統計量を導出するためのパラメータを求める（ステップＳ４５）。この場合のパラメータとして、温度センサ２０の各箇所における新規モデリング期間における平均値や不偏分散共分散行列および逆行列などが含まれる。

次に、異常判定部３３は、ステップＳ４１で収集した現在のデータセットを用いて以下のステップＳ４７，Ｓ４８の外れ値検定を開始する（ステップＳ４６）。まず、異常判定部３３は、いずれかの箇所のマハラノビス平方距離が許容値１を超え、かつ当該箇所のマハラノビス平方距離の積算値が許容値２を超えたか否かを判定する（ステップＳ４７）。ステップＳ４７で「Ｎｏ」と判定された場合、異常判定部３３は、推定有効時間が所定の時間（例えば３０秒間隔での測定で１０分）未満であるか否かを判定する（ステップＳ４８）。ステップＳ４８で「Ｙｅｓ」と判定された場合、異常判定部３３は、異常に係る信号を出力する（ステップＳ４９）。

ステップＳ４７で「Ｎｏ」と判定された場合、異常判定部３３は、推定有効時間を、所定の時間よりも大きな固定値にセットする（ステップＳ５０）。その後、ステップＳ４１から再度実行される。ステップＳ４８で「Ｎｏ」と判定された場合、異常判定部３３は、推定有効時間をゼロにリセットする（ステップＳ５１）。その後、ステップＳ４１から再度実行される。

図２０の例では、図１０と異なり、モデリングの後に新しいデータセットを用いて外れ値検定を行うのではなく、現在のデータセットで再度外れ値検定を行い、その結果により推定有効時間をリセットするか否かを決めている。これは、モデリングが行われずに一定程度時間が経過している状態である時点で新しいデータセットで許容値１または許容値２を超えた場合に、それが実際に突発的な事象の発生に起因したものなのか、それとも、単純に炭種変更や出力指令変更に起因したものかを切り分けるためである。この技術思想は、実施形態１にも適用することができる。もしモデリング後に再度許容値１または許容値２を超過したとすればそれは前者に相当すると考えてよく、推定有効時間を評価できる状態にすることで、加速度的な状況に入っているのか否かを判断できるようになる。この方法を用いる際には、データセットの次元に関わらず、図１８に開示した方法、もしくは上述の通り、標準偏差で割り算するのではなく、平均値で割り算する方法を用いて無次元化を行うのがより好適である。なお、図２０のステップＳ４２およびステップＳ４７において、「ａｎｄ」判定を「ｏｒ」判定としてもよい。

図２０で例示する第３実施形態と、特許第５３０８５０１号などでも使われている標準的なムービングウィンドウという、逐次的に外れ値検定の前段階で必ずモデリングを行う図２１の比較例とを比較し、第３実施形態の効果を明確化する。図２１の比較例では、セットされたデータ群を用いて（ステップＳ６１）、いずれかの箇所のマハラノビス距離が許容値１および許容値２を超えた場合（ステップＳ６２）に異常に係る信号が出力される（ステップＳ６３）。ステップＳ６２で「Ｎｏ」の場合にマハラノビス距離が許容値３を超えた場合（ステップＳ６４）にも異常に係る信号が出力される（ステップＳ６５）。ステップＳ６４で「Ｎｏ」の場合に、その時点から過去データ（例えば３０秒間隔の測定で１時間）を用いて再度マハラノビス平方距離やＭＳＤ法におけるマハラノビス平方距離に相当する統計量を導出するためのパラメータが求められる（ステップＳ６６）。

図２０および図２１のいずれの例においても、図１８で例示した無次元化処理を行った。但し、ここでは上述に示したもう一つの方法である、標準偏差で割り算する方法ではなく平均値で割り算する方法を用いた。対象として第１実施形態に示した石炭の燃焼サイクルを利用する火力発電設備を例にする。図６の壁面温度１〜３の代わりに、図９（ｂ）、図１３、図１７で示される、多数の捲回部を持つシート形状に加工したものを図１６のごとく設備の壁面に設置し、６４点の測定データを抽出し、それを測定毎のデータセットとして扱う。

図２２に比較例の結果を例示する。瞬間値とは、新しいデータセットを用いて計算した外れ値距離であり、２点平均とは、前回計算した外れ値距離と今回計算した外れ値距離の２点平均値のことである。許容値１はモデル更新期間の瞬間値の平均値＋瞬間値の標準偏差の３倍（３σ）値とし、許容値２は２点平均値の平均値＋瞬間値の標準偏差の３倍（３σ）値とした。許容値３は瞬間値の平均値＋瞬間値の標準偏差の８倍（８σ）値とした。異常判定は、図２１で異常でない場合には０を出力し、異常確定の場合に１を出力するものである。図２２の例では、システムが停止した１８：１８よりも４分はやい１８：１４に異常判定が頻発している。しかしながら、図２２では午前０時から３回誤検知が発生しており、システムとしての信頼性が損なわれていることがわかる。

図２３に図２０に記載の方法の結果を例示する。瞬間値、２点平均、異常判定、許容値１〜３の定義は図２２と同じである。推定有効時間の閾値は４分とした。システムが停止した１８：１８よりも２分早いが図２１の方法よりも２分遅い１８：１６で異常判定が頻発している。しかしながら、図２２と異なり、誤検知がなく、異常判定の信頼度が保たれていることがわかる。

図２２および図２３に対して、許容値１を瞬間値の平均値＋瞬間値の標準偏差の２倍（２σ）値とし、許容値２を２点平均の平均値＋瞬間値の標準偏差の２倍（２σ）値に変更したものを図２４および図２５に示す。いずれも異常判定が頻発し始める時刻は図２２および２３と変化はない。図２１に記載の方法ではさらに誤検知が頻発しているが、図２０の方法では誤検知が発生しておらず、異常判定の信頼度が保たれていることがわかる。

すなわち、特許第５３０８５０１号にみられるようなムービングウィンドウ型の外れ値検定等とは異なる本実施形態を用いることで、許容値１，２について、平均値＋標準偏差の２倍（２σ）〜標準偏差の３倍（３σ）というあいまいな幅で設定しても予兆検知を行うことが可能となるといえる。

実施例４では、温度だけでなく、第１実施形態で例示した石炭の燃焼サイクルを利用する火力発電設備の様々なセンシングデータに対して図２６で例示する異常予兆の検知を実行する。具体的には、図６に示す石炭火力発電設備の目的変数１〜３および説明変数群を合わせた１６個のセンシングデータの時系列データに対して、図２６で例示する処理を実施する。図１８で例示した方法を用いてまずは無次元化を行う。但し、ここでは、実施例３と同様に、標準偏差で割り算する方法ではなく平均値で割り算する方法を用いる。

図２６の処理と図２０の処理との違いは、逐次的処理をしていく中で、一定の長期の時間を超えたら強制的に再度モデリングを行う点にある。具体的には、ステップＳ４１の実行後かつステップＳ４２の実行前に、異常判定部３３が、現時点が強制更新時間以内であるか否かを判定する（ステップＳ５２）点が異なる。ステップＳ５２で「Ｙｅｓ」と判定された場合にはステップＳ４２が実行される。ステップ５２で「Ｎｏ」と判定された場合にはステップＳ４５が実行される。ステップＳ５２を実行するのは、外れ値としては問題ない値で推移したとしても、一定時間が経過すると外れ値に寄与するデータが明確に変化するため、それを修正することがより好適であるためである。これを強制更新時間とした。強制更新時間は推定有効時間よりも数倍以上の長い時間で設定する。

本実施例では、実施例３と同様に許容値１はモデル更新間の瞬間値の平均値＋瞬間値の標準偏差の３倍（３σ）値とし、許容値２は２点平均値の平均値＋瞬間値の標準偏差の３倍（３σ）値とした。許容値３は瞬間値の平均値＋瞬間値の標準偏差の８倍（８σ）値とした。推定有効時間は５分間、強制更新時間４０分とした。なお、各センシングデータは２分毎に収集している。

図２７に結果を例示する。図２７において、システムが停止した１８:１８よりも１時間近く早い１７:２４の時点で異常が判定されている。また、その後一旦平常に戻り、再度１８：１０に異常が頻発している。そこで、各センシングデータを、同一日の０時〜２時の各平均値を用いて規格化したものを図２８に例示する。このように固定値で規格化したのは、モデリング等を行っていないゆえに、同じ指標で比較したかったためである。

図２８を参照すると、１７：２２から給炭量が減少し、それに伴い複数のセンシングデータが異なる推移を示し始めていることがわかる。現在の推定有効時間が５分であるため、１７：２７以降には異常が確定する可能性が高いが、給炭量を減らし始めたその次の測定の１７:２４の時点で異常となった。これは、給炭量を含むいくつかのデータはさらに減少が進んだが、その他のデータは明確な変化がなかったために平均値＋瞬間値の標準偏差の８倍（８σ）の値を超過したためである。給炭量の減少が緩和されて一旦落ち着いたが、再度１８：１０から給炭量を含む多くのパラメータが変化し始めるが、変化が急峻なので、その時点ですぐに判定結果は平均値＋瞬間値の標準偏差の８倍（８σ）の値を超過したため、異常判定がなされている。

以上から、本発明を用いた予兆検知によれば、温度だけでなく、様々な次元をもつセンシングデータ同士の相関関係から、一定程度あいまいさを含む閾値を設定しても即応性が高くて正確な異常の予兆検知を行うことが可能であることがわかる。なお、図２６のステップＳ４２およびステップＳ４７において、「ａｎｄ」判定を「ｏｒ」判定としてもよい。

以上、本発明の実施例について詳述したが、本発明は係る特定の実施例に限定されるものではなく、特許請求の範囲に記載された本発明の要旨の範囲内において、種々の変形・変更が可能である。例えば、上記各例においては推定有効時間に対する閾値を１つとしたが、当該第１の閾値よりも大きい第２の閾値をさらに設定してもよい。この場合、推定有効時間が第１の閾値よりも長く第２の閾値未満である場合に、異常警報よりも安全度が１段階高いが注意を促すものとして注意情報警報を出力するようにしてもよい。

１０ 説明変数取得部
１１ レーザ
１２ ビームスプリッタ
１３ 光スイッチ
１４ フィルタ
１５ａ，１５ｂ 検出器
１６ 演算部
２０ 温度センサ
２１ センサ部
２２ 測定器
３０ 判定部
３１ モデル作成部
３２ 閾値設定部
３３ 異常判定部
３４ 出力部
４０ ミル中間ハウジング部
４１ 石炭
４２ 貯留部
４３ 粉砕リング
４４ ローラ
４５ 微粉炭
４６ 一次空気室
１００ 判定装置

Claims (17)

1. 複数の第１センサ検出値の基準モデルを作成するモデル作成部と、
   所定の時点から、第２センサ検出値と作成された前記基準モデルとから求められた推定誤差が第１閾値を超過するか否かを判定し、前記第２センサ検出値と作成された前記基準モデルとから求められた前記推定誤差の積算値が第２閾値を超過するか否かを判定し、前記推定誤差が前記第１閾値を超えるまでの第１時間が第１の所定時間よりも短いか否かを判定し、前記推定誤差の積算値が前記第２閾値を超えるまでの第２時間が第２の所定時間より短いか否かを判定する判定部と、
   前記判定部が、前記推定誤差が前記第１閾値を超過したと判定し、且つ、前記第１時間が前記第１の所定時間より短いと判定した場合、及び前記積算値が前記第２閾値を超過したと判定し、且つ、前記第２時間が前記第２の所定時間より短いと判定した場合、異常に係る信号を出力する出力部と、を備えることを特徴とする判定装置。
2. 前記判定部が、前記推定誤差が前記第１閾値を超過したと判定し、前記積算値が前記第２閾値を超過したと判定し、前記第１時間が前記第１の所定時間よりも短いと判定し、且つ、前記第２時間が前記第２の所定時間より短いと判定した場合、前記出力部は、異常に係る信号を出力することを特徴とする請求項１に記載の判定装置。
3. 前記モデル作成部は、前記複数の第１センサ検出値を用いて回帰分析により前記基準モデルを作成することを特徴とする請求項１又は２に記載の判定装置。
4. 前記判定部が、前記推定誤差が前記第１閾値を超過したと判定し、且つ前記第１時間が前記第１の所定時間より短いと判定した場合、及び前記積算値が前記第２閾値を超過したと判定し、且つ、前記第２時間が前記第２の所定時間より短いと判定した場合、前記モデル作成部は、その時点から一定の過去時間における前記複数の第１センサ検出値を用いて前記基準モデルを再作成することを特徴とする請求項１〜３のいずれか一項に記載の判定装置。
5. 前記複数の第１センサ検出値及び／又は前記第２センサ検出値は、同一光ファイバの異なる長さ位置の後方散乱光によって得られる結果であることを特徴とする請求項１〜４のいずれか一項に記載の判定装置。
6. 前記第１閾値は、前記基準モデルの作成後の一定期間の推定誤差を基に決定され、前記第２閾値は、前記基準モデルの作成後の一定期間の推定誤差の積算値を基に決定されることを特徴とする請求項１〜５のいずれか一項に記載の判定装置。
7. 前記判定部は、更に前記第２センサ検出値と作成された前記基準モデルとから求められた推定誤差が前記第１閾値及び前記第２閾値よりも大きい第３閾値を超えたことを判定し、超えたと判定された場合は、前記出力部は、異常に係る信号を出力することを特徴とする請求項１〜６のいずれか一項に記載の判定装置。
8. 前記判定部が更に新しい第２検出値を取得した時点が強制更新時間以内かどうかを判定し、強制更新時間以内でないと判断した場合は、前記モデル作成部は、前記時点からの過去データを用いて再度基準モデルを作成し、前記時点の第２検出値を用いて、前記判定部が異常判定を行うことを特徴とする請求項１〜７のいずれか一項に記載の判定装置。
9. 複数の第１センサ検出値を用いて基準モデルを作成するモデル作成部と、
   所定の時点から、複数の第２センサ検出値と作成された前記基準モデルとから求められたマハラノビス平方距離が第１閾値を超過するか否かを判定し、複数の前記第２センサ検出値と作成された前記基準モデルとから求められたマハラノビス平方距離の平均値が第２閾値を超過するか否かを判定し、前記マハラノビス平方距離が前記第１閾値及び前記第２閾値よりも大きい第３閾値を超過したか否かを判定し、前記マハラノビス平方距離が前記第１閾値を超えるまでの第１時間が第１の所定時間よりも短いか否かを判定し、前記マハラノビス平方距離の平均値が前記第２閾値を超えるまでの第２時間が第２の所定時間より短いか否かを判定する判定部と、
   前記判定部が、前記マハラノビス平方距離が前記第１閾値と前記第３閾値とを超過したと判定し、且つ、前記第１時間が前記第１の所定時間より短いと判定した場合、及び前記平均値が前記第２閾値を超過したと判定し、前記マハラノビス平方距離が前記第３閾値を超過したと判定し、且つ、前記第２時間が前記第２の所定時間より短いと判定した場合、異常に係る信号を出力する出力部と、を備えることを特徴とする判定装置。
10. 前記出力部は、前記判定部が、前記マハラノビス平方距離が前記第１閾値を超過したと判定し、前記平均値が前記第２閾値を超過したと判定し、前記第１時間が前記第１の所定時間よりも短いと判定し、且つ、前記第２時間が前記第２の所定時間より短いと判定した場合、異常に係る信号を出力することを特徴とする請求項９に記載の判定装置。
11. 前記モデル作成部は、前記複数の第１センサ検出値の平均値及び分散共分散行列を用いて前記基準モデルを作成することを特徴とする請求項９または１０に記載の判定装置。
12. 前記判定部が、更に新しい第２センサ検出値を取得した時点が強制更新時間以内かどうかを判定し、強制更新時間以内でないと判断した場合、前記モデル作成部は、前記時点からの過去データを用いて再度基準モデルを作成し、前記時点の第２センサ検出値を用いて、前記判定部は再度異常判定を行うことを特徴とする請求項９〜１１のいずれか一項に記載の判定装置。
13. 前記複数の第１センサ検出値及び／又は前記複数の第２センサ検出値は、同一光ファイバの異なる長さ位置の後方散乱光によって得られる結果であることを特徴とする請求項９〜１２のいずれか一項に記載の判定装置。
14. 複数の第１センサ検出値の基準モデルを作成し、
    所定の時点から、第２センサ検出値と作成された前記基準モデルとから求められた推定誤差が第１閾値を超過するか否かを判定し、
    前記第２センサ検出値と作成された前記基準モデルとから求められた推定誤差の積算値が第２閾値を超過するか否かを判定し、
    前記推定誤差が前記第１閾値を超えるまでの第１時間が第１の所定時間よりも短いか否かを判定し、
    前記推定誤差の積算値が前記第２閾値を超えるまでの第２時間が第２の所定時間より短いか否かを判定し、
    前記推定誤差が前記第１閾値を超過したと判定され、且つ、前記第１時間が前記第１の所定時間より短いと判定された場合、及び前記積算値が前記第２閾値を超過したと判定され、且つ、前記第２時間が前記第２の所定時間より短いと判定された場合、異常に係る信号を出力する、ことを特徴とする判定方法。
15. コンピューターに、
    複数の第１センサ検出値の基準モデルを作成する処理と、
    所定の時点から、第２センサ検出値と作成された前記基準モデルとから求められた推定誤差が第１閾値を超過するか否かを判定し、前記第２センサ検出値と作成された前記基準モデルとから求められた前記推定誤差の積算値が第２閾値を超過するか否かを判定し、前記推定誤差が前記第１閾値を超えるまでの第１時間が第１の所定時間よりも短いか否かを判定し、前記推定誤差の積算値が前記第２閾値を超えるまでの第２時間が第２の所定時間より短いか否かを判定する処理と、
    前記判定する処理において、前記推定誤差が前記第１閾値を超過したと判定し、且つ、前記第１時間が前記第１の所定時間より短いと判定した場合、及び前記積算値が前記第２閾値を超過したと判定し、且つ、前記第２時間が前記第２の所定時間より短いと判定した場合、異常に係る信号を出力する処理と、を実行させることを特徴とする判定プログラム。
16. 複数の第１センサ検出値の基準モデルを作成し、
    所定の時点から、複数の第２センサ検出値と作成された前記基準モデルとから求められたマハラノビス平方距離が第１閾値を超過するか否かを判定し、
    複数の前記第２センサ検出値と作成された前記基準モデルとから求められたマハラノビス平方距離の平均値が第２閾値を超過するか否かを判定し、
    前記マハラノビス平方距離が前記第１閾値及び前記第２閾値よりも大きい第３閾値を超えたことを判定し、
    前記マハラノビス平方距離が前記第１閾値を超えるまでの第１時間が第１の所定時間よりも短いか否かを判定し、
    前記マハラノビス平方距離の平均値が前記第２閾値を超えるまでの第２時間が第２の所定時間より短いか否かを判定し、
    前記マハラノビス平方距離が前記第１閾値と前記第３閾値とを超過したと判定され、且つ、前記第１時間が前記第１の所定時間よりも短いと判定された場合、及び前記平均値が前記第２閾値を超過したと判定され、前記マハラノビス平方距離が前記第３閾値を超過したと判定され、且つ、前記第２時間が前記第２の所定時間より短いと判定された場合、異常に係る信号を出力することを特徴とする判定方法。
17. コンピューターに、
    複数の第１センサ検出値の基準モデルを作成する処理と、
    所定の時点から、複数の第２センサ検出値と作成された前記基準モデルとから求められたマハラノビス平方距離が第１閾値を超過するか否かを判定し、複数の前記第２センサ検出値と作成された前記基準モデルとから求められたマハラノビス平方距離の平均値が第２閾値を超過するか否かを判定し、前記マハラノビス平方距離が前記第１閾値及び前記第２閾値よりも大きい第３閾値を超過したか否かを判定し、前記マハラノビス平方距離が前記第１閾値を超えるまでの第１時間が第１の所定時間よりも短いか否かを判定し、前記マハラノビス平方距離の平均値が前記第２閾値を超えるまでの第２時間が第２の所定時間より短いか否かを判定する処理と、
    前記マハラノビス平方距離が前記第１閾値と前記第３閾値とを超過したと判定され、且つ、前記第１時間が前記第１の所定時間よりも短いと判定された場合、及び前記平均値が前記第２閾値を超過したと判定され、前記マハラノビス平方距離が前記第３閾値を超過したと判定され、且つ、前記第２時間が前記第２の所定時間より短いと判定された場合、異常に係る信号を出力する処理と、を実行させることを特徴とする判定プログラム。

Images