JP2017138929A - 時系列データ波形分析装置、方法、及びプログラム - Google Patents

Abstract

【課題】計算量が低減された時系列データ波形分析装置、方法、及びプログラムを提供する。【解決手段】一実施形態に係る時系列データ波形分析装置は、特徴波形追加部と、特徴量計算部と、パラメータ更新部と、特徴波形除去部と、を備える。特徴波形追加部は、ラベル付き時系列データに含まれる部分時系列の一部を、特徴波形集合に特徴波形として追加する。特徴量計算部は、１つ又は複数のラベル付き時系列データをランダムに抽出し、抽出されたラベル付き時系列データに対する、特徴波形の特徴量をＴＳＳ法に基づいて計算する。パラメータ更新部は、特徴量に基づいて、特徴波形と、特徴波形に対応する重み係数と、を含むパラメータを確率的勾配降下法に基づいて更新する。特徴波形除去部は、対応する重み係数が０である特徴波形を、特徴波形集合から除去する。評価関数生成部は、特徴波形と、重み係数と、に基づいて、評価関数を生成する。【選択図】図１

Description

本発明の実施形態は、時系列データ波形分析装置、方法、及びプログラムに関する。

近年、経済時系列分析や工場内製造プロセスのセンサ監視などの様々な分野において、時系列データの波形を分析する技術が重要となっている。この技術は、時系列データの分類、予測、及びスコア付けなど（以下、「分類等」という）に利用される。

時系列データの波形を分析する方法として、Time Series Shapelets法（以下、「ＴＳＳ法」という）が注目されている。ＴＳＳ法を利用することにより、高精度な分析が可能であることが報告されている。ＴＳＳ法では、学習用の複数の時系列データに基づいて、特徴波形（shapelet）の特徴量と、特徴量の重み係数と、を含む評価関数を生成し、この評価関数を利用して、時系列データの分類等を行う。

従来のＴＳＳ法では、特徴波形の全候補数をＱ、時系列データの数をＮ、各時系列データの長さ（時系列長）をＬとすると、評価関数を生成するための計算量は、Ｏ（Ｑ×Ｎ×Ｌ×Ｌ）となった。全候補数Ｑは、時系列データから生成可能な部分時系列の総数であるため膨大である。また、計算量は、時系列長の２乗で増大する。結果として、従来のＴＳＳ法では、計算量が大きく、実時間での計算が困難になる恐れがあった。

計算量を低減するために、予めＫ個（Ｋ≪Ｑ）の特徴波形候補を用意する方法も提案されている。この方法によれば、計算量は、Ｏ（Ｋ×Ｎ×Ｌ×Ｌ）となる。しかしながら、この方法では、用意された特徴波形候補が適当でない場合、評価関数を求めるための目的関数の解が収束せず、評価関数が得られないという問題もある。

"Learning Time-Series Shapelets", KDD '14 Proceedings of the 20th ACM SIGKDD international conference on Knowledge discovery and data miningPages 392-401／Josif Grabocka et.al／TSS法に基づく効率的な特徴波形および線形評価関数の最適化に関する研究

計算量が低減された時系列データ波形分析装置、方法、及びプログラムを提供する。

一実施形態に係る時系列データ波形分析装置は、特徴波形追加部と、特徴量計算部と、パラメータ更新部と、特徴波形除去部と、を備える。特徴波形追加部は、ラベル付き時系列データに含まれる部分時系列の一部を、特徴波形集合に特徴波形として追加する。特徴量計算部は、１つ又は複数のラベル付き時系列データをランダムに抽出し、抽出されたラベル付き時系列データに対する、特徴波形の特徴量をＴＳＳ法に基づいて計算する。パラメータ更新部は、特徴量に基づいて、特徴波形と、特徴波形に対応する重み係数と、を含むパラメータを確率的勾配降下法に基づいて更新する。特徴波形除去部は、対応する重み係数が０である特徴波形を、特徴波形集合から除去する。評価関数生成部は、特徴波形と、重み係数と、に基づいて、評価関数を生成する。

第１実施形態に係る時系列データ波形分析装置の機能構成の一例を示す図。 ラベル付きデータセットの一例を示す図。 特徴量の計算方法の具体例を示す図。 コンピュータの一例を示す図。 第１実施形態に係る評価関数の生成処理の一例を示すフローチャート。 第２実施形態に係る時系列データ波形分析装置の機能構成の一例を示す図。 第２実施形態に係る部分時系列の追加処理の一例を示すフローチャート。

以下、本発明の実施形態について図面を参照して説明する。

（第１実施形態）
第１実施形態に係る時系列データ波形分析装置（以下、「分析装置」という）について、図１〜図５を参照して説明する。まず、本実施形態に係る分析装置による、時系列データ波形分析方法（以下、「分析方法」という）の概要について説明する。

本実施形態では、まず、分析装置の入力として、学習用の時系列データと、部分時系列集合と、目的関数パラメータと、が用意される。分析装置は、特徴波形として、部分時系列集合からＫ個の部分時系列をランダムに抽出する。分析装置は、抽出した特徴波形と、学習用の時系列データと、に基づいて、各時系列データにおける各特徴波形の特徴量を計算する。分析装置は、得られた特徴量に基づいて、所望の制約条件を有する目的関数（最適化問題）を解き、最適な特徴波形と、最適な重み係数と、を計算する。本実施形態において、目的関数の最適解は、確率的勾配降下法により求められる。分析装置は、こうして得られた最適な特徴波形及び重み係数を含む評価関数を生成し、生成した評価関数により、分析用の時系列データの分類等を行う。

次に、本実施形態に係る分析装置の機能構成について説明する。図１は、本実施形態に係る分析装置の機能構成の一例を示す図である。図１の分析装置は、学習用データ記憶部１と、部分時系列記憶部２と、分析用データ記憶部３と、パラメータ記憶部４と、特徴波形追加部５と、特徴量計算部６と、パラメータ更新部７と、更新終了判定部８と、評価関数生成部９と、時系列データ波形分析部１０と、特徴波形除去部１１と、を備える。

学習用データ記憶部１は、評価関数を生成するための学習用データとして、ラベル付きデータセットＴを記憶する。図２は、ラベル付きデータセットＴの一例を示す図である。図２に示すように、ラベル付きデータセットＴは、Ｎ個のラベル付き時系列データｔ_ｉ（ｉ＝１〜Ｎ）（以下、「時系列データｔ_ｉ」という）の集合である。各時系列データｔ_ｉは、時系列長がＬ_ｉの時系列データに、ラベルが付与されたものである。ラベル付きデータセットＴに含まれる各時系列データｔ_ｉの時系列長Ｌ_ｉは、それぞれ同一であってもよいし、異なってもよい。

図２の例では、各時系列データｔ_ｉに、ラベルとして２値の整数（＋１又は−１）が付与されている。以下では、ラベル＋１が付与された時系列データｔ_ｉの個数はａ個、ラベル−１が付与された時系列データｔ_ｉの個数はｂ個であるものとする（ａ＋ｂ＝Ｎ）。

本実施形態において、分析装置が生成する評価関数は、例えば、ラベルが未知の時系列データに対して、推定されるラベルに応じたスコア（評価値）を算出する関数となる。具体的には、評価関数は、ラベルが＋１と推定される時系列データに対して高い（又は低い）スコアを算出し、ラベルが−１と推定される時系列データに対して低い（又は高い）スコアを算出する関数となる。

また、評価関数は、算出されたスコアに基づいて推定される、時系列データのラベルを出力する関数であってもよい。この場合、評価関数は、時系列データのスコアを閾値と比較し、比較結果に応じたラベルを出力すればよい。

なお、本実施形態において、時系列データｔ_ｉに付与されるラベルは、３種類以上であってもよいし、任意の実数値であってもよい。以下では、学習用データ記憶部１に、図２のラベル付きデータセットＴが記憶されている場合を例に説明する。

部分時系列記憶部２は、部分時系列集合Ｇを記憶する。部分時系列集合Ｇは、Ｑ個の部分時系列ｇ_ｉ（ｉ＝１〜Ｑ）の集合である。各部分時系列ｇ_ｉは、ラベル付きデータセットＴに含まれる各時系列データｔ_ｉから抽出された、時系列長ｒ_ｉを有する任意の部分時系列である。部分時系列集合Ｇに含まれる各部分時系列ｇ_ｉの時系列長ｒ_ｉは、それぞれ同一であってもよいし、異なってもよい。本実施形態において、部分時系列集合Ｇは、ＴＳＳ法における特徴波形の候補として、予め用意される。例えば、部分時系列集合Ｇは学習用データの全部分時系列集合でもよいし、以前本分析装置で似た学習データを分析済みであれば、その際出力された特徴波形を、今回の学習データに対する特徴波形の候補として用意してもよい。

分析用データ記憶部３は、分析用データとして、評価関数による分析対象となるラベル無し時系列データを記憶する。分析装置は、評価関数を利用して、ラベル無し時系列データのスコアやラベルを推定する。これにより、分析装置は、ラベル無し時系列データを、分類等することができる。

パラメータ記憶部４は、評価関数を生成するための各種のパラメータを記憶する。パラメータ記憶部４が記憶するパラメータについて、詳しくは後述する。

特徴波形追加部５（以下、「追加部５」という）は、部分時系列集合ＧからＫ（≪Ｑ）個の部分時系列ｇ_ｊをランダムに抽出し、特徴波形集合Ｓに追加する。特徴波形集合Ｓは、Ｋ個の特徴波形ｓ_ｉ（ｉ＝１〜Ｋ）の集合である。特徴波形集合Ｓに追加された部分時系列ｇは、以降の処理において、特徴波形ｓとして扱われる。

特徴量計算部６は、ラベル付きデータセットＴから時系列データｔをランダムに抽出し、抽出した各時系列データｔの特徴量ベクトルＸを計算する。特徴量ベクトルＸは、時系列データｔに対する、各特徴波形ｓ_ｊの特徴量ｘ_ｓｊを要素として含むベクトルである。

ここで、特徴量ｘの計算方法について説明する。本実施形態における特徴量ｘとは、ＴＳＳ法における、時系列データｔと特徴波形ｓとの間の距離に応じた値のことである。具体的には、特徴量ｘは、時系列長ｒを有する特徴波形ｓと、時系列データｔに含まれる部分時系列ｇ_ｓと、の間の距離を、時系列長ｒで除算した値である。部分時系列ｇ_ｓとは、時系列データｔに含まれる部分時系列のうち、時系列長ｒを有し、かつ、特徴波形ｓとマッチする（特徴波形ｓとの間の平均距離が最も小さい）部分時系列のことである。特徴量ｘは、以下の式で表される。

式（１）におけるｇは時系列データｔに含まれる部分時系列であり、Ｇは時系列データｔの部分時系列集合である。式（１）からわかるように、特徴量ｘは、特徴波形ｓと、時系列データｔに含まれる部分時系列ｇ_ｓと、の間の平均距離に相当する。言い換えると、特徴量ｘは、特徴波形ｓと、時系列データｔと、の間の平均距離の最小値に相当する。

ここで、図３は、特徴量ｘの計算方法の具体例を示す図である。図３の例では、２つの特徴波形ｓ_１，ｓ_２と、時系列データｔと、が用意されている。特徴波形ｓ_１は、３つのデータ（０．１，１．０，０．１：時系列長３）からなり、特徴波形ｓ_２は、４つのデータ（０．１，０．５，０．５，０．１：時系列長４）からなり、時系列データｔは、９つのデータ（０．１，１．１，０．１，０．１，０．１，０．３，０．３，０．１，０．１）からなる。

このとき、時系列データｔの特徴波形ｓ_１に対する部分時系列ｇ_ｓ１は、時系列データｔにおける１番目から３番目のまでのデータからなる部分時系列であり、特徴波形ｓ_１の特徴量ｘ_ｓ１は、０．００３（＝１／３×｛（０．１−０．１）^２＋（１．１−１．０）^２＋（０．１−０．１）^２｝）となる。

また、時系列データｔの特徴波形ｓ_２に対する部分時系列ｇ_ｓ２は、時系列データｔにおける５番目から８番目のまでのデータからなる部分時系列であり、特徴波形ｓ_２の特徴量ｘ_ｓ２は、０．０２（＝１／４×｛（０．１−０．１）^２＋（０．３−０．５）^２＋（０．３−０．５）^２＋（０．１−０．１）^２｝）となる。

本実施形態において、特徴量計算部６は、時系列データｔに対して、上記のようにＫ個の特徴波形ｓ_ｊ（ｊ＝１〜Ｋ）の特徴量ｘ_ｓｊを計算する。したがって、時系列データｔの特徴量ベクトルＸは、Ｋ個の特徴量ｘ_ｓｊを要素として含むベクトルとなる。

パラメータ更新部７は、パラメータ記憶部４に記憶されたパラメータと、時系列データｔの特徴量ベクトルＸと、に基づいて、評価関数ｆの各種のパラメータを更新する。パラメータ更新部７は、更新した各パラメータを、パラメータ記憶部４に記憶させる。

本実施形態において、評価関数ｆは、各特徴波形の特徴量ｘ_ｓｉ（ｉ＝１〜Ｋ）と、各特徴量ｘ_ｓｉの重み係数ｗ_ｓｉと、を含む線形関数で表される。評価関数ｆは、例えば、特徴量ベクトルＸと、重みベクトルＷと、の内積で表される（ｆ＝Ｘ・Ｗ）。重みベクトルＷは、重み係数ｗ_ｓｉ（ｉ＝１〜Ｋ）を要素として含むベクトルのことである。

より詳細には、パラメータ更新部７は、確率的勾配降下法を利用して、重みベクトルＷ（重み係数ｗ）及び特徴波形集合Ｓ（特徴波形ｓ）を更新し、目的関数の解を求める。評価関数ｆは、パラメータ更新部７により更新された重みベクトルＷ及び特徴波形集合Ｓにより構成される。パラメータ更新部７は、特徴波形更新部７１と、重み係数更新部７２と、重み係数正則化部７３と、他変数更新部７４と、を備える。

特徴波形更新部７１は、各特徴波形ｓの勾配を計算し、得られた勾配に基づいて、特徴波形ｓを更新する。これにより、特徴波形集合Ｓが更新される。特徴波形ｓの更新方法について、詳しくは後述する。

重み係数更新部７２は、重み係数ｗの勾配を計算し、得られた勾配に基づいて、重み係数ｗを更新する。これにより、重みベクトルＷが更新される。重み係数ｗの更新方法について、詳しくは後述する。

重み係数正則化部７３（以下、「正則化部７３」という）は、重み係数更新部７２が更新した重みベクトルＷを、目的関数に含まれる正則化条件に基づいて正則化する。例えば、正則化部７３は、Ｌ１空間写像アルゴリズムを利用して、重みベクトルＷを、||Ｗ||_１≦λの空間に写像することができる。このとき、正則化の計算量は、Ｏ（Ｋ）程度となる。Ｌ１空間写像アルゴリズムについて、詳しくは以下の文献に記載の通りである。
参考：J.Duchi et al., Efficient projections onto the l1-ball for learning in high dimensions

なお、重みベクトルＷの正則化方法は上記のものに限られず、任意に選択可能である。分析装置のユーザは、所望の正則化条件を設定することができる。

他変数更新部７４は、目的関数に含まれる他変数の勾配を計算し、得られた勾配に基づいて、他変数を更新する。他変数とは、目的関数に含まれる、特徴波形集合Ｓ及び重みベクトルＷを除く最適化対象のことである。目的関数に他変数が含まれない場合、他変数更新部７４は不要である。他変数の更新方法について、詳しくは後述する。

更新終了判定部８（以下、「判定部８」という）は、パラメータの更新を終了するか判定する。具体的には、判定部８は、更新終了条件が満たされたか判定する。更新終了条件は、例えば、更新回数により設定される。この場合、判定部８は、パラメータ更新部７による更新回数が所定回数に達すると、更新を終了すると判定する。このように、更新終了条件を更新回数により設定することにより、評価関数ｆの生成処理に要する時間を所望の範囲内に設定することができる。

また、更新終了条件は、得られた評価関数ｆの予測精度により設定されてもよい。この場合、判定部８は、学習用データ記憶部１から複数の時系列データｔを取得し、パラメータ更新部７が更新した特徴波形集合Ｓ及び重みベクトルＷにより構成される評価関数ｆにより、取得した時系列データｔのラベルを予測する。判定部８は、予測されたラベルの正解率が所定値以上の場合、更新を終了すると判定する。このように、更新終了条件を予測精度により設定することにより、得られる評価関数ｆの予測精度を保証することができる。

評価関数生成部９は、パラメータ記憶部４に記憶された各パラメータに基づいて、評価関数ｆを生成する。評価関数ｆは、例えば、重みベクトルＷと、特徴波形集合Ｓに対応する特徴量ベクトルＸと、の内積で表される（ｆ＝Ｗ・Ｘ＝ｗ_ｉｘ_ｉ（ｉ＝１〜Ｋ））。

時系列データ波形分析部１０（以下、「分析部１０」という）は、評価関数生成部９が生成した評価関数ｆに基づいて、分析用データを分析する。すなわち、分析部１０は、ラベル無し時系列データの分類等を行う。

特徴波形除去部１１（以下、「除去部１１」という）は、判定部８が、更新を終了しない（更新を継続する）と判定した場合に、対応する重み係数ｗが０である特徴波形ｓを、特徴波形集合Ｓから除去する。除去部１１が除去する特徴波形ｓの数をｋ個とする。除去部１１がｋ個の特徴波形ｓを除去すると、特徴波形集合Ｓに含まれる特徴波形ｓは、（Ｋ−ｋ）個となる。

除去部１１がｋ個の特徴波形ｓを除去すると、追加部５は、部分時系列記憶部２からｋ個の部分時系列ｇを新たに抽出し、特徴波形ｓとして特徴波形集合Ｓに追加する。以降、新たな特徴波形集合Ｓに基づいて、パラメータの更新処理が行われる。

次に、本実施形態に係る分析装置のハードウェア構成について説明する。本実施形態に係る分析装置は、コンピュータ１００により構成される。コンピュータ１００には、サーバ、クライアント、マイコン、及び汎用コンピュータなどが含まれる。

図４は、コンピュータ１００の一例を示す図である。図４のコンピュータ１００は、プロセッサ１０１と、入力装置１０２と、表示装置１０３と、通信装置１０４と、記憶装置１０５と、を備える。プロセッサ１０１、入力装置１０２、表示装置１０３、通信装置１０４、及び記憶装置１０５は、バス１０６により相互に接続されている。

プロセッサ１０１は、コンピュータ１００の制御装置及び演算装置を含む電子回路である。プロセッサ１０１として、例えば、汎用目的プロセッサ、中央処理装置（ＣＰＵ）、マイクロプロセッサ、デジタル信号プロセッサ（ＤＳＰ）、コントローラ、マイクロコントローラ、状態マシン、特定用途向け集積回路、フィールドプログラマブルゲートアレイ（ＦＰＧＡ）、プログラム可能論理回路（ＰＬＤ）、及びこれらの組合せを用いることができる。

プロセッサ１０１は、バス１０６を介して接続された各装置（例えば、入力装置１０２、通信装置１０４、記憶装置１０５）から入力されたデータやプログラムに基づいて演算処理を行い、演算結果や制御信号を、バス１０６を介して接続された各装置（例えば、表示装置１０３、通信装置１０４、記憶装置１０５）に出力する。具体的には、プロセッサ１０１は、コンピュータ１００のＯＳ（オペレーティングシステム）や、時系列データ波形分析プログラム（以下、「分析プログラム」という）などを実行し、コンピュータ１００を構成する各装置を制御する。

分析プログラムとは、コンピュータ１００に、分析装置の上述の各機能構成を実現させるプログラムである。分析プログラムは、一時的でない有形のコンピュータ読み取り可能な記憶媒体 に記憶される。上記の記憶媒体は、例えば、光ディスク、光磁気ディスク、磁気ディスク、磁気テープ、フラッシュメモリ、半導体メモリであるが、これに限られない。プロセッサ１０１が分析プログラムを実行 することにより、コンピュータ１００が分析装置として機能する。

入力装置１０２は、コンピュータ１００に情報を入力するための装置である。入力装置１０２は、例えば、キーボード、マウス、及びタッチパネルであるが、これに限られない。

表示装置１０３は、画像や映像を表示するための装置である。表示装置１０３は、例えば、ＬＣＤ（液晶ディスプレイ）、ＣＲＴ（ブラウン管）、及びＰＤＰ（プラズマディスプレイ）であるが、これに限られない。表示装置１０３は、学習用データ、分析用データ、特徴波形ｓ、評価関数ｆ、及び分析結果など、分析装置により記憶又は生成された任意の情報を表示することができる。

通信装置１０４は、コンピュータ１００が外部装置と無線又は有線で通信するための装置である。通信装置１０４は、例えば、モデム、ハブ、及びルータであるが、これに限られない。

記憶装置１０５は、コンピュータ１００のＯＳや、分析プログラム、分析プログラムの実行に必要なデータ、及び分析プログラムの実行により生成されたデータなどを記憶する記憶媒体である。記憶装置１０５には、主記憶装置と外部記憶装置とが含まれる。主記憶装置は、例えば、ＲＡＭ、ＤＲＡＭ、ＳＲＡＭであるが、これに限られない。また、外部記憶装置は、例えば、ハードディスク、光ディスク、フラッシュメモリ、及び磁気テープであるが、これに限られない。学習用データ記憶部１、部分時系列記憶部２、分析用データ記憶部３、パラメータ記憶部４は、記憶装置１０５上に構築されてもよいし、通信装置１０４を介して接続された外部のサーバ上に構築されてもよい。

なお、コンピュータ１００は、プロセッサ１０１、入力装置１０２、表示装置１０３、通信装置１０４、及び記憶装置１０５を、それぞれ１つ又は複数備えてもよいし、プリンタやスキャナなどの周辺機器を接続されていてもよい。

また、分析装置は、単一のコンピュータ１００により構成されてもよいし、相互に接続された複数のコンピュータ１００からなるシステムとして構成されてもよい。

さらに、分析プログラムは、コンピュータ１００の記憶装置１０５に予め記憶されていてもよいし、コンピュータ１００の外部の記憶媒体に記憶されていてもよいし、インターネット上にアップロードされていてもよい。いずれの場合も、分析プログラムをコンピュータ１００にインストールして実行することにより、分析装置の機能が実現される。

次に、本実施形態に係る分析装置の動作について具体的に説明する。以下では、動作の開始時点で、学習用データ（ラベル付きデータセットＴ）、部分波形集合Ｇ、分析用データ（ラベル無し時系列データ）、目的関数、及び各種のパラメータは、予め用意されているものとする。また、目的関数は、以下の式で表されるものとする。

式（２）において、λは、重みベクトルＷを正則化する空間を規定するパラメータである。λは、０ではない重み係数ｗの個数に影響を与える。||Ｗ||_１≦λは、正則化条件である。また、νは、汎化性能に関連する調整パラメータであり、０以上１以下に設定される。パラメータλ，νは、学習用データに対するスコア付けの精度が最大になるように決定されるのが好ましい。パラメータλ，νは、その設定値がパラメータ記憶部４に記憶される。

また、ρ（＞０）は、最適化されるスカラー変数であり、大きいほどスコア付けの汎用性能が高いことを示すパラメータである。パラメータρは、上述の他変数に相当し、重みベクトルＷ及び特徴波形集合Ｓとともに、目的関数による最適化対象となる。パラメータρ及び重みベクトルＷは、その初期値がパラメータ記憶部４に記憶される。

また、ｌは、損失関数であり、以下の式で表される。

式（３）の損失関数ｌは、ペアワイズヒンジ損失関数である。ペアワイズヒンジ損失関数を利用することにより、２値のラベルに応じたスコアを精度良く算出する評価関数ｆを生成することができる。

また、Ｘ^＋はラベル＋１を付与された時系列データｔ^＋から得られる特徴量ベクトルＸであり、Ｘ⁻はラベル−１を付与された時系列データｔ⁻から得られる特徴量ベクトルＸである。式（２）の目的関数は、できるだけ多くの特徴量ベクトルＸ^＋，Ｘ⁻のペアに対して、ｆ（Ｘ^＋）−ｆ（Ｘ⁻）＞０となるようなｆ（Ｘ）＝Ｗ・Ｘを求めるためのものである。

また、パラメータ記憶部４には、確率的勾配降下法のパラメータである学習率ηの設定値が記憶される。なお、学習率ηは、変数毎に異なる値が設定されてもよい。また、学習率ηは、変数毎の更新回数に応じ、値を変えていっても構わない。例えば、変数の更新回数をｃとしたとき、各変数の更新の際、一般にη=η／ｃ、η=η／ｃ^１／２などとすることで、最適値に収束しやすいことが知られている。

なお、目的関数は、上記のものに限られず、評価関数ｆの目的に応じて、任意に設定可能である。例えば、損失関数ｌとして、ヒンジ損失関数、ロジスティック損失関数、又はε−インセンシティブ損失関数を利用することで、ラベルを予測する評価関数ｆを生成することができる。また、ロジスティック損失関数を利用することにより、ロジスティック回帰問題に適用可能な評価関数ｆを生成することができる。また、ε−インセンシティブ損失関数を利用することにより、時系列データの値を予測する回帰問題に、適用可能な評価関数ｆを生成することができる。これらの損失関数について、詳しくは以下の文献に記載の通りである。
参考：S. Shalev-Shwartz, Y. Singer, and N. Srebro. Pegasos: Primal estimated sub-gradient solver for svm. In Proceedings of the 24th International Conference on Machine Learning, ICML '07, pages 807-814, New York, NY, USA, 2007. ACM.

図５は、本実施形態に係る分析装置による評価関数の生成処理の一例を示すフローチャートである。図５に示すように、生成処理が開始されると、まず、追加部５が、特徴波形集合Ｓを空集合に初期化し、部分時系列記憶部２からＫ個の部分時系列ｇをランダムに抽出し、抽出したＫ個の部分時系列ｇを特徴波形ｓとして特徴波形集合Ｓに追加する（ステップＳ１）。特徴波形集合Ｓには、特徴波形ｓ_ｉ（ｉ＝１〜Ｋ）が含まれる。

次に、特徴量計算部６が、学習用データ記憶部１から、ラベル＋１を付与された時系列データｔ^＋と、ラベル−１を付与された時系列データｔ⁻と、のペアをランダムに抽出する（ステップＳ２）。

特徴量計算部６は、抽出した時系列データｔ^＋，ｔ⁻に対して、それぞれ特徴量ベクトルＸ^＋，Ｘ⁻を計算する（ステップＳ３）。これにより、特徴量ベクトルＸ^＋，Ｘ⁻のペアが得られる。

なお、損失関数ｌが、ヒンジ損失関数、ロジスティック損失関数、又はε−インセンシティブ損失関数である場合には、特徴量計算部６は、時系列データｔを１つランダムに抽出し、その特徴量ベクトルＸを計算すればよい。

続いて、特徴波形更新部７１が、特徴量ベクトルＸ^＋，Ｘ⁻と、パラメータ記憶部４に記憶された各種のパラメータと、に基づいて、特徴波形集合Ｓを更新する（ステップＳ４）。具体的には、特徴波形更新部７１は、特徴波形ｓ_ｉの勾配を、以下の式により計算する。

式（４）において、ｗ_ｉは特徴波形ｓ_ｉに対応する重み係数であり、ｘ^＋は時系列データｔ^＋に対する特徴波形ｓｉの特徴量であり、ｘ⁻は時系列データｔ⁻に対する特徴波形ｓｉの特徴量である。また、Ｉは、指示関数であり、入力が真の場合に１、偽の場合に０を返す。

特徴波形更新部７１は、式（４）により得られた勾配に基づいて、特徴波形ｓ_ｉを以下のように更新する。

特徴波形更新部７１は、以上のように、特徴波形集合Ｓに含まれるＫ個の特徴波形ｓ_ｉ（ｉ＝１〜Ｋ）をそれぞれ更新する。これにより、特徴波形集合Ｓが更新される。

また、重み係数更新部７２が、特徴量ベクトルＸ^＋，Ｘ⁻と、パラメータ記憶部４に記憶された各種のパラメータと、基づいて、重みベクトルＷを更新する（ステップＳ５）。具体的には、重み係数更新部７２は、重みベクトルＷの勾配を、以下の式により計算する。

重み係数更新部７２は、式（６）により得られた勾配に基づいて、重みベクトルＷを以下のように更新する。

そして、正則化部７３が、更新された重みベクトルＷを、正則化条件に従って正則化する。式（２）の例では、正則化部７２は、重みベクトルＷをＬ１正則化する。

また、他変数更新部７４が、特徴量ベクトルＸ^＋，Ｘ⁻と、パラメータ記憶部４に記憶された各種のパラメータと、基づいて、他変数を更新する（ステップＳ６）。式（２）の例では、目的関数には、他変数としてパラメータρが含まれる。そこで、他変数更新部７４は、パラメータρの勾配を、以下の式により計算する。

他変数更新部７４は、式（８）により得られた勾配に基づいて、パラメータρを以下のように更新する。

こうして更新された特徴波形集合Ｓ、重みベクトルＷ、及びパラメータρは、パラメータ記憶部４に記憶される。なお、上記のステップＳ４〜Ｓ６の順番は任意である。また、目的関数に他変数が含まれない場合、ステップＳ６は省略される。

次に、判定部８が、更新終了条件が満たされたか判定する（ステップＳ７）。更新終了条件については、上述の通りである。

更新終了条件が満たされていない場合（ステップＳ７のＮＯ）、除去部１１が、対応する重み係数ｗが０であるｋ個の特徴波形ｓを、特徴波形集合Ｓから除去する（ステップＳ８）。

そして、追加部５が、部分時系列記憶部２からｋ個の部分時系列ｇをランダムに抽出し、抽出したｋ個の部分時系列ｇを、新たな特徴波形ｓとして特徴波形集合Ｓに追加する（ステップＳ９）。

その後、評価関数ｆの生成処理は、ステップＳ２に戻る。以降、更新終了条件が満たされるまで、ステップＳ２〜Ｓ９の処理が繰り返される。

一方、更新終了条件が満たされた場合（ステップＳ７のＹＥＳ）、評価関数生成部９が、その時点でパラメータ記憶部４に記憶されている特徴波形集合Ｓ及び重みベクトルＷに基づいて、評価関数ｆを生成する（ステップＳ１０）。評価関数ｆは、例えば、重みベクトルＷと、特徴波形集合Ｓに対応する特徴量ベクトルＸと、の内積で表される。

上述の通り、本実施形態では、更新処理毎にランダムに選択された時系列データｔに基づいて、各パラメータが更新され、目的関数の解が求められる。この最適化アルゴリズムは、いわゆる、確率的勾配降下法に相当する。すなわち、本実施形態では、確率的勾配降下法を利用して、評価関数ｆが生成される。

分析部１０は、こうして生成された評価関数ｆを利用して、分析用データの分析を行う。これにより、分析部１０は、ラベルが＋１と推定される時系列データに対して高い（又は低い）スコアを付与し、ラベルが−１と推定される時系列データに対して低い（又は高い）スコアを付与することができる。

以上説明した通り、本実施形態に係る分析装置は、Ｋ個の特徴波形ｓに対して複数回の更新処理を実行することにより、評価関数ｆを生成することができる。更新処理の回数をｎ（≪Ｎ）とすると、評価関数ｆを生成するために利用される時系列データｔの数はｎ個である。また、各時系列データの時系列長をＬとする。このとき、本実施形態における、評価関数ｆを生成するための計算量は、Ｏ（Ｋ×ｎ×Ｌ）となる。したがって、本実施形態により、評価関数ｆを生成するための計算量を大幅に低減することが可能となり、計算量が低減された分析装置を実現することができる。

また、本実施形態に係る分析装置は、ランダムに選択されたＫ個の部分時系列ｇを特徴波形ｓとして利用して評価関数ｆを生成する。したがって、Ｋ個の適当な特徴波形候補を用意することなく、任意のラベル付きデータセットＴに基づいて、評価関数ｆを生成することができる。

また、本実施形態に係る分析装置は、対応する重み係数ｗが０である特徴波形ｓが除去された特徴波形集合Ｓに基づいて、評価関数ｆを生成する。これにより、特徴波形集合Ｓに含まれる各特徴波形ｓの独立性を高め、予測精度が高い評価関数ｆを生成することができる。

本実施形態に係る分析装置が分析する時系列データとして、例えば、製造装置等に搭載されている様々なセンサの時系列データが考えられる。分析装置によりセンサの時系列データを分析することにより、製造装置等の異常を検知するための特徴波形及び評価関数を生成することができる。この特徴波形及び評価関数を用いて、現在のセンサの時系列データを分析することにより、製造装置等の異常をリアルタイムに検知することができる。本実施形態に係る分析装置を利用することにより、時系列データの時系列長が非常に長い場合（例えば、サンプリングレートが高い場合や、処理時間が長い場合）であっても、高速に分析することができる。

また、エネルギー需要量等も時系列データで表される。そこで、この時系列データに、高需要状況、低需要状況、平常状況などの、需要状況を示すラベルを付与した学習用データを用意し、この学習用データを分析装置で分析することで、エネルギー需要量等の需要状況を推定するための特徴波形及び評価関数を生成することができる。この特徴波形及び評価関数を用いて、現在の時系列データにラベルを付与することにより、現在の需要状況を推定し、電力供給量の調整に活用することができる。本実施形態に係る分析装置を利用することにより、需要量等の測定スパン（１日、１週間など）が長い場合であっても、高速に分析することができる。

また、加速度計などの生体センサの時系列データに、立っている、歩いている、眠っているなどの、人間の状態や行動を示すラベルを付与した学習用データを用意し、この学習用データを分析装置で分析することで、人間の状態や行動を推定するための特徴波形及び評価関数を生成することができる。この特徴波形及び評価関数を用いて、現在の時系列データにラベルを付与することにより、現在の状態や行動を推定し、医療行為などに活用することができる。本実施形態に係る分析装置を利用することにより、サンプリングレートが非常に高い生体センサであっても、高速に分析することができる。

なお、本実施形態において、パラメータ更新部７による各種のパラメータの更新は、バッチ処理により実行されてもよい。この場合、ステップＳ２において、特徴量計算部６は、複数のペアを抽出し、抽出された各ペアに対応する特徴量ベクトルＸ^＋，Ｘ⁻のペアを計算する。そして、ステップＳ４〜Ｓ６において、パラメータ更新部７は、特徴量ベクトルＸ^＋，Ｘ⁻の各ペアに基づいて、特徴波形集合Ｓ、重みベクトルＷ、及び他変数をそれぞれ更新する。すなわち、１回の更新処理で、各パラメータがそれぞれ複数回更新される。これにより、更新処理の回数を低減することができる。

（第２実施形態）
第２実施形態に係る分析装置について、図６及び図７を参照して説明する。本実施形態では、特徴波形集合Ｓは、優先度に応じたｋ個の部分時系列ｇを、特徴波形ｓとして追加される。図６は、本実施形態に係る分析装置の機能構成の一例を示す図である。図６の分析装置は、優先度計算部１２を備える。図６の分析装置の他の構成は、図１と同様である。また、本実施形態に係る分析装置のハードウェア構成は、第１実施形態と同様である。

優先度計算部１２は、部分時系列ｇの優先度を計算する。部分時系列ｇの優先度とは、除去部１１が特徴波形集合Ｓから除去した特徴波形ｓに対する、部分時系列ｇの非類似度を示す値である。優先度計算部１２は、非類似度が高いほど優先度を高く計算する。言い換えると、部分時系列ｇの優先度は、除去された特徴波形ｓに対する非類似度が高いほど高くなる。

非類似度として、例えば、特徴量を利用することができる。特徴量が大きいほど、非類似度は高くなる。この場合、優先度計算部１２は、部分時系列ｇに対する、除去された特徴波形ｓの特徴量を計算し、得られた特徴量が大きい程、優先度を高く計算すればよい。

除去された特徴波形ｓが複数ある場合、優先度計算部１２は、部分時系列ｇの優先度として、各特徴波形ｓに対する部分時系列ｇの優先度の最大値、最小値、又は平均値などを計算すればよい。

本実施形態において、追加部５は、優先度計算部１２が計算した優先度が高い部分時系列ｇを、特徴波形ｓとして特徴波形集合Ｓに追加する。追加方法について、詳しくは後述する。

次に、本実施形態に係る分析装置の動作について説明する。図７は、本実施形態に係る分析装置による部分時系列ｇの追加処理の一例を示すフローチャートである。図７のフローチャートは、図５のステップＳ９の内部処理に相当する。なお、本実施形態における評価関数ｆの生成処理における、ステップ９以外の処理は図５と同様である。

本実施形態では、除去部１１がｋ個の特徴波形ｓを特徴波形集合Ｓから除去すると（ステップＳ８）、図７に示すように、追加部５が、部分時系列記憶部２から１つの部分時系列ｇをランダムに抽出する（ステップＳ９１）。次に、優先度計算部１２が、抽出された部分時系列ｇの優先度を計算する（ステップＳ９２）。続いて、追加部５が、部分時系列ｇの優先度が予め設定された閾値以上であるか判定する（ステップＳ９３）。

優先度が閾値未満の場合（ステップＳ９３のＮＯ）、処理はステップＳ９１に戻る。

一方、優先度が閾値以上である場合（ステップＳ９３のＹＥＳ）、追加部５は、ステップＳ９１で抽出した部分時系列ｇを、特徴波形ｓとして特徴波形集合Ｓに追加する（ステップＳ９４）。

その後、追加部５は、特徴波形集合Ｓにｋ個の部分時系列ｇを追加したか判定する（ステップＳ９５）。ｋ個の部分時系列ｇが追加されていない場合（ステップＳ９５のＮＯ）、処理はステップＳ９１に戻る。一方、ｋ個の部分時系列ｇが追加されている場合（ステップＳ９５のＹＥＳ）、部分時系列ｇの追加処理が終了する。その後、優先度が高い部分時系列ｓが追加された特徴波形集合Ｓに基づいて、評価関数ｆが生成される。

以上説明した通り、本実施形態に係る分析装置は、対応する重み係数ｗが０である特徴波形ｓの代わりに、優先度が高い特徴波形ｓ（部分時系列ｇ）が追加された特徴波形集合Ｓに基づいて、評価関数ｆを生成する。これにより、特徴波形集合Ｓに含まれる各特徴波形ｓの独立性をより高め、予測精度が高い評価関数ｆを生成することができる。

なお、部分時系列ｇの抽出及び優先度の計算は、バッチ処理により実行されてもよい。この場合、ステップＳ９１において、追加部５が、複数の部分時系列ｇを抽出しステップＳ９２において、優先度計算部１２が、抽出された複数の部分時系列ｇの優先度を計算する。

また、本実施形態において、優先度計算部１２は、部分時系列記憶部２に記憶された部分時系列ｇの一部又は全部の優先度を、予め計算してもよい。この場合、追加部５は、優先度が計算された部分時系列ｇの中から、優先度が高い順にｋ個の部分時系列ｇを特徴波形集合Ｓに追加すればよい。このように部分時系列ｇを追加することにより、優先度が高い順に部分時系列ｇを特徴波形集合Ｓに追加することができる。

なお、本発明は上記各実施形態そのままに限定されるものではなく、実施段階ではその要旨を逸脱しない範囲で構成要素を変形して具体化できる。また、上記各実施形態に開示されている複数の構成要素を適宜組み合わせることによって種々の発明を形成できる。また例えば、各実施形態に示される全構成要素からいくつかの構成要素を削除した構成も考えられる。さらに、異なる実施形態に記載した構成要素を適宜組み合わせてもよい。

１：学習用データ記憶部、２：部分時系列記憶部、３：分析用データ記憶部、４：パラメータ記憶部、５：特徴波形追加部、６：特徴量計算部、７：パラメータ更新部、８：更新終了判定部、９：評価関数生成部、１０：時系列データ波形分析部、１１：特徴波形除去部、１２：優先度計算部、７１：特徴波形更新部、７２：重み係数更新部、７３：重み係数正則化部、７４：他変数更新部、１００：コンピュータ、１０１：プロセッサ、１０２：入力装置、１０３：表示装置、１０４：通信装置、１０５：記憶装置、１０６：バス

Claims (9)

1. ラベル付き時系列データに含まれる部分時系列の一部を、特徴波形集合に特徴波形として追加する特徴波形追加部と、
   １つ又は複数のラベル付き時系列データをランダムに抽出し、抽出された前記ラベル付き時系列データに対する、前記特徴波形の特徴量をＴＳＳ法に基づいて計算する特徴量計算部と、
   前記特徴量に基づいて、前記特徴波形と、前記特徴波形に対応する重み係数と、を含むパラメータを確率的勾配降下法に基づいて更新するパラメータ更新部と、
   対応する前記重み係数が０である前記特徴波形を、前記特徴波形集合から除去する特徴波形除去部と、
   前記特徴波形と、前記重み係数と、に基づいて、評価関数を生成する評価関数生成部と、
   を備える時系列データ波形分析装置。
2. 前記パラメータ更新部は、前記特徴量に基づいて、前記パラメータの勾配を計算し、前記勾配に基づいて、前記パラメータを更新する
   請求項１に記載の時系列データ波形分析装置。
3. 前記特徴量は、前記特徴波形と、前記ラベル付き時系列データと、の間の平均距離の最小値である
   請求項１又は請求項２に記載の時系列データ波形分析装置。
4. 前記評価関数に基づいて、ラベル無し時系列データを分析する時系列データ波形分析部を備える
   請求項１乃至請求項３のいずれか１項に記載の時系列データ波形分析装置。
5. 前記パラメータ更新部は、所定の正則化条件に従って、前記重み係数を正則化する
   請求項１乃至請求項４のいずれか１項に記載の時系列データ波形分析装置。
6. 更新回数又は前記評価関数の精度に基づいて、前記パラメータの更新を終了するか判定する更新終了判定部を備える
   請求項１乃至請求項５のいずれか１項に記載の時系列データ波形分析装置。
7. 対応する前記重み係数が０である前記特徴波形と、前記ラベル付き時系列データに含まれる部分時系列と、の非類似度に応じた優先度を計算する優先度計算部を備える
   請求項１乃至請求項６のいずれか１項に記載の時系列データ波形分析装置。
8. ラベル付き時系列データに含まれる部分時系列の一部を、特徴波形集合に特徴波形として追加する工程と、
   １つ又は複数のラベル付き時系列データをランダムに抽出し、抽出された前記ラベル付き時系列データに対する、前記特徴波形の特徴量をＴＳＳ法に基づいて計算する工程と、
   前記特徴量に基づいて、前記特徴波形と、前記特徴波形に対応する重み係数と、を含むパラメータを確率的勾配降下法に基づいて更新する工程と、
   対応する前記重み係数が０である前記特徴波形を、前記特徴波形集合から除去する工程と、
   前記特徴波形と、前記重み係数と、に基づいて、評価関数を生成する工程と、
   を含む時系列データ波形分析方法。
9. ラベル付き時系列データに含まれる部分時系列の一部を、特徴波形集合に特徴波形として追加する工程と、
   １つ又は複数のラベル付き時系列データをランダムに抽出し、抽出された前記ラベル付き時系列データに対する、前記特徴波形の特徴量をＴＳＳ法に基づいて計算する工程と、
   前記特徴量に基づいて、前記特徴波形と、前記特徴波形に対応する重み係数と、を含むパラメータを確率的勾配降下法に基づいて更新する工程と、
   対応する前記重み係数が０である前記特徴波形を、前記特徴波形集合から除去する工程と、
   前記特徴波形と、前記重み係数と、に基づいて、評価関数を生成する工程と、
   をコンピュータに実行させる時系列データ波形分析プログラム。

Images