(57)【要約】【課題】時系列データを高精度にクラス分類するのに有効なクラス分類モデルのモデルパラメータ及び特徴波形を学習する。【解決手段】本発明の実施形態としての時系列データ分析方法は、第１クラスに属する時系列データである複数の第１時系列データと、複数の特徴波形との距離に基づき、前記複数の特徴波形の特徴量を含む複数の第１特徴ベクトルを生成し、第２クラスに属する時系列データである複数の第２時系列データと、前記複数の特徴波形との距離に基づき、前記複数の特徴波形の特徴量を含む複数の第２特徴ベクトルを生成する、特徴ベクトル生成ステップと、前記複数の第１特徴ベクトルと、前記複数の第２特徴ベクトルと、クラス分類モデルの性能指標に関するパラメータである性能指標パラメータとに基づき、前記複数の特徴波形の重みを含むモデルパラメータと、前記複数の特徴波形とを更新する更新処理ステップとを備える。【選択図】図１
2020170214.tif 000002
【特許請求の範囲】
【請求項１】
第１クラスに属する時系列データである複数の第１時系列データと、複数の特徴波形との距離に基づき、前記複数の特徴波形の特徴量を含む複数の第１特徴ベクトルを生成し、第２クラスに属する時系列データである複数の第２時系列データと、前記複数の特徴波形との距離に基づき、前記複数の特徴波形の特徴量を含む複数の第２特徴ベクトルを生成する、特徴ベクトル生成ステップと、
前記複数の第１特徴ベクトルと、前記複数の第２特徴ベクトルと、クラス分類モデルの性能指標に関するパラメータである性能指標パラメータとに基づき、前記複数の特徴波形の重みを含むモデルパラメータと、前記複数の特徴波形とを更新する更新処理ステップと
を備えた時系列データ分析方法。
【請求項２】
前記特徴ベクトル生成ステップは、テスト対象となる時系列データと、更新された前記複数の特徴波形との距離に基づき、更新された前記複数の特徴波形の特徴量を含む第３特徴ベクトルを生成し、
前記第３特徴ベクトルと、更新された前記モデルパラメータとに基づいて、前記テスト対象となる時系列データのクラスを前記第１クラス及び前記第２クラスのいずれかに予測する予測ステップを備えた
請求項１に記載の時系列データ分析方法。
【請求項３】
前記予測ステップは、前記第３特徴ベクトルと、更新された前記モデルパラメータとに基づきスコアを算出し、前記スコアを閾値と比較することにより、前記テスト対象となる時系列データのクラスを予測する
請求項２に記載の時系列データ分析方法。
【請求項４】
前記更新処理ステップは、
前記複数の第１時系列データに対して、前記複数の第２特徴ベクトルと前記モデルパラメータとに基づき、前記第１クラスに属する可能性の大きさに関する複数のスコアを算出し、
前記複数の第２時系列データに対して、前記複数の第２特徴ベクトルと前記モデルパラメータとに基づき、前記第１クラスに属する可能性の大きさに関する複数のスコアを算出し、
前記複数の第１時系列データから少なくとも１つの前記第１時系列データを選択し、
前記第１時系列データのスコアよりも前記第２時系列データのスコアの方が前記第１クラスに属する可能性が大きくなる、前記第１時系列データと前記第２時系列データの組み合わせを特定し、
特定した組み合わせと、選択した少なくとも１つの前記第１時系列データの前記第１特徴ベクトルに基づき、前記複数の特徴波形を更新する
請求項１〜３のいずれか一項に記載の時系列データ分析方法。
【請求項５】
前記更新処理ステップは、少なくとも１つの前記第１時系列データをランダムに選択する
請求項４に記載の時系列データ分析方法。
【請求項６】
前記更新処理ステップは、
前記複数の第１時系列データに対して、前記複数の第１特徴ベクトルと前記モデルパラメータとに基づき、前記第１クラスに属する可能性の大きさに関する複数のスコアを算出し、
前記複数の第２時系列データに対して、前記複数の第２特徴ベクトルと前記モデルパラメータとに基づき、前記第１クラスに属する可能性の大きさに関する複数のスコアを算出し、
前記複数の第２時系列データの前記複数のスコアに基づき、少なくとも１つの前記第２時系列データを選択し、
前記第１時系列データのスコアよりも前記第２時系列データのスコアの方が前記第１クラスに属する可能性が大きくなる、前記第１時系列データと前記第２時系列データの組み合わせを特定し、
特定した組み合わせと、選択した少なくとも１つの前記第２時系列データの前記第２特徴ベクトルに基づき、前記複数の特徴波形を更新する
請求項１〜３のいずれか一項に記載の時系列データ分析方法。
【請求項７】
前記更新処理ステップは、前記複数の第２時系列データの中から、前記スコアが大きい上位の所定数の前記第２時系列データを選択し、前記所定数の第２時系列データから少なくとも１つの前記第２時系列データを選択する
請求項６に記載の時系列データ分析方法。
【請求項８】
前記更新処理ステップは、
前記複数の第２時系列データに対して、前記複数の第２特徴ベクトルと前記モデルパラメータとに基づき、前記第１クラスに属する可能性の大きさに関する複数のスコアを算出し、
前記複数の第２時系列データの前記スコアに基づき少なくとも１つの前記第２特徴ベクトルを選択し、選択した第２特徴ベクトルと、前記複数の第１特徴ベクトルとに基づき、前記モデルパラメータを更新する
請求項１〜７のいずれか一項に記載の時系列データ分析方法。
【請求項９】
前記更新処理ステップは、前記モデルパラメータに関するＬ１ノルム正則化付の目的関数の最小化を行うことにより、前記複数の特徴波形を更新し、
前記目的関数は、前記重みを表す変数と、前記複数の特徴波形を表す変数とを含む
請求項１〜８のいずれか一項に記載の時系列データ分析方法。
【請求項１０】
前記性能指標パラメータの値は調整可能であり、
前記性能指標パラメータの値を指定する情報を受信し、前記情報に基づき前記性能指標パラメータを設定する入力設定ステップ
を備えた請求項１〜９のいずれか一項に記載の時系列データ分析方法。
【請求項１１】
前記性能指標は、ｐＡＵＣ（partial Area Under the Curve）であり、
前記性能指標パラメータは、偽陽性率の範囲に関するパラメータである
請求項１〜１０のいずれか一項に記載の時系列データ分析方法。
【請求項１２】
前記クラス分類モデルは、サポートベクタマシン（ＳＶＭ）モデルである
請求項１〜１１のいずれか一項に記載の時系列データ分析方法。
【請求項１３】
前記特徴量は、前記特徴波形と前記時系列データとの最小距離である
請求項１〜１２のいずれか一項に記載の時系列データ分析方法。
【請求項１４】
前記第１クラスは正常クラス、前記第２クラスは異常クラスである
請求項１〜１３のいずれか一項に記載の時系列データ分析方法。
【請求項１５】
前記更新処理ステップにより更新された前記特徴波形を表示する表示ステップ
を備えた請求項１〜１４のいずれか一項に記載の時系列データ分析方法。
【請求項１６】
第１クラスに属する時系列データである複数の第１時系列データと、複数の特徴波形との距離に基づき、前記複数の特徴波形の特徴量を含む複数の第１特徴ベクトルを生成し、第２クラスに属する時系列データである複数の第２時系列データと、前記複数の特徴波形との距離に基づき、前記複数の特徴波形の特徴量を含む複数の第２特徴ベクトルを生成する、特徴ベクトル生成部と、
前記複数の第１特徴ベクトルと、前記複数の第２特徴ベクトルと、クラス分類モデルの性能指標に関するパラメータである性能指標パラメータとに基づき、前記複数の特徴波形の重みを含むモデルパラメータと、前記複数の特徴波形とを更新する更新処理部と
を備えた時系列データ分析装置。
【請求項１７】
第１クラスに属する時系列データである複数の第１時系列データと、複数の特徴波形との距離に基づき、前記複数の特徴波形の特徴量を含む複数の第１特徴ベクトルを生成し、第２クラスに属する時系列データである複数の第２時系列データと、前記複数の特徴波形との距離に基づき、前記複数の特徴波形の特徴量を含む複数の第２特徴ベクトルを生成する、特徴ベクトル生成ステップと、
前記複数の第１特徴ベクトルと、前記複数の第２特徴ベクトルと、クラス分類モデルの性能指標に関するパラメータである性能指標パラメータとに基づき、前記複数の特徴波形の重みを含むモデルパラメータと、前記複数の特徴波形とを更新する更新処理ステップと
をコンピュータに実行させるためのコンピュータプログラム。
【請求項１８】
クラス分類モデルの性能指標に関するパラメータである性能指標パラメータをユーザ入力情報に基づき調整するステップと、
調整された前記性能指標パラメータと、第１クラスに属する時系列データである複数の第１時系列データと、第２クラスに属する時系列データである複数の第２時系列データとに基づき、時系列データのクラス分類モデルにおける複数の特徴波形の重みを含むモデルパラメータと、前記複数の特徴波形とを更新する更新処理ステップと、
更新された前記複数の特徴波形と、更新された前記モデルパラメータとに基づき、前記第１時系列データ、前記第２時系列データ、又はテスト対象となる時系列データのスコアを算出するスコア算出ステップと、
更新された前記複数の特徴波形と、前記スコアに基づく情報とを表示する表示ステップと
を備えた時系列データ分析方法。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明の実施形態は、時系列データ分析方法、時系列データ分析装置及びコンピュータプログラムに関する。
【背景技術】
【０００２】
センサデータなどの時系列データを正常及び異常の２クラスのいずれかに分類する２クラス分類では、予測の性能に加えて、予測の根拠を明確にする必要もある。
【０００３】
時系列データに対して予測の根拠を明確にする２クラス分類技術として、クラス分類モデルに加えて予測に有効な部分波形パターンである特徴波形（shapelets）を同時に学習する方法が近年盛んに研究されている。この方法は、特徴波形学習法と呼ばれている。特徴波形は、異常有無の予測に用いられるとともに、予測の根拠として用いられる。
【０００４】
一方、クラス分類性能や、ランキングなどの予測性能指標としてＡＵＣ（Area Under the Curve）がよく用いられる。ＡＵＣは横軸を偽陽性率（false positive rate）、縦軸を真陽性率（true positive rate）とした曲線であるＲＯＣ曲線（ＲＯＣ（Receiver Operating Characteristics）curve）と、横軸との間の面積である。
【０００５】
偽陽性率を狭い範囲に絞り、この範囲内での面積であるｐＡＵＣ（partial AUC）を評価指標として用いる場合も多い。例えば、トラブルの見逃しを十分低く保ったもとで、正しくトラブルを予測したい、誤診断を低く抑えたもとで、正しく診断事例を予測したい、少数の上位のランキングの予測精度を正確にしたい等の場面では、ｐＡＵＣが性能指標として用いられる。
【０００６】
しかし、時系列データのクラス分類問題に対して、ｐＡＵＣ又はＡＵＣといった性能指標を最適化するようにクラス分類モデルを学習する方法は無かった。また、性能指標を最適化すると同時に、予測の根拠となる特徴波形を同時に学習する方法も無かった。
【先行技術文献】
【非特許文献】
【０００７】
【非特許文献１】KDD '14 Proceedings of the 20th ACM SIGKDD international conference on Knowledge discovery and data mining， Pages 392-401／Josif Grabocka et al.
【非特許文献２】KDD '13 Proceedings of the 19th ACM SIGKDD international conference on Knowledge discovery and data mining, SVMpAUCtight: a new support vector method for optimizing partial AUC based on a tight convex upper bound／Harikrishna Narasimhan et al.
【非特許文献３】Efficient Projections onto the l1-Ball for Learning in High Dimensions, ICML2008
【発明の概要】
【発明が解決しようとする課題】
【０００８】
本発明の実施形態は、時系列データを高精度にクラス分類するのに有効なクラス分類モデルのモデルパラメータ及び特徴波形を学習することが可能な時系列データ分析方法、時系列データ分析装置及びコンピュータプログラムを提供する。
【課題を解決するための手段】
【０００９】
本発明の実施形態としての時系列データ分析方法は、第１クラスに属する時系列データである複数の第１時系列データと、複数の特徴波形との距離に基づき、前記複数の特徴波形の特徴量を含む複数の第１特徴ベクトルを生成し、第２クラスに属する時系列データである複数の第２時系列データと、前記複数の特徴波形との距離に基づき、前記複数の特徴波形の特徴量を含む複数の第２特徴ベクトルを生成する、特徴ベクトル生成ステップと、前記複数の第１特徴ベクトルと、前記複数の第２特徴ベクトルと、クラス分類モデルの性能指標に関するパラメータである性能指標パラメータとに基づき、前記複数の特徴波形の重みを含むモデルパラメータと、前記複数の特徴波形とを更新する更新処理ステップとを備える。
【図面の簡単な説明】
【００１０】
【図１】第１の実施形態に係る時系列データ分析装置を表すブロック図。
【図２】学習用データ記憶部に格納されている時系列データ集合の例を示す図。
【図３】２つの特徴波形を含む特徴波形集合の例を示す図。
【図４】ＲＯＣ曲線、ＡＵＣ、ｐＡＵＣを説明するための図。
【図５】ＧＵＩの画面の例を示す図。
【図６】Ｌ２ボールに重みベクトルを射影する例を示す図。
【図７】特徴波形更新部の動作の一例を示すフローチャート。
【図８】表示部に表示された学習結果データの例を示す図。
【図９】ＲＯＣ曲線の例を示す図。
【図１０】ＧＵＩの画面の例を示す図。
【図１１】学習用の時系列データと、学習された特徴波形との表示例を示す図。
【図１２】ＲＯＣ曲線の例を示す図。
【図１３】学習されたモデルパラメータ（重みベクトル）により表される識別境界を模式的に示す図。
【図１４】学習フェーズの動作のフローチャート。
【図１５】表示部に表示された評価結果データの一例を示す図。
【図１６】テストフェーズの動作のフローチャート。
【図１７】本実施形態に係る時系列データ分析装置のハードウェア構成を示す図。
【図１８】第２の実施形態に係る時系列データ解析装置のブロック図。
【図１９】Ｌ１ボールに重みベクトルを射影する例を示す図。
【図２０】第２の実施形態に係る学習フェーズの動作のフローチャート。
【図２１】第３の実施形態に係る時系列データ分析システムを示す図。
【発明を実施するための形態】
【００１１】
以下、図面を参照しながら、本発明の実施形態について説明する。
（第１実施形態）
図１は、本発明の実施形態に係る時系列データ分析装置を表すブロック図である。図１の時系列データ分析装置は、入力設定部１０、特徴ベクトル生成部１１、更新処理部１２、予測部１６、表示部１７、学習用データ記憶部１８、テスト用データ記憶部１９、及び出力情報記憶部２０を備える。更新処理部１２は、重み更新部１３、特徴波形更新部１４、及び更新終了判定部１５を備える。
【００１２】
本時系列データ分析装置は、学習フェーズと、テストフェーズとを備える。学習フェーズでは、学習用の時系列データと、クラス分類モデルの性能指標に関するパラメータである性能指標パラメータとに基づき、クラス分類モデルのモデルパラメータと、複数の特徴波形とを学習する。モデルパラメータは複数の特徴波形に対する重みを含む。テストフェーズでは、学習フェーズで学習したモデルパラメータと複数の特徴波形とを用いて、テスト対象となる時系列データのクラスを予測することにより、当該テスト対象となる時系列データに異常があるかを判断する。
【００１３】
以下、学習フェーズとテストフェーズに分けて、本装置について詳細に説明する。
【００１４】
＜学習フェーズ＞
【００１５】
学習用データ記憶部１８は、学習用の入力データを記憶している。学習用の入力データとして、２値ラベル付きの学習用の時系列データのセット、時系列データに関するパラメータ情報、特徴波形（shapelet）に関するパラメータ情報、クラス分類モデルの性能指標の情報、当該性能指標に関するパラメータ（性能指標パラメータ）の情報とを記憶している。
【００１６】
学習用データ記憶部１８は、２値ラベル付きの学習用の時系列データのセットと、時系列データのパラメータ情報を記憶している。時系列データは、一例として、分析対象装置に設置されたセンサの検出値に基づく時系列データである。時系列データは、センサの検出値そのものでもよいし、検出値の統計値（平均、最大、最小、標準偏差など）でもよいし、複数のセンサの検出値の演算値（例えば電流と電圧とを乗算した電力）でもよい。
【００１７】
学習用の時系列データは、教師付き時系列データであり、正常又は異常を表す２値ラベルが付与されている。一例として、正常ラベルのラベル値は１、異常ラベルのラベル値は−１であるが、これに限定されない。正常ラベルを付与された時系列データ（第１時系列データ）は正常クラスに属する。異常ラベルを付与された時系列データは異常クラスに属する。
【００１８】
一例として正常クラスは第１クラス、異常クラスは第２クラスに対応する。第１クラスに属する時系列データは第１時系列データに対応する。第２クラスに属する時系列データは第２時系列データに対応する。
【００１９】
時系列データのパラメータ情報は、時系列データの個数、及び時系列データの長さの情報である。以下の説明で、時系列データの集合を“Ｔ”、時系列データの個数を“Ｉ”とする。また、各時系列データの長さを“Ｑ”とする。すなわち、各時系列データは、Ｑ個の点からなるデータである。時系列データセットＴは、Ｉ×Ｑの行列で表すことができる。
【００２０】
図２に、学習用データ記憶部１８に格納されている時系列データ集合Ｔの例を示す。集合ＴにはＩ個の時系列データが含まれる。各時系列データの長さは同じＱである。すなわち、各時系列データは、Ｑ個の点を含む。図では、Ｑ個の点を線でつないだ例が示される。個々の時系列データをＴｉ（ｉ＝１，２，…,Ｉ）によって表している。任意の時系列データは時系列データｉと表現する。本実施形態では、各時系列データの長さは同じＱであるが、長さが異なる場合への拡張も可能である。各時系列データは等間隔でサンプリングされており、データ欠損はないものとする。データ欠損がある場合は、補間処理によりデータ補間を行えばよい。
【００２１】
また、学習用データ記憶部１８は、特徴波形（shapelet）に関するパラメータ情報として、特徴波形の個数と、特徴波形の長さを表す値を記憶している。特徴波形の個数を“Ｋ”、特徴波形の長さを“Ｌ”とする。Ｌは、時系列データの長さＱよりも小さい値である。
【００２２】
特徴波形は、Ｌ個の点からなるデータである。特徴波形の集合をＳとすると、ＳはＫ×Ｌの行列である。特徴波形は、Ｔｉｍｅ  Ｓｅｒｉｅｓ  Ｓｈａｐｅｌｅｔｓ法（ＴＳＳ法）でｓｈａｐｅｌｅｔと呼ばれるものに相当する。後述するように、特徴波形は、学習フェーズの開始時に初期の形状が決定された後は、繰り返し更新されていく。
【００２３】
図３に、２つ（Ｋ＝２）の特徴波形を含む特徴波形集合Ｓの例を示す。各特徴波形の長さはＬである。各特徴波形をＳ１、Ｓ２で表している。本実施形態では、各特徴波形の長さは同じＬであるが、長さが異なる場合への拡張も可能である。
【００２４】
ここで、時系列データｉと特徴波形ｋとの距離について説明する。時系列データｉと特徴波形ｋとの距離は、時系列データｉにおける長さＫの各区間の部分時系列と、特徴波形ｋとの距離のうち、最も小さい距離として定義される。具体的には、時系列データｉにおいて波形の開始位置（先頭）からの長さであるオフセットを波形の末尾方向に順次動かす。各オフセットから長さＬの区間の部分時系列と、特徴波形ｋとの距離を計算する。そして、最も小さくなる距離を、時系列データｉと特徴波形ｋとの距離とする。距離が小さいほど、特徴波形ｋは時系列データと類似（フィット）している。距離にはユークリッド距離を用いる。但し、波形間の類似度を評価可能な距離であれば、どのような種類の距離でもよい。
【００２５】
時系列データｉにおけるオフセットｊから長さＬの区間の部分時系列と、特徴波形ｋとの距離は、以下の式（１）で計算される。
【００２６】
【数１】
2020170214.tif 000003
【００２７】
Ｔｉ，ｊ＋ｌ−１は、時系列データ集合Ｔに含まれる時系列データｉにおけるオフセットｊの位置から数えてｌ−１番目の位置の値を表す。Ｓｋ，ｌは、特徴波形集合Ｓに含まれる特徴波形ｋの先頭から数えてｌ番目の位置の値を表す。Ｄｉ，ｋ，ｊは、時系列データｉにおけるオフセットｊから長さＬの区間の部分時系列（部分波形）と、特徴波形ｋとの間の平均距離に相当する。
【００２８】
時系列データｉと特徴波形ｋとの距離は、上記の式（１）に基づき、以下の式（２）で計算される。
【００２９】
【数２】
2020170214.tif 000004
【００３０】
本学習フェーズでは、クラス分類モデルのモデルパラメータの学習と、特徴波形集合Ｓの学習を行う。クラス分類モデルとして、Ｓｕｐｐｏｒｔ  Ｖｅｃｔｏｒ  Ｍａｃｈｉｎｅ（ＳＶＭ）モデルを想定する。この場合、モデルパラメータは、識別境界の重みベクトルＷに対応する。重みベクトルＷはＫ次元のベクトルであり、Ｋ個の特徴波形の重みを含む。特徴波形集合ＳはＫ×Ｌの行列である。前述した通り、Ｋは特徴波形の個数を表し、Ｌは特徴波形の長さを表す。
【００３１】
学習用データ記憶部１８は、学習においてクラス分類モデルの性能を評価するための性能指標と、及び性能指標に関するパラメータ（性能指標パラメータ）とを記憶する。ここでは、性能指標の例としてｐＡＵＣ（partial Area Under the ROC curve）、性能指標パラメータとして偽陽性率（false positive rate）の範囲を指定するパラメータを用いる。偽陽性率とは、正解が負ラベルのデータのラベルを、間違えて正ラベルと予測する割合である。すなわち、異常クラスに属するデータのクラスを、正常クラスと予測する割合である。
【００３２】
ここで図４を用いて、ＲＯＣ曲線、ＡＵＣ、ｐＡＵＣについて説明する。
【００３３】
図４は、ＲＯＣ曲線、ＡＵＣ、ｐＡＵＣを模式的に示す。ＲＯＣ曲線は、縦軸を真陽性率（true positive rate）、横軸を偽陽性率（false positive rate）とした座標系に描かれたグラフである。真陽性率は、正解が正ラベルのデータを、正しく正ラベルと予測する割合である。すなわち、正常クラスに属するデータのクラスを、正常クラスと予測する割合である。予測は、クラス分類モデルの出力値（スコア）を閾値と比較することで行う。閾値以上であれば正常と判断（正ラベルを予測）、閾値未満であれば異常と判断（負ラベルを予測）する。閾値を例えばスコアの最大値から最小値まで変更させると、ＲＯＣ曲線が描かれる。
【００３４】
正解が負ラベルのデータを、間違えて正ラベルと予測した回数をＦＰ、正解が正ラベルのデータを、正しく正ラベルと予測した回数をＴＰ、正解が正ラベルのデータを、間違えて負ラベルと予測した回数をＦＮ、正解が負ラベルのデータを、正しく負ラベルと予測する回数をＴＮとする。このとき真陽性率はＴＰ／（ＴＰ＋ＦＮ）、偽陽性率はＦＰ／（ＦＰ＋ＴＮ）により計算できる。
【００３５】
ＡＵＣは、ＲＯＣ曲線の下の面積、すなわち、ＲＯＣ曲線と偽陽性率の軸とによって囲まれた領域の面積である。
【００３６】
ｐＡＵＣは、横軸の偽陽性率における特定の範囲と、ＲＯＣ曲線とによって囲まれる領域の面積である。横軸の偽陽性率の範囲は、０以上１以下であり、上記の特定の範囲は、性能指標パラメータで指定される。例えば、０以上０．１以下とする。但し、この範囲は一例であり、０以上０．０５以下、又は、０以上０．０１以下など、他の範囲でもよい。
【００３７】
ここでは性能指標としてｐＡＵＣ、性能指標パラメータとして偽陽性率の範囲を用いるが、他の性能指標でもよい。例えば、偽陰性率（false negative rate）に基づいた性能指標及び性能指標パラメータを定義し、それを用いてもよい。偽陰性率（false negative rate）は、正解が正ラベルのデータを、間違えて負ラベルと予測する割合である。
【００３８】
本装置の操作者（ユーザ）は、ＧＵＩ（Graphical User Interface）を介して、学習用データ記憶部１８に学習用の入力データを設定してもよい。この場合、ＧＵＩは表示部１７の画面に表示させる。
【００３９】
図５は、ＧＵＩの画面の例を示す。性能指標としてｐＡＵＣを指定している。性能指標パラメータとして調整つまみ（バー）３１の位置を調整することで、偽陽性率の範囲を指定している。偽陽性率の範囲の下限値はゼロとし、ここでは範囲の上限値（β）を指定している。図示の例では、上限値（β）を０．０１としている。なお、偽陽性率の範囲の下限値をゼロより大きい値とすることを排除しない。調整つまみ３１を左右に移動させることで、上限値β（偽陽性率の範囲）を容易に調整可能になっている。また、学習用の時系列データセットをファイルパスにより指定している。この場合、ファイルパスで指定されたファイルを読み出す。ファイルには時系列データの長さ及び個数の情報が格納されていてもよい。また、図の画面で、特徴波形の長さとして５、特徴波形の個数として２を指定している。
【００４０】
入力設定部１０は、学習用データ記憶部１８から学習用の入力データとして、２値ラベル付きの学習用の時系列データのセット、時系列データのパラメータ情報（時系列データの個数、時系列データの長さ）、特徴波形のパラメータ情報（特徴波形の個数、長さ）、クラス分類モデルの性能指標の情報（ここではｐＡＵＣ）、性能指標パラメータの情報（ここでは偽陽性率の範囲）を読み出す。入力設定部１０は、読み出したデータを特徴ベクトル生成部１１に入力する。学習用の入力データの一部又は全部を、本装置の操作者又は管理者であるユーザが、入力装置を用いて入力してもよい。入力装置は、キーボード、マウス、タッチパッド又はスマートフォンなど、本装置に各種データ又は指示を入力するための装置である。この場合、入力装置から受信したデータが、特徴ベクトル生成部１１に入力される。
【００４１】
ここで、特徴波形のパラメータ情報を特徴ベクトル生成部１１に入力しない構成も可能である。この場合、特徴ベクトル生成部１１では、特徴波形の個数及び長さとして、デフォルト値を用いればよい。例えば、特徴波形の個数（最大個数）Ｋを２、特徴波形の長さＬをＱ×０．１などとする。
【００４２】
また、入力設定部１０は、特徴波形集合Ｓと、モデルパラメータ（重みベクトル）Ｗとを初期化する。重みベクトルＷは、各特徴波形の重みを含む。
【００４３】
重みベクトルＷの初期化は、例えば、Ｋ個の全ての要素を０とする。
【００４４】
特徴波形集合Ｓの初期化は、例えば、以下のように行う。各時系列データの先頭から長さＬの窓フレームを一定間隔ずつシフトしながら、窓フレームに含まれる長さＬの部分波形（セグメント）を抽出する。これらのセグメントに対して、ｋ−ｍｅａｎｓ法などのクラスタリングを行うことにより、Ｋ個のクラスタを生成する。Ｋ個のクラスタのセントロイド（重心）を計算する。重心は、例えばクラスタに属するすべてのセグメントの平均を計算することで取得する。Ｋ個のクラスタから計算されたＫ個の重心を、初期化した特徴波形集合Ｓとする。
【００４５】
特徴波形集合Ｓの初期化は、他の方法でもよい。例えば乱数を用いて長さＬのＫ個の特徴波形を生成し、これらの特徴波形を特徴波形集合Ｓとすることも排除されない。
【００４６】
特徴ベクトル生成部１１が、各時系列データに対してＫ次元の特徴ベクトルを生成する。具体的には、各時系列データを対象として、対象となる時系列データとＫ個の特徴波形との距離を特徴量として計算する。計算したＫ個の特徴量（距離）を、Ｋ個の特徴波形に対応する要素に格納したＫ次元の特徴ベクトルを生成する。この特徴ベクトルは学習用の特徴ベクトルである。第１時系列データ（例えば正常クラスの時系列データ）に基づき生成された特徴ベクトルは第１特徴ベクトル、第２時系列データ（例えば異常クラスの時系列データ）に基づき生成された特徴ベクトルは第２特徴ベクトルに対応する。
【００４７】
ｉ番目の時系列データに対する特徴ベクトルを“Ｘｉ”とする。特徴ベクトルＸｉのｋ番目の要素（ｋ番目の特徴波形との距離）は、前述した式（２）で定義される通り、Ｘｉ，ｋである。したがって、例えば、特徴波形セットＳが特徴波形１，２，・・・，Ｋを含み、特徴波形１、２、・・・、Ｋと時系列データｉとの距離をＸｉ，１，Ｘｉ，２，・・・，Ｘｉ，ｋとすると、特徴ベクトルＸｉ＝（Ｘｉ，１，Ｘｉ，２，・・・，Ｘｉ，ｋ）である。
【００４８】
更新処理部１２は、複数の時系列データに対して生成された複数の特徴ベクトルと、性能指標パラメータとに基づき、クラス分類モデルにおけるモデルパラメータ（重みベクトル）Ｗと、複数の特徴波形（特徴波形セットＳ）とを更新する。以下、更新処理部１２が備える重み更新部１３、特徴波形更新部１４及び更新終了判定部１５について説明する。
【００４９】
重み更新部１３は、機械学習により、クラス分類モデルのモデルパラメータ（重みベクトル）の学習と、特徴波形集合の学習とを同時に行う。ここでは、クラス分類モデルとして、Ｓｕｐｐｏｒｔ  Ｖｅｃｔｏｒ  Ｍａｃｈｉｎｅ（ＳＶＭ）を用いる。ＳＶＭは、特徴空間において、正常と異常を判別する識別境界を学習するアルゴリズム、または当該識別境界に基づき判定を行うクラス分類モデルである。特徴空間は、Ｘｉ，ｋ（ｋ＝１，２，…，Ｋ）を軸とするＫ次元の空間である。特徴波形の個数Ｋが２であれば、特徴空間は、Ｘｉ，１とＸｉ，２を軸とする２次元の空間である。モデルパラメータ（重みベクトル）は、当該識別境界に対応する。モデルパラメータ（重みベクトル）Ｗは、各特徴波形に対応するパラメータ（重み）ｗ１，ｗ２，・・・，ｗｋを含む。識別境界は、線形とするが、非線形でもよい。
【００５０】
識別境界が非線形の場合、モデルパラメータ（重みベクトル）は無限次元のベクトルとなるため、代わりに、識別境界のモデルパラメータ（重みベクトル）Ｗとして、サポートベクトル集合Ｓｖと、集合Ｓｖに属するサポートベクトルの寄与率の集合Ｓａとを用いる。サポートベクトルは，識別境界の決定に寄与する特徴ベクトルである。寄与率は，そのサポートベクトルが、識別境界の決定にどの程度寄与するかを表しており，寄与率の絶対値が大きいほど、その決定に大きく寄与する（寄与率が０の場合は，識別境界の決定に寄与せず、それに対応する特徴ベクトルはサポートベクトルではない）。ＳＶＭでは，カーネル（内積を拡張した関数）と、サポートベクトルとその寄与率とを用いて，非線形の識別境界を表現できる。
【００５１】
本実施形態においてＳＶＭによるモデルパラメータ（重みベクトル）の学習を、特徴波形集合の学習と同時に行う。これらの学習は、以下の最適化問題として定式化される。
【００５２】
【数３】
2020170214.tif 000005
【００５３】
上記の最適化問題において、Ｘ＋は、正ラベルを持つ特徴ベクトルの集合を表す。Ｘ−は、負ラベルを持つ特徴ベクトルの集合を表す。正ラベルを持つ特徴ベクトルを正常特徴ベクトル、負ラベルをもつ特徴ベクトルを異常特徴ベクトルと称する場合がある。
【００５４】
式（３）は、関数Ｈ（Ｘ，ｚ，π，Ｗ）の最大値を、制約条件（Subject to）の下で最小化することにより、モデルパラメータ（重みベクトル）Ｗと、特徴波形集合Ｓとを求めることを規定している。
【００５５】
βは、性能指標パラメータとして指定されたｐＡＵＣの偽陽性率の範囲の上限値である。例えば偽陽性率の範囲が０以上０．１以下として指定された場合、βは０．１である。
【００５６】
Ｉ＋は、時系列データセットにおける正ラベルの時系列データ（正例）の個数を表す。Ｉ−は、時系列データセットにおける正ラベルの時系列データ（負例）の個数を表す。全時系列データの個数は前述した通りＩであるから、Ｉ＝Ｉ＋＋Ｉ−である。
【００５７】
ｚβは、負例の特徴ベクトル（異常特徴ベクトル）の集合から任意にｉβ−個（少なくとも１つ）を取り出した場合の異常特徴ベクトルの部分集合を表す。ｉβ−は、式（５）で定義されており、βＩ−を引数として床関数を計算した値である。
2020170214.tif 000006
は床関数であり、実数である引数以下の最大の整数である。例えば、引数が４．６５であれば、床関数の出力は４である。なお、βＩ−が整数にならない場合や、偽陽性率の範囲の下限値が０以外の値を取る場合は、文献（A Structural SVM Based Approach for Optimizing Partial AUC, JMLR2013）と同様の方法で厳密に定式化することもできる。
【００５８】
ｚは、式（６）で定義される通り、スコアの大きい順に選択したｉβ−個の異常特徴ベクトルを格納することを定義している。Ｘ−の添字であるａ１、ａ２、・・・はインデックスを表している。スコアはクラス分類モデルの出力値であり、後述するように、モデルパラメータ（重みベクトル）と特徴ベクトルの内積である。クラス分類モデルの出力値をＹとすると、クラス分類モデルはＹ＝ＷＴ・Ｘと表せる。本実施形態では、スコアが大きいほど正ラベルの可能性が高いことを意味する。但し、スコアが小さいほど正ラベルの可能性が高くなるように変形してもよい。
【００５９】
式（３）のΠは、時系列データセットから選択したＩ＋個の正例とｉβ−個の負例とに基づくオーダーリング行列（ｏｒｄｅｒｉｎｇ  ｍａｔｒｉｘ）の集合を表す。オーダーリング行列は、Ｉ＋行ｉβ−列の行列である。オーダーリング行列は、Ｉ＋個の正例と、ｉβ−個の負例とを１つずつ組み合わせたペア（計Ｉ＋×ｉβ−個のペア）に対応する要素に、正例のスコアと負例のスコアとの大小関係に応じた値を格納する。正例のスコアが負例のスコア以上であれば、対応する要素の値は０、正例のスコアが負例のスコアより小さければ、そのペアに対応する要素の値は１である。すなわち、各ペアについて、正例よりも正ラベルの可能性が高いスコアを負例が有する場合に１、それ以外の場合に０を格納する。換言すれば、スコアの大小関係が本来の関係と逆転している場合に１、それ以外は０を格納する。
【００６０】
πは上記オーダーリング行列であり、Πに帰属している。
【００６１】
式（４）は、Ｌｐノルムがλ以下であることを定めている（Ｌｐノルム正則化）。λはハイパーパラメータであり、予め与えておく。ｐはノルムの次数を表す。ｐ＝１はＬ１ノルム、ｐ＝２はＬ２ノルムを表す。
ｐ＝１の場合は｜｜Ｗ｜｜１＝｜ｗ１｜＋｜ｗ２｜＋・・・＋｜ｗＫ｜である。
ｐ＝２の場合は｜｜Ｗ｜｜２＝｜ｗ１｜２＋｜ｗ２｜２＋・・・＋｜ｗＫ｜２である。
【００６２】
ｐの値は予め決めておく。本実施形態ではパラメータ情報で指定した特徴波形の個数と同数の特徴波形を学習することを想定し、例えばｐ＝２とする（Ｌ２ノルム）。ｐ＝２とすることで、学習される特徴波形の個数が、指定した個数と同数になりやすい（重みｗ１，ｗ２，・・・，ｗｋのいずれもゼロになりにくい）。後述する第２の実施形態では、スパースモデリングにより、パラメータ情報で多数の個数を指定し、学習される特徴波形の個数を、指定した個数よりも小さい値に絞り込むため、ｐ＝１とする（Ｌ１ノルム）。この場合、ｗ１，ｗ２，・・・，ｗｋのうち多くのｗがゼロになり、結果として学習される特徴波形の個数が絞り込まれる。但し、これらに限定されず、ｐの値を任意に定めることができる。ｐを３以上としてもよい。
【００６３】
式（７）は、目的関数Ｈ（Ｘ，ｚ，π，Ｗ）を定義している。変数として、特徴ベクトルＸ、重みベクトルＷ、ｚ、πを含んでいる。
【００６４】
目的関数Ｈに含まれるΔβ（π＊，π）は、式（９）で定義されている。Δβ（π＊，π）は、行列πに含まれるすべての要素を合計し、行列πの要素数（Ｉ＋×βＩ−）で除算することより、行列πに含まれる１の個数の割合（誤り率）を計算している。πの添字ｉ＋は行、ｉ−は列を表す。π＊は、上記オーダーリング行列πと同じサイズの行列であり、すべての要素をゼロ（０）とした行列である。
【００６５】
目的関数Ｈに含まれるφｚ（Ｘ，π）は式（８）で定義されている。φｚ（Ｘ，π＊）は、式（８）のπをπ＊に置き換えればよい。φｚ（Ｘ，π＊）−φｚ（Ｘ，π）は、πにおいて０になっている要素に対応する成分はゼロになり、１になっている要素に対応する成分が基本的に非ゼロになる（もしπ＝π＊であれば、φｚ（Ｘ，π＊）−φｚ（Ｘ，π）はゼロベクトルである）。したがって、πにおいて０である要素に対応する正例と負例の特徴ベクトル
2020170214.tif 000007
が大きくなる（つまり、それらの特徴ベクトルに対して正しく重み付けされる）ようになる。
【００６６】
式（１０）は、前述した通り、時系列データｉと特徴波形ｋとの距離を定義している。式（１１）は、前述した通り、時系列データｉの特徴ベクトルを定義している。
【００６７】
上述の最適化問題は、確率的勾配法を用いて効率的に計算することが可能である。これにより、重みベクトルＷと特徴波形集合Ｓを効率的に算出できる。より詳細には、重み更新部１３は、確率的勾配法に基づき、重みベクトルＷを更新する。次に、特徴波形更新部１４は、確率的勾配法に基づき、特徴波形集合Ｓを更新する。特徴ベクトル生成部１１、重み更新部１３及び特徴波形更新部１４の処理を繰り返すことで、重みベクトルＷ及び特徴波形集合Ｓを学習する。ここでは確率的勾配法を用いたが、最急降下法など、他の種類の勾配法を用いてもよい。以下、重み更新部１３及び特徴波形更新部１４について詳細に説明する。
【００６８】
重み更新部１３が、モデルパラメータ（重みベクトル）Ｗを射影勾配降下法（確率的勾配法の一例である）により更新する。具体的には、目的関数Ｈ（Ｘ，ｚ，π，Ｗ）をモデルパラメータ（重みベクトル）Ｗで偏微分する。この際、式（３）のmaxに該当するｚとπを見つける。ｚは前述したようにスコアの大きい順に選択したｉβ−個の異常特徴ベクトルである。πは後述するように、ｚと、正常特徴ベクトルとに基づき特定すればよい。この特定したｚとπを目的関数Ｈ（Ｘ，ｚ，π，Ｗ）に代入してＷで直接微分することで、勾配∂Ｈ／∂Ｗの値を計算する。重み更新部１３は、計算した値（偏微分値）に基づき、重みベクトルＷを更新する。
【００６９】
例えば重みベクトルＷから、偏微分値（Ｗと同じ次元のベクトル）を減算する。つまり、Ｗの値を、偏微分値とは逆方向に動かす。これによりＨの値が小さくなる。ここでは偏微分値を減算したが、偏微分値に一定の係数をかけた値を減算してもよい。減算されたＷが前述の正則化の制約（式（４））を満たしているかを判断する。正則化の制約を満たしている場合は、減算されたＷを、更新されたＷとする。一方、減算されたＷが正則化の制約を満たしていない場合は、減算されたＷの値を、Ｌ２ボールに原点方向（すなわちＬ２ボールとのユークリッド距離が最小になる方向）に射影し、Ｌ２ボールに射影された点の値を、更新されたＷとする。Ｌ２ボールは、Ｌ２ノルムの距離（ユークリッド距離）を表すグラフである。ここでは半径λの球である。Ｌ２ボールは、Ｗの値の取り得る範囲を示している。一例としてλは１である。但し、λは１より大きくても、１より小さくてもよい。
【００７０】
図６は、λが１のときＬ２ボールに重みベクトルの値を射影する例を示す。Ｌ２ボールは半径１の円である。減算された重みベクトルＷａがＬ２ボールの外側に位置している。つまり、減算された重みベクトルＷａは正則化の制約を満たしていない。重みベクトルＷａの位置と、原点とを結ぶ直線とを算出する。算出した直線と、Ｌ２ボールとの交点Ｃａを算出する。交点Ｃａの値を、更新されたＷとする。減算された重みベクトルがＷｂの場合についても同様にして、重みベクトルＷｂの位置と原点とを結ぶ直線とを算出する。算出した直線と、Ｌ２ボールとの交点Ｃｂを算出する。交点Ｃｂの値を、更新されたＷとする。
【００７１】
上記の最適化問題では特徴ベクトルＸまたは特徴波形集合Ｓに依存する部分は、モデルパラメータ（重みベクトル）Ｗに依存しない。例えば式（１０）は、Ｗに依存しない。
【００７２】
特徴波形更新部１４が、特徴波形集合Ｓを確率的勾配降下法（確率的勾配法の一例である）により更新する。
【００７３】
図７は、特徴波形更新部１４の動作の一例を示すフローチャートである。まず正ラベル及び負ラベルの一方をランダムに選択する（Ａ０１）。ランダムに選択することで、正ラベルと負ラベルが均等に（同等の確率で）選択される。ランダムでなく、正ラベルと負ラベルとを交互に選択するなど、他の方法を用いてもよい。
【００７４】
選択したラベルが正ラベルか負ラベルかを判断する（Ａ０２）。
【００７５】
正ラベルが選択された場合、正ラベルの時系列データの特徴ベクトル（正常特徴ベクトル）のセットから、ランダムに１つの正常特徴ベクトルを選択する（Ａ０３）。正常特徴ベクトルは、一例として第１特徴ベクトルに対応する。
【００７６】
選択した正常特徴ベクトルに対応する時系列データのスコアを算出する。正常特徴ベクトルに対応する時系列データとは、正常特徴ベクトルが生成される元となった時系列データのことである。また、負ラベルの時系列データの特徴ベクトル（異常特徴ベクトル）のセットにおいて、各異常特徴ベクトルに対応する時系列データのスコアを算出する（Ａ０４）。異常特徴ベクトルに対応する時系列データとは、異常特徴ベクトルが生成される元となった時系列データのことである。異常特徴ベクトルは、一例として第２特徴ベクトルに対応する。
【００７７】
以下、正常特徴ベクトルに対応する時系列データのスコアを単に正常特徴ベクトルのスコア、異常特徴ベクトルに対応する時系列データのスコアを単に異常特徴ベクトルのスコアと称する。
【００７８】
スコアは、重み更新部１３で直近に更新されたモデルパラメータ（重みベクトル）と、特徴ベクトルとの内積により算出される。一例として、重みベクトルＷが（ｗ１，ｗ２）であり、時系列データｉの特徴ベクトルＸｉが（Ｘｉ，１，Ｘｉ，２）であれば、スコアはＷＴ・Ｘｉ＝ｗ１Ｘｉ，１，＋ｗ２Ｘｉ，２である。なおＴは転置を意味する。
【００７９】
スコアは、時系列データが正例（正常の時系列データ）及び負例（異常の時系列データ）をいずれに分類するかを判定するためのクラス分類モデルの出力値（予測値）に対応する。ここでは、スコアが大きいほど、その時系列データは正常である可能性が大きいことを意味する。但し、スコアが小さいほど、その時系列データは正常である可能性が大きいことを意味するようにしてもよい。
【００８０】
ステップＡ０３で選択した正常特徴ベクトルを対象特徴ベクトルとする（Ａ０４）。また、式（３）のmaxに該当するｚとπを見つける（同Ａ０４）。これに該当するｚは、異常特徴ベクトルＸｉをＷＴ・Ｘｉが大きい順にｉβ−個取り出した異常特徴ベクトルらである。これに該当するπを見つけるために、異常特徴ベクトルのスコアが正常特徴ベクトルのスコアよりも高くなってしまっている、異常特徴ベクトルと正常特徴ベクトルの組み合わせを見つける。この組み合わせは、性能指標が悪化する、異常特徴ベクトルと正常特徴ベクトルの組み合わせである。換言すれば、正常特徴ベクトルのスコアよりも正常クラスに属する可能性が高いスコアの異常特徴ベクトルを見つける。この組み合わせから、該当するπを見つける。つまり各異常特徴ベクトルをπの各列へ割り当てる割り当てパターンは複数存在するが、この中からmaxを満たすパターンを見つけ、見つけたパターンによりπを特定する。文献（A Structural SVM Based Approach for Optimizing Partial AUC, JMLR2013）のAlgorithm 2（Find Most-Violated Constraint）に本手法の詳細が記載されている（同Ａ０４）。
【００８１】
このように、本ステップでは、勾配降下法に必要な情報である、式（３）のmaxに該当するｚとπとの組を１つ見つける。
【００８２】
選択した対象特徴ベクトルと、maxに該当するｚとπとを用いて、確率的勾配降下法に基づき、特徴波形集合Ｓを更新する（Ａ０５）。
【００８３】
具体的には、目的関数Ｈ（Ｘ，ｚ，π，Ｗ）を特徴波形集合Ｓで偏微分した勾配∂Ｈ／∂Ｓを計算する。計算の一例を示す。微分公式の連鎖律（ｃｈａｉｎ  ｒｕｌｅ）を用いて、∂Ｈ／∂Ｓを以下のように変形できる。Ｘｉは対象特徴ベクトルを表す変数である。
∂Ｈ／∂Ｓ＝∂Ｈ／∂Ｘｉ・∂Ｘｉ／∂Ｓ      ・・・（１２）
【００８４】
ｚとπとに基づき∂Ｈ／∂Ｘｉの式を導出し、当該式における変数Ｘｉに上記対象特徴ベクトルを入力し、当該式における変数Ｗに、更新された重みベクトルＷの値を入力することで、∂Ｈ／∂Ｘｉの値を計算する。また、∂Ｘｉ／∂Ｓの式における変数Ｓに、現在の特徴波形（初期の特徴波形もしくは前回更新された特徴波形）を入力することで、∂Ｘｉ／∂Ｓの値を計算する。なお、∂Ｘｉ／∂Ｓの式は式（１０）から導けばよい。∂Ｈ／∂Ｘｉの値と、∂Ｘｉ／∂Ｓの値とを乗算することで、∂Ｈ／∂Ｓの値（偏微分値）を計算する。
【００８５】
∂Ｈ／∂Ｓの値（偏微分値）に基づき、特徴波形集合Ｓを更新する。例えば特徴波形集合Ｓから∂Ｈ／∂Ｓの値を減算する。つまり、Ｓの値を、偏微分値とは逆方向に動かす（Ｈの値を小さくする方向に動かす）。ここでは偏微分値を減算したが、偏微分値に係数をかけた値を減算してもよい。これにより、更新された特徴波形集合Ｓを得る。これは、前述したスコアの大小関係を訂正するように、特徴波形集合Ｓを更新することに相当する。
【００８６】
ステップＡ０２で負ラベルが選択されたと判断された場合、負ラベルの時系列データの特徴ベクトル（異常特徴ベクトル）から、ｐＡＵＣに影響する１つの異常特徴ベクトルを選択する（Ａ０６）。具体的には、まず、各異常特徴ベクトルのスコアを算出する。スコアは、重み更新部１３で直近に更新されたモデルパラメータ（重みベクトル）と、異常特徴ベクトルとの内積により算出する。これらの異常特徴ベクトルをスコアの大きい順にソートする（但し小さい順でもよい）。スコアが大きい順に、Ｕ（Ｕは１以上の整数）個の異常特徴ベクトルを特定する。Ｕ個は、性能指標パラメータに応じた決まる値であり、一例として、β×Ｉ個である。Ｕ個の異常特徴ベクトルの中から、１つの異常特徴ベクトルを選択する。例えば、ランダムに選択してもよい、最もスコアが大きい異常特徴ベクトルを選択してもよいし、その他の方法で選択してもよい。ここでは１つの異常特徴ベクトルを選択したが、２以上の少数（所定数）の異常特徴ベクトルを選択してもよい。ここではスコアの大きい順にＵ個の異常特徴ベクトルを選択したが、閾値以上のスコアを有する異常特徴ベクトルを選択してもよい。閾値は予め決めても良いし、任意に決めてもよい。
【００８７】
スコアが大きい異常特徴ベクトルは、偽陽性になる可能性が高い特徴ベクトルであるといえる。換言すれば、パラメータ情報として指定した偽陽性率の範囲に影響する特徴ベクトル、すなわちｐＡＵＣに影響する特徴ベクトルである。このような特徴ベクトルを選択することで、偽陽性率の指定した範囲でｐＡＵＣを最大化するという偽陽性率の条件に応じた特徴波形を効率的に学習することができる。
【００８８】
偽陽性率の範囲（β）は一例として０．０５又は０．０１以下など、通常、小さい範囲に指定される。したがって、単純に全ての異常特徴ベクトルの中から選択してしまうと、偽陽性率の範囲に影響する異常特徴ベクトルは殆ど選ばれず、その結果、偽陽性率の範囲のｐＡＵＣを最適化（最大化）する条件に適した特徴波形が学習されにくい。このため、本ステップＡ０６では、偽陽性率の範囲に影響する時系列データの特徴ベクトル（異常特徴ベクトル）をスコアに基づき特定し、その中から１つ（又は少数）の異常特徴ベクトルを選択している。この選択された特徴ベクトルを対象特徴ベクトルとする。
【００８９】
ステップＡ０６で選択した異常特徴ベクトルを対象特徴ベクトルとする（Ａ０７）。また、式（３）のmaxに該当するｚとπを見つける（同Ａ０７）。これに該当するｚは、異常特徴ベクトルＸｉをＷＴ・Ｘｉが大きい順にｉβ−個取り出した異常特徴ベクトルらである。これに該当するπを見つけるために、異常特徴ベクトルのスコアが正常特徴ベクトルのスコアよりも高くなってしまっている異常特徴ベクトルと正常特徴ベクトルの組み合わせを見つける。この組み合わせは、性能指標が悪化する、異常特徴ベクトルと正常特徴ベクトルの組み合わせである。換言すれば、対象特徴ベクトルのスコアよりも正常クラスに属する可能性が低いスコアの正常特徴ベクトルを見つける。この組み合わせから、上述したステップＡ０４と同様の方法で、該当するπを見つける（文献（A Structural SVM Based Approach for Optimizing Partial AUC, JMLR2013）のAlgorithm 2（Find Most-Violated Constraint）参照）。
【００９０】
このように、本ステップでは、勾配降下法に必要な情報である、式（３）のmaxに該当するｚとπとの組を１つ見つける。
【００９１】
選択した対象特徴ベクトルと、maxに該当するｚとπとを用いて、確率的勾配降下法に基づき、特徴波形集合Ｓを更新する（Ａ０５）。本ステップＡ０５の詳細は前述したため、省略する。
【００９２】
更新終了判定部１５は、モデルパラメータ（重みベクトル）Ｗ及び特徴波形集合Ｓの更新を終了するか判定する。具体的には、更新終了条件が満たされたか判定する。更新終了条件は、例えば特徴ベクトル生成部１１と、重み更新部１３と、特徴波形更新部１４との一連の処理の繰り返し回数により規定される。例えば繰り返し回数が予め定めた回数（１００００回など）に達したかを判定する。予め定めた回数に達した場合は、重みベクトルＷと特徴波形集合Ｓの学習が十分に行われたと判断し、処理を終了する。予め定めた回数に達していない場合は、予め定めた回数に達するまで、上記一連の処理を繰り返す。更新終了条件を更新回数により規定することにより、学習に要する時間を所望の範囲内に設定することができる。
【００９３】
更新終了判定部１５は、学習の結果から、予測用（判定用）の閾値を決定してもよい。例えば偽陽性率が一定値以下になるよう閾値を決定してもよい。あるいは、閾値は予め定め与えてもよい。
【００９４】
出力情報記憶部２０は、学習により得られた特徴波形集合Ｓのデータと重みベクトルＷのデータとを内部に格納する。
【００９５】
表示部１７は、学習により得られた特徴波形集合Ｓ及び重みベクトルＷを含む学習結果データを画面に表示する。
【００９６】
図８に表示部１７に表示された学習結果データの例を示す。図８（Ａ）は、学習用の正ラベルの時系列データ（サンプル３）と、学習された２つの特徴波形（それぞれＳ１、Ｓ２とする）が表示されている。特徴波形Ｓ１、Ｓ２はいずれも異常（負ラベル）を検出するのに有効な特徴波形である。ここでは図５に示した条件で学習を行った例を示す（偽陽性率の範囲の上限値（β）の値を０．０１としている）。横軸は時間、縦軸は振幅である。サンプル３は、学習用の全時系列データのうち３番目の時系列データである。この時系列データ（サンプル３）と、特徴波形Ｓ１、Ｓ２と、学習された重みベクトルＷから計算されるスコアは８である。８は閾値以上であり、サンプル３は正しく分類されている。特徴波形Ｓ１、Ｓ２のグラフは、サンプル３に対してそれぞれ最も距離が近い位置に配置されている。特徴波形Ｓ１、Ｓ２は、サンプル３に対して最も距離が近い位置でも、サンプル３との距離が離れている。すなわち、特徴波形Ｓ１、Ｓ２は、サンプル３とフィットしていない。正常を示す情報を表示してもよい。当該情報及びスコアは、スコアに基づく情報の一例である。ここではサンプル３を表示したが、他の時系列データも同様にして表示されてよい。
【００９７】
図８（Ｂ）は、学習用の負ラベルの時系列データ（サンプル１１、４）と、学習された特徴波形Ｓ１及び特徴波形Ｓ２が表示されている。特徴波形Ｓ１及び特徴波形Ｓ２は図８（Ａ）と同じである。この時系列データ（サンプル１１）と、特徴波形Ｓ１、Ｓ２と、学習された重みベクトルＷから計算されるスコアは−５である。−５は閾値未満であり、したがって、サンプル１１は正しく分類されている。同様に、時系列データ（サンプル４）と、特徴波形Ｓ１、ｓ２と、学習された重みベクトルＷから計算されるスコアも−５である。−５は閾値未満であり、したがって、サンプル４も正しく分類されている。サンプル１１、４に対して、異常を示す情報をとして表示してもよい。当該情報及びスコアは、スコアに基づく情報の一例である。
【００９８】
サンプル１１に対して特徴波形Ｓ１、Ｓ２のグラフがそれぞれ最も距離が近い位置に配置されている。特徴波形Ｓ２はサンプル１１にフィットしていない（サンプル１１との距離が大きい）ものの、特徴波形Ｓ１が、サンプル１１の部分波形にフィットしている（サンプル１１との距離が近い）。特徴波形Ｓ１はこのフィットしている部分波形を、異常に特有な形状として検出するのに有効である。
【００９９】
サンプル４に対して特徴波形Ｓ１、Ｓ２のグラフがそれぞれ最も距離が近い位置に配置されている。特徴波形Ｓ１はサンプル４にフィットしていないが、特徴波形Ｓ２が、サンプル４の部分波形にフィットしている。特徴波形Ｓ２はこのフィットしている部分波形を異常に特有な形状として検出するのに有効である。学習用の負ラベルの時系列データのセットにおいて、このような形状を持つ時系列データの個数がたとえ少なくても（例えば負ラベルの時系列データの多くは、サンプル１１において特徴波形Ｓ１にフィットしている部分波形に類似の部分波形を有する場合でも）、当該時系列データを負ラベルに分類すべきものとして有効に検出できる。
【０１００】
図９は、学習された特徴波形Ｓ１、Ｓ２とモデルパラメータ（重みベクトル）Ｗとを用いて、学習に用いた時系列データに対して正ラベル又は負ラベルの予測（判定）を行い、その結果に基づき作成したＲＯＣ曲線の例を示す。判定は、時系列データに対してスコアを計算し、スコアを閾値と比較することで行う。スコアが閾値以上であれば正ラベル、閾値未満であれば負ラベルに分類する。偽陽性率が０以上０．０１以下の範囲でＲＯＣ曲線と囲まれた部分の面積がｐＡＵＣである。ｐＡＵＣが大きくなっている。このため、偽陽性率が０．０１以下の範囲でも、高精度に時系列データに対する正常・異常の判定を行うことができる。
【０１０１】
図８では偽陽性率の範囲の上限値（β）を０．０１としたが、以下、βを１．０として学習を行った場合の例を図１０及び図１１を用いて説明する。
【０１０２】
図１０は、この場合のＧＵＩの画面の例を示す。ユーザは、ＧＵＩを解して、性能指標パラメータの範囲を調整つまみ３１の操作（ユーザ入力情報）により変更し、最大値の１．０に設定する。それ以外の条件は図５と同じである。変更後の偽陽性率の範囲で学習を行う（つまりＡＵＣを最大化するように学習する）。このように性能指標パラメータを変更することで、学習される特徴波形、重みベクトルＷが変わり、また、計算されるスコア、判定されるラベル、偽陽性率の範囲（β）における予測性能も変わる。変更後の偽陽性率の範囲で学習を行った結果得られる学習結果データの例を、図１１を用いて説明する。
【０１０３】
図１１（Ａ）は、学習用の正ラベルの時系列データ（サンプル２）と、学習された特徴波形（Ｓ３とする）が表示されている。特徴波形は１つのみ学習されている（これは、重み更新処理部の処理で他の特徴波形の重みが０になったため、当該他の特徴波形が出力されなかった場合である）。この時系列データ（サンプル２）と、特徴波形Ｓ３と、学習された重みベクトルＷから計算されるスコアは１０である。１０は閾値以上であり、サンプル２は正しく分類されている。特徴波形Ｓ３のグラフは、サンプル２に対してそれぞれ最も距離が近い位置に配置されている。特徴波形Ｓ３は、サンプル２に対して最も距離が近い位置でも、サンプル２との距離が離れている。すなわち、特徴波形Ｓ３は、サンプル２とフィットしていない。
【０１０４】
図１１（Ｂ）は、学習用の負ラベルの時系列データ（サンプル８）と、学習された特徴波形Ｓ３が表示されている。特徴波形Ｓ３は図１１（Ａ）と同じである。この時系列データ（サンプル８）と、特徴波形Ｓ３と、学習された重みベクトルＷから計算されるスコアは−５である。−５は閾値未満であり、サンプル８は正しく分類されている。
【０１０５】
図１２は、学習された特徴波形Ｓ３と重みベクトルＷとに基づき、学習に用いた時系列データに対して正ラベル又は負ラベルの判定を行い、その結果に基づき作成したＲＯＣ曲線の例を示す。偽陽性率が０以上０．０１以下の範囲でＲＯＣ曲線と囲まれた部分の面積（ｐＡＵＣ）は、図９よりも小さいことが理解される。
【０１０６】
つまり、β＝１．０として、（ｐＡＵＣではなく）ＡＵＣを最適化（最大化）するように学習された結果、図８で示した特徴波形Ｓ２のような、発生回数が少ない形状を異常として検出するのに有効な特徴波形を学習できない。このため、このような形状をもつ時系列データを異常として検出できない。このことから、図９のようにβの小さい範囲でｐＡＵＣを最適化するよう学習を行うことで、発生回数の少ない形状の異常をも検出可能な特徴波形及び重みベクトルＷを学習できる。また、ＲＯＣ曲線は単調増加関数であるため、ｐＡＵＣを最適化した結果、ＡＵＣも、図１２に比べて大きくなる。よって、図９のようにβを小さい範囲に指定して学習を行ったクラス分類モデル（重みベクトルＷ）の方が、全体としての予測性能も高い。但し、β＝１．０としても、異常検出に有効な特徴波形Ｓ３を学習できるため、本実施形態はこの場合も、異常検出のための学習に有効である。
【０１０７】
図１３は、学習されたモデルパラメータ（重みベクトル）により表される識別境界を模式的に示す図である。図１３（Ａ）は、線形の識別境界の例、図１３（Ｂ）が、非線形の識別境界の例を示す。いずれもこの例では、特徴空間は、２次元である。図１３（Ａ）に示すように、線形の識別境界の場合、識別境界は直線によって表される。この例では、直線に対して、上側が正常領域（上側）、下側が異常領域である。黒丸は特徴ベクトルを表す。図１３（Ｂ）に示すように、非線形の識別境界の場合、識別境界は、複雑な形状になっている。この例では、識別境界の内側が正常領域、外側が異常領域である。
【０１０８】
図１４は、学習フェーズの動作のフローチャートである。
【０１０９】
ステップＡ１１において、入力設定部１０は、複数の特徴波形（特徴波形集合）と、クラス分類モデルの重みベクトルＷとを初期化する。クラス分類モデルは、クラス分類器又は予測器などとも称される。
【０１１０】
ステップＡ１２において、特徴ベクトル生成部１１は、学習用の各時系列データに対して複数の特徴波形との距離を計算し、計算した複数の距離を要素として含む特徴ベクトルを生成する。
【０１１１】
ステップＡ１３において、重み更新部１３は、性能指標パラメータに基づく目的関数Ｈ（式（３）参照）と、各時系列データの特徴ベクトルとに基づき、確率的勾配法に基づき、重みベクトルＷを更新する。具体的には、目的関数をＷで偏微分（δＨ／δＷ）した式における変数Ｘに各特徴ベクトルを入力して、偏微分値（Ｗの勾配）を計算する。現在の重みベクトルＷから、偏微分値又はこれに係数を乗じた値を減算する。減算された重みベクトルが、Ｌ２ノルムの正則化の制約を満たすかを判断する（式（４）参照）。制約を満たしていない場合には、Ｌ２ボールへ重みベクトルを射影することにより（図６参照）、重みベクトルを更新する。Ｌ２ノルムの正則化の制約を満たす場合は、減算された重みベクトルを、更新された重みベクトルとする。
【０１１２】
ステップＡ１４において、特徴波形更新部１４は、性能指標パラメータに基づき例えば１つの時系列データｉを対象特徴ベクトルＸｉとして選択する。選択した対象特徴ベクトルＸｉに基づき、各特徴波形を更新する。例えば目的関数Ｈを特徴波形集合Ｓで偏微分し（δＨ／δＳ）、これをδＨ／δＸｉ・δＸｉ／δＳに変形する。δＨ／δＸｉ及びδＸｉ／δＳの値をそれぞれ計算し、計算した値を乗算することで、δＨ／δＳ（Ｓの勾配）を計算する。特徴波形集合Ｓから、Ｓの勾配の値を減算することにより、特徴波形集合Ｓを更新する。
【０１１３】
ステップＡ１５において、更新終了判定部１５が、更新終了条件が満たされたか判定する。更新終了条件は、例えばステップＡ１２〜Ａ１４の繰り返し回数（更新回数）が閾値に達することである。更新終了条件が満たされない間は（ＮＯ）、ステップＡ１２〜Ａ１４を繰り返す。更新終了条件が満たされた場合は（ＹＥＳ）、ステップＡ１６に進む。
【０１１４】
ステップＡ１６において、更新された複数の特徴波形のデータと、更新された重みベクトルのデータとを出力し、出力情報記憶部２０に格納する。なお、重みが０の特徴波形が存在する場合は、その特徴波形を出力しなくてよい。また、当該重み０の要素を重みベクトルから除去する（削除した要素の個数だけ、重みベクトルが短くなる）。
【０１１５】
＜テストフェーズ＞
テストフェーズでは、学習した特徴波形の集合及び重みベクトルを入力として与えて、テスト用の時系列データのスコアを算出し、算出したスコアに基づき、時系列データに対して正ラベル又は負ラベルを決定（正常か異常かを決定）する。以下、詳細に説明する。
【０１１６】
テストフェーズでは、入力設定部１０と、テスト用データ記憶部１９と、出力情報記憶部２０と、特徴ベクトル生成部１１と、予測部１６と、表示部１７とを用いる。
【０１１７】
出力情報記憶部２０には、学習フェーズで最終的に得られた更新後の特徴波形集合Ｓ（Ｋ個の更新後の特徴波形を含む）と、重みベクトル（モデルパラメータ）Ｗとが記憶されている。
【０１１８】
テスト用データ記憶部１９は、テスト対象となる時系列データを記憶している。この時系列データは、テスト対象となる分析対象装置に設置されたセンサの検出値に基づくものである。
【０１１９】
入力設定部１０は、テスト対象となる時系列データをテスト用データ記憶部１９から読み出して、特徴ベクトル生成部１１に入力する。
【０１２０】
特徴ベクトル生成部１１は、テスト対象の時系列データ（ｔとする）を読み出し、当該時系列データと、特徴波形集合Ｓとの距離に基づき、Ｋ次元の特徴ベクトル（Ｘｔとする）を計算する。
【０１２１】
予測部１６は、特徴ベクトルＸｔと、重みベクトルＷとに基づきスコアを算出する。具体的には、クラス分類モデルの式（スコアの算出式）を
Ｙ＝Ｘｔ・Ｗ      ・・・（１３）
とする。このとき、特徴ベクトルＸｔと、重みベクトルＷとの内積を計算することで、スコアを算出する。例えば、Ｋ＝２であり、特徴ベクトルＸｔが（Ｘｔ，１，Ｘｔ，２）、重みベクトルＷが（ｗ１，ｗ２）であれば、Ｘｔ，１・ｗ１＋Ｘｔ，２・ｗ２によりスコアを算出する。
【０１２２】
予測部１６は、算出したスコアを閾値と比較する。閾値は、前述したように、更新終了判定部１５が決定した値でもよいし、予め定めた値でもよい。スコアが閾値以上であれば、テスト対象となる時系列データは正常であると判定し、正ラベルを付与する。閾値未満であれば、テスト対象となる時系列データは異常であると判定し、負ラベルを付与する。
【０１２３】
表示部１７は、予測部１６の評価結果に基づくデータ（評価結果データ）を画面に表示する。評価結果データは、テスト対象となる時系列データ、付与されたラベルの情報（正常又は異常の判定結果）、当該時系列データに最も近い位置に配置された特徴波形、算出されたスコアを含む。評価結果データは、上述した学習結果データと同様の形式を有する。判定結果が異常の場合のみ、評価結果データを表示してもよい。出力情報記憶部２０は評価結果データを内部に記憶してもよい。
【０１２４】
図１５は、表示部１７に表示された評価結果データの一例を示す。テスト対象となる時系列データ（テストデータＮｏ．１１８）が表示されている。特徴波形Ｓ１、Ｓ２がそれぞれ時系列データに対して最も近い位置に配置されている。特徴波形Ｓ２が時系列データの部分波形にフィットしている。スコアは−４．９であり、負ラベルが付与されている。すなわち、時系列データは異常である（分析対象装置に異常が発生している）と判定されている。異常と判定されたことを示す情報が表示されている。当該情報及びスコアは、スコアに基づく情報の一例である。ユーザは、評価結果データを確認することで、異常の発生を判断できる。また、ユーザは、どの特徴波形が時系列データにフィットしているかに応じて、異常と判断された根拠を特定できる。例えば特徴波形Ｓ２がフィットしていれば、発生頻度の低い異常が発生したと判断できる。
【０１２５】
なお、重みベクトル（モデルパラメータ）Ｗとして、サポートベクトルの集合Ｓｖと、寄与率の集合Ｓａを用いる場合、クラス分類モデルを、以下のように生成する。（Ｓａ，Ｓｖ）は識別境界のモデルパラメータに対応し、Ｘは入力変数（特徴ベクトル）に対応する。Ｙは出力である。Ｙに−１を掛けた“−Ｙ”を異常度と定義する。Ｋはカーネル関数であり、ＳｖはサポートベクトルＳ’ｖの集合である。ＳａはＳｖに属するサポートベクトルの寄与率Ｓ’ａの集合である。
【数４】
2020170214.tif 000008
【０１２６】
予測部１６は、計算された異常度“−Ｙ”が閾値以上であれば、テスト対象となる時系列データは異常であると判定し、負ラベルを付与する。異常度“−Ｙ”が閾値未満であれば、テスト対象となる時系列データは正常であると判定し、正ラベルを付与する。
【０１２７】
図１６は、テストフェーズの動作のフローチャートである。
【０１２８】
ステップＡ２１において、特徴ベクトル生成部１１が、テスト対象の時系列データと、学習された各特徴波形との距離を計算し、計算した距離を各特徴波形に対応する要素に格納した特徴ベクトルを生成する。
【０１２９】
ステップＡ２２において、予測部１６が、特徴ベクトルと、学習された重みベクトルとの内積により、スコア（評価値）を算出する。
【０１３０】
ステップＡ２３において、予測部１６が、算出したスコアを閾値と比較する。閾値以上であれば、テスト対象となる時系列データは正常であると判定し、正ラベルを付与する。閾値未満であれば、テスト対象となる時系列データは異常（分析対象装置に異常が発生している）であると判定し、負ラベルを付与する。
【０１３１】
ステップＡ２４において、表示部１７は、評価結果データを画面に表示する。評価結果データは、一例として、テスト対象となる時系列データ、付与されたラベルの情報（正常又は異常の判定結果）、当該時系列データに最も近い位置に配置された特徴波形、及び算出されたスコアを含む。
【０１３２】
本実施形態では性能指標パラメータとして偽陽性率の範囲を直接、ＧＵＩから指定する例を示したが、ユーザに調整可能として提示する性能指標パラメータの名称は偽陽性率でなくてもよい。例えば、誤検知と見逃しを重視する度合いを調整するパラメータをユーザに調整可能に提示してもよい。または、上位ランキングの何位まで（スコアが最も大きい側又は最も小さい側から数えた時系列データの個数）の判定の正しさを重視するかを調整するパラメータをユーザに調整可能に提示してもよい。いずれの場合にも、本装置の学習において、内部的に、ユーザが指定したパラメータの値を、偽陽性率の範囲に置き換える。例えばユーザの指定したパラメータの値と、偽陽性率の範囲とを対応付けた情報（例えば関数又はテーブル）を予め学習用データ記憶部１８に格納しておくことで、このような置き換えを実現する。この後の処理は、前述した本実施形態と同様である。
【０１３３】
図１７に、本実施形態に係る時系列データ分析装置のハードウェア構成を示す。本実施形態に係る時系列データ分析装置は、コンピュータ装置１００により構成される。コンピュータ装置１００は、ＣＰＵ１０１と、入力インターフェース１０２と、表示装置１０３と、通信装置１０４と、主記憶装置１０５と、外部記憶装置１０６とを備え、これらはバス１０７により相互に接続されている。
【０１３４】
ＣＰＵ（中央演算装置）１０１は、主記憶装置１０５上で、コンピュータプログラムである分析プログラムを実行する。分析プログラムは、時系列データ分析装置の上述の各機能構成を実現するプログラムのことである。ＣＰＵ１０１が、分析プログラムを実行することにより、各機能構成は実現される。
【０１３５】
入力インターフェース１０２は、キーボード、マウス、及びタッチパネルなどの入力装置からの操作信号を、時系列データ分析装置に入力するための回路である。
【０１３６】
表示装置１０３は、時系列データ分析装置から出力されるデータまたは情報を表示する。表示装置１０３は、例えば、ＬＣＤ（液晶ディスプレイ）、ＣＲＴ（ブラウン管）、及びＰＤＰ（プラズマディスプレイ）であるが、これに限られない。学習用データ記憶部１８、テスト用データ記憶部１９、出力情報記憶部２０に記憶されたデータまたは情報は、この表示装置１０３により表示することができる。
【０１３７】
通信装置１０４は、時系列データ分析装置が外部装置と無線又は有線で通信するための回路である。学習用の入力データまたはテスト用の時系列データなどのデータは、通信装置１０４を介して外部装置から入力することができる。外部装置から入力したデータを、学習用データ記憶部１８またはテスト用データ記憶部１９に格納することができる。
【０１３８】
主記憶装置１０５は、分析プログラム、分析プログラムの実行に必要なデータ、及び分析プログラムの実行により生成されたデータなどを記憶する。分析プログラムは、主記憶装置１０５上で展開され、実行される。主記憶装置１０５は、例えば、ＲＡＭ、ＤＲＡＭ、ＳＲＡＭであるが、これに限られない。学習用データ記憶部１８、テスト用データ記憶部１９及び出力情報記憶部２０は、主記憶装置１０５上に構築されてもよい。
【０１３９】
外部記憶装置１０６は、分析プログラム、分析プログラムの実行に必要なデータ、及び分析プログラムの実行により生成されたデータなどを記憶する。これらのプログラムやデータは、分析プログラムの実行の際に、主記憶装置１０５に読み出される。外部記憶装置１０６は、例えば、ハードディスク、光ディスク、フラッシュメモリ、及び磁気テープであるが、これに限られない。学習用データ記憶部１８、テスト用データ記憶部１９及び出力情報記憶部２０は、外部記憶装置１０６上に構築されてもよい。
【０１４０】
なお、分析プログラムは、コンピュータ装置１００に予めインストールされていてもよいし、ＣＤ−ＲＯＭなどの記憶媒体に記憶されていてもよい。また、分析プログラムは、インターネット上にアップロードされていてもよい。
【０１４１】
本実施形態では、時系列データ分析装置が、学習フェーズと、テストフェーズとの両方を行う構成を備えていたが、いずれか一方のみを行う構成でもよい。つまり、学習フェーズを行う装置と、テストフェーズを行う装置を別々に構成してもよい。
【０１４２】
本実施形態ではクラス分類モデルとしてＳＶＭを用いたが、ロジスティック回帰モデルなど、他のモデルを用いてもよい。
【０１４３】
以上、本実施形態によれば、ｐＡＵＣ等の性能指標を最適化するように時系列データに対するクラス分類モデルの重みを学習する。したがって、性能指標の条件（例えば偽陽性率の範囲）に合わせてクラス分類モデルを生成することで、その条件の元で高い予測性能を達成することができる。また、本実施形態では、当該条件の元で予測に有効な特徴波形を学習できる。
【０１４４】
また、本実施形態によれば、時系列データにおいて、ｐＡＵＣ等の性能指標を直接最適化しつつ、予測に有効な特徴波形を学習できる。時系列データのラベル（異常有無）の予測では、学習された重みに基づき高い精度での予測が可能になるとともに、予測の根拠（いずれの特徴波形が時系列データにフィットしたか）も提示できる。
【０１４５】
ｐＡＵＣを性能指標として最適化することで、例えば、（１）トラブルの見逃しを十分低く保ったもとで、正しくトラブルを予測する、（２）誤診断を低く抑えたもとで、正しく診断事例を予測する、（３）少数の上位のランキングの予測精度を正確するなどの要求に応えることができる。
【０１４６】
（第２の実施形態）
本実施形態では、最初に多数の特徴波形を学習用に指定し、これらの特徴波形の中から、予測に有効な特徴波形をスパースモデリングにより絞り込むことを特徴とする。
【０１４７】
図１８は、第２の実施形態に係る時系列データ解析装置のブロック図である。特徴波形絞り込み部２１が追加されている。それ以外のブロックは図１と同じであるため、拡張又は変更された機能を除き、説明を省略する。
【０１４８】
第１の実施形態では特徴波形のパラメータ情報として特徴波形の個数を例えば２など、所望数もしくはそれに近い値を指定した。これに対して、本実施形態では、入力設定部１０は、第１の実施形態よりも十分大きな個数、例えば最大個数（例えば５００）を指定する。
【０１４９】
重み更新部１３は、最適化問題における正則化の制約（式（４）参照）としてＬ１ノルム（Ｌａｓｓｏ）正則化を用いる。つまりｐ＝１とする。第１の実施形態ではＬ２ノルム正則化（ｐ＝２）であったが、本実施形態では、Ｌ１ノルム正則化（ｐ＝１）を用いる。
【０１５０】
この変更により、重み更新部１３では、Ｌ１ノルム正則化による射影勾配降下法を用いる。第１の実施形態ではＬ２ノルム正則化の制約のため射影用にＬ２ボールを用いたが、本実施形態ではＬ１ノルム正則化の制約のため、射影用にＬ１ノルムの距離（マンハッタン距離）を表すＬ１ボールを用いる。Ｌ１ボールは、原点から各頂点までの距離がλである正方形であり、重みベクトルＷの取り得る範囲を示している。一例としてλは１である。但し、λは１より大きくても、１より小さくてもよい。δＨ／δＷの値だけ減算された重みベクトルＷがＬ１ノルム正則化の制約を満たさない場合は、上記減算された重みベクトルＷをＬ１ボールに対して射影する。射影は、ユークリッド距離が最小になる辺上の位置に対して行う。射影された位置の値を、更新された重みベクトルＷとする。
【０１５１】
図１９は、λが１のときのＬ１ボールと、Ｌ１ボールに重みベクトルＷを射影する例を示す。減算された重みベクトルＷｄがＬ１ボールの外側に位置している。重みベクトルＷｄの位置からユークリッド距離が最小となる辺上の位置に射影する。射影された点Ｃｄの値を、更新されたＷとする。減算された重みベクトルがＷｅ又はＷｆの場合についても同様にして、Ｌ１ボールに対してユークリッド距離が最小となる辺上の位置に射影する。Ｌ１ボールとの交点Ｃｆ又はＣｇを算出し、算出した交点Ｃｂ又はＣｇの値を、更新されたＷとする。
【０１５２】
Ｌ１ボールでは、Ｌ２ボール（図６参照）と比べて各軸との交点が尖っているため、その交点に射影される場合が多くなり、その結果、重みベクトルのいくつかの成分が０となるもしくは０になり易い。
【０１５３】
特徴波形更新部１４の動作は、第１の実施形態と同じである。
【０１５４】
特徴波形絞り込み部２１は、重み更新部１３で更新された重みベクトルに基づき、重みが０の要素が存在するかを判断し、重みが０の要素が存在する場合は、当該要素に対応する特徴波形を削除する。すなわち、現在存在する複数の特徴波形のうち、重みが０でない特徴波形に絞り込む。削除した要素の個数分だけ、重みベクトルの長さは短くなる。
【０１５５】
更新終了判定部１５は、更新終了条件が満たされたかを判断し、更新終了条件が満たされるまで、特徴ベクトル生成部１１、重み更新部１３、特徴波形更新部１４、及び特徴波形絞り込み部２１の動作を繰り返す。更新終了条件は第１の実施形態と同様である。
【０１５６】
Ｌ１ノルム正則化では、最終的に、多くの成分の重みが０になる。更新終了判定部１５は、最終的に重みの値が０となった特徴波形を学習結果として出力しない。重みが０でない特徴波形を、学習された特徴波形として出力し、出力情報記憶部２０に格納する。これにより、第１の実施形態よりも、予測精度の高いクラス分類モデル（重みベクトルＷ）と、検出精度の高い必要最小限の個数の特徴波形とを学習することができる。重みが０でない特徴波形を出力したが、重みが所定値以上の特徴波形のみを出力し、所定値未満の特徴波形を出力しないことも可能である。
【０１５７】
図２０は、第２の実施形態に係る学習フェーズの動作のフローチャートである。
【０１５８】
ステップＡ１１、Ａ１２は、第１の実施形態の図１４のステップＡ１１、Ａ１２と同じである。
【０１５９】
ステップＡ１７は、Ｌ１ノルム正則化の制約を満たすための処理を除き、第１の実施形態の図１４のステップＡ１３と同じである。第１の実施形態のステップＡ１３ではＬ２ノルム正則化を用いたが、本実施形態のステップＡ１３ではＬ１ノルム正則化を用いる。すなわち、δＨ／δＷの値だけ減算された重みベクトルが、Ｌ１ノルムの正則化の制約を満たすかを判断する（式（４）参照）。制約を満たしていない場合には、Ｌ１ボールへ、減算された重みベクトルを射影することにより（図１９参照）、重みベクトルを更新する。Ｌ１ノルム正則化の制約を満たす場合は、減算された重みベクトルを、更新された重みベクトルとする。
【０１６０】
ステップＡ１４は、第１の実施形態の図１４のステップＡ１４と同じである。
【０１６１】
ステップＡ１８において、特徴波形絞り込み部２１は、更新された重みベクトルに基づき、重みが０の要素が存在するかを判断し、重みが０の要素が存在する場合は、当該要素に対応する特徴波形を削除する。すなわち、現在存在する複数の特徴波形のうち、重みが０でない特徴波形に絞り込む。
【０１６２】
ステップＡ１５は、第１の実施形態の図１４のステップＡ１５と同じである。
【０１６３】
ステップＡ１９において、更新された１つ又は複数の特徴波形（最終的に残った１つ又は複数の特徴波形、すなわち重みが０でない特徴波形）のデータと、更新された重みベクトルのデータとを出力し、出力情報記憶部２０に格納する。
【０１６４】
以上、本実施形態によれば、スパースモデリングを用いることで多数の初期の特徴波形から絞り込みをかけつつ、予測に有効な特徴波形を学習できる。また、精度の高いクラス分類モデル（重み）を学習できる。
【０１６５】
（第３の実施形態）
第３の実施形態では、時系列データ分析装置が、通信ネットワークを介して、分析対象装置に接続された時系列データ分析システムの実施形態を示す。
【０１６６】
図２１に、第３の実施形態に係る時系列データ分析システムを示す。時系列データ分析装置４１は、第１〜第２の実施形態のいずれかに係る時系列データ分析装置に相当する。時系列データ分析装置４１は、通信ネットワーク４２を介して、複数の分析対象装置４３に接続されている。分析対象装置４３には、物理量を検出するセンサが搭載されている。分析対象装置４３は、センサの検出値に基づく時系列データを生成し、生成した時系列データを、通信ネットワーク４２を介して、時系列データ分析装置４１に送信する。時系列データ分析装置４１は、学習フェーズ用に時系列データを収集する場合、各分析対象装置４３が事前に正常状態及び異常状態のいずれにあるかを確認しておく。時系列データ分析装置４１は、正常状態にある分析対象装置４３から受信した時系列データに正ラベルを付し、異常状態にある分析対象装置４３から受信した時系列データに負ラベルを付して、学習用データ記憶部１８に格納する。また、時系列データ分析装置４１は、テストフェーズ用に時系列データを収集する場合は、受信した時系列データをテスト用データ記憶部１９に格納する。これにより、リアルタイムに分析対象装置４３の異常有無をテストできる。
【０１６７】
なお、本発明は上記各実施形態そのままに限定されるものではなく、実施段階ではその要旨を逸脱しない範囲で構成要素を変形して具体化できる。また、上記各実施形態に開示されている複数の構成要素を適宜組み合わせることによって種々の発明を形成できる。また例えば、各実施形態に示される全構成要素からいくつかの構成要素を削除した構成も考えられる。さらに、異なる実施形態に記載した構成要素を適宜組み合わせてもよい。
【符号の説明】
【０１６８】
１０：入力設定部
１１：特徴ベクトル生成部
１２：更新処理部
１３：重み更新部
１４：特徴波形更新部
１５：更新終了判定部
１６：予測部
１７：表示部
１８：学習用データ記憶部
１９：テスト用データ記憶部
２０：出力情報記憶部
２１：特徴波形絞り込み部
【図１】
2020170214.tif 000009
【図２】
2020170214.tif 000010
【図３】
2020170214.tif 000011
【図４】
2020170214.tif 000012
【図５】
2020170214.tif 000013
【図６】
2020170214.tif 000014
【図７】
2020170214.tif 000015
【図８】
2020170214.tif 000016
【図９】
2020170214.tif 000017
【図１０】
2020170214.tif 000018
【図１１】
2020170214.tif 000019
【図１２】
2020170214.tif 000020
【図１３】
2020170214.tif 000021
【図１４】
2020170214.tif 000022
【図１５】
2020170214.tif 000023
【図１６】
2020170214.tif 000024
【図１７】
2020170214.tif 000025
【図１８】
2020170214.tif 000026
【図１９】
2020170214.tif 000027
【図２０】
2020170214.tif 000028
【図２１】
2020170214.tif 000029