(57)【特許請求の範囲】
【請求項１】
行列として表される観測データに関するヘッセ行列の行列式の近似値を計算する計算部と、
前記行列式の近似値を利用して隠れ変数の変分確率を計算する変分確率計算部と、
変分分布に基づいて隠れ状態を削除する隠れ状態削除部と、
周辺化対数尤度関数を完全変数に対する推定量に関してラプラス近似した近似量の下界として定義される基準値に対してパラメータを最適化するとともに、基準値を計算するパラメータ最適化部と、
前記基準値が収束したか否かを判定する収束判定部とを
備えることを特徴とするモデル推定装置。
【請求項２】
計算部がヘッセ行列の行列式の近似値を計算し、変分確率計算部が隠れ変数の変分確率を計算し、隠れ状態削除部が隠れ状態を削除し、パラメータ最適化部がパラメータを最適化し、計算部がヘッセ行列の行列式の近似値を計算し、パラメータ最適化部が基準値を計算し、収束判定部が、前記基準値が収束したか否かを判定するループ処理を、収束判定部が、前記基準値が収束したと判定するまで繰り返す
請求項１に記載のモデル推定装置。
【請求項３】
モデル推定装置は、隠れ属性モデル推定装置である請求項１または請求項２に記載のモデル推定装置。
【請求項４】
コンピュータが、
行列として表される観測データに関するヘッセ行列の行列式の近似値を計算し、
前記行列式の近似値を利用して隠れ変数の変分確率を計算し、
変分分布に基づいて隠れ状態を削除し、
周辺化対数尤度関数を完全変数に対する推定量に関してラプラス近似した近似量の下界として定義される基準値に対してパラメータを最適化し、
ヘッセ行列の行列式の近似値を計算し、
前記基準値を計算し、
前記基準値が収束したか否かを判定する
ことを特徴とするモデル推定方法。
【請求項５】
コンピュータが、
ヘッセ行列の行列式の近似値を計算し、隠れ変数の変分確率を計算し、隠れ状態を削除し、パラメータを最適化し、ヘッセ行列の行列式の近似値を計算し、前記基準値を計算し、前記基準値が収束したか否かを判定するループ処理を、前記基準値が収束するまで繰り返す
請求項４に記載のモデル推定方法。
【請求項６】
モデル推定方法は、隠れ属性モデル推定方法である請求項４または請求項５に記載のモデル推定方法。
【請求項７】
コンピュータに、
行列として表される観測データに関するヘッセ行列の行列式の近似値を計算する計算処理、
前記行列式の近似値を利用して隠れ変数の変分確率を計算する変分確率計算処理、
変分分布に基づいて隠れ状態を削除する隠れ状態削除処理、
周辺化対数尤度関数を完全変数に対する推定量に関してラプラス近似した近似量の下界として定義される基準値に対してパラメータを最適化するパラメータ最適化処理、
前記基準値を計算する基準値計算処理、および、
前記基準値が収束したか否かを判定する収束判定処理
を実行させるためのモデル推定プログラム。
【請求項８】
コンピュータに、
計算処理、変分確率計算処理、隠れ状態削除処理、パラメータ最適化処理、計算処理、基準値計算処理、および収束判定処理のループ処理を、基準値が収束したと判定されるまで繰り返し実行させる
請求項７に記載のモデル推定プログラム。
【請求項９】
モデル推定プログラムは、隠れ属性モデル推定プログラムである請求項７または請求項８に記載のモデル推定プログラム。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、多変量データの隠れ属性モデルを推定するモデル推定装置、モデル推定方法およびモデル推定プログラムに関し、特に、モデル事後確率を近似し、その下界を最大化する事によって多変量データの隠れ属性モデルを推定するモデル推定装置、モデル推定方法およびモデル推定プログラムに関する。
【背景技術】
【０００２】
自動車から取得されるセンサデータ、健康診断の検査値、電力需要履歴などに代表されるデータには、背後に、観測されない状態（例えば、自動車の故障状態、生活習慣、次の日の気象条件など）が存在する。このようなデータを分析するためには、観測されない変数の存在を仮定した隠れ変数モデルが重要な役割を果たす。隠れ変数は、上記のような観測に重大な影響を与える因子を表現する。隠れ変数モデルを利用したデータ分析は、産業上重要な多くの分野に適用される。例えば、自動車から取得されるセンサデータの分析によって、自動車の故障原因を解析して素早い修理を実現することが考えられる。また、健康診断の検査値の分析によって、疾患のリスクの推定および疾患の予防を実現することが考えられる。また、電力需要履歴の分析によって、電力の需要を予測して過不足に備えられるようにすることが考えられる。
【０００３】
隠れ変数モデルの最も代表的な例として、混合分布モデルが挙げられる。混合分布モデルは、観測データが複数の性質をもったグループから独立に観測されていると仮定し、そのグループ構造を隠れ変数として表現するモデルである。混合分布モデルでは、各グループが独立であるという仮定に基づく。しかし、現実のデータは複数の要因が絡まりあって観測されることが多い。そこで、混合分布モデルを拡張した、隠れ属性モデルが提案されている（例えば、非特許文献１参照）。このモデルは、各観測データの背後に複数の因子（属性）の存在を仮定し、その組み合わせによって観測が得られるという仮定に基づく。
【０００４】
隠れ属性モデルを学習するためには、隠れ状態数、観測確率分布の種類、及び分布パラメータを決定する必要がある。特に、隠れ状態数や観測確率の種類を決定する問題は、一般的に「モデル選択問題」や「システム同定問題」と呼ばれ、信頼性のあるモデルを構築するために極めて重要な問題である。そのための技術が種々提案されている。
【０００５】
例えば、非特許文献１では、隠れ状態を決定する方法として、変分ベイズ法によって、変分自由エネルギーを最大化する方法が提案されている。以下、この方法を第１の公知技術と記す。
【０００６】
また、例えば、非特許文献１では、隠れ状態を決定する方法として、階層Dirichlet過程事前分布を用いたノンパラメトリックベイズ法が提案されている。以下、この方法を第２の公知技術と記す。
【０００７】
混合モデルでは、隠れ変数が独立であり、隠れ変数に対してパラメータが独立である。また、隠れマルコフモデルでは、隠れ変数が時間依存性を持ち、隠れ変数に対してパラメータが独立である。そして、非特許文献２、非特許文献３では、混合モデルおよび隠れマルコフモデルに適用される技術として、因子化漸近ベイズ推論（Factorized Asymptotic Bayesian Inference）という技術が提案されている。この技術は、変分ベイズ法やノンパラメトリックベイズ法に対して、速度・精度の面で優位性を持つ。
【０００８】
また、完全周辺尤度関数を近似し、その下界を最大化することが、非特許文献２、非特許文献３に記載されている。
【先行技術文献】
【非特許文献】
【０００９】
【非特許文献１】Thomas L. Griffiths  and Zoubin Ghahramani, "Infinite Latent Feature Models and the Indian Buffet Process", Technical Report 2005-001, Gatsby Computational Neuroscience Unit, 2005.
【非特許文献２】Ryohei Fujimaki, Satoshi Morinaga: "Factorized Asymptotic Bayesian Inference for Mixture Modeling". Proceedings of the the fifteenth international conference on Artificial Intelligence and Statistics (AISTATS), 2012
【非特許文献３】Ryohei Fujimaki, Kohei Hayashi: "Factorized Asymptotic Bayesian Hidden Markov Model". Proceedings of the 25th international conference on machine learning (ICML), 2012
【発明の概要】
【発明が解決しようとする課題】
【００１０】
第１の公知技術は、周辺化尤度関数の下界を最大化する際に、変分分布上における隠れ状態と分布パラメータの独立性を仮定する。そのため、第１の公知技術では、周辺化尤度の近似精度が悪くなるという問題がある。
【００１１】
第２の公知技術は、モデルが複雑になるため計算量が非常に大きくなるという問題や、入力パラメータによって結果が大きく変わるという問題がある。
【００１２】
また、非特許文献２、非特許文献３等に記載された技術では、本質的に隠れ変数に関するパラメータの独立性が重要である。そのため、隠れ属性モデルのように隠れ変数に対してパラメータが依存関係をもつモデルに直接、因子化漸近ベイズ推論を適用することができなかった。
【００１３】
本発明は、因子化漸近ベイズ推論に基づいて、隠れ属性モデルに関するモデル選択問題を解決する隠れ属性モデル推定装置、隠れ属性モデル推定方法および隠れ属性モデル推定プログラムを提供することを目的とする。
【課題を解決するための手段】
【００１４】
本発明によるモデル推定装置は、行列として表される観測データに関するヘッセ行列の行列式の近似値を計算する計算部と、その行列式の近似値を利用して隠れ変数の変分確率を計算する変分確率計算部と、変分分布に基づいて隠れ状態を削除する隠れ状態削除部と、周辺化対数尤度関数を完全変数に対する推定量に関してラプラス近似した近似量の下界として定義される基準値に対してパラメータを最適化するとともに、基準値を計算するパラメータ最適化部と、基準値が収束したか否かを判定する収束判定部とを備えることを特徴とする。
【００１５】
また、本発明によるモデル推定方法は、コンピュータが、行列として表される観測データに関するヘッセ行列の行列式の近似値を計算し、その行列式の近似値を利用して隠れ変数の変分確率を計算し、変分分布に基づいて隠れ状態を削除し、周辺化対数尤度関数を完全変数に対する推定量に関してラプラス近似した近似量の下界として定義される基準値に対してパラメータを最適化し、ヘッセ行列の行列式の近似値を計算し、基準値を計算し、基準値が収束したか否かを判定することを特徴とする。
【００１６】
また、本発明によるモデル推定プログラムを記録したコンピュータ読み取り可能な記録媒体は、コンピュータに、行列として表される観測データに関するヘッセ行列の行列式の近似値を計算する計算処理、その行列式の近似値を利用して隠れ変数の変分確率を計算する変分確率計算処理、変分分布に基づいて隠れ状態を削除する隠れ状態削除処理、周辺化対数尤度関数を完全変数に対する推定量に関してラプラス近似した近似量の下界として定義される基準値に対してパラメータを最適化するパラメータ最適化処理、基準値を計算する基準値計算処理、および、基準値が収束したか否かを判定する収束判定処理を実行させるためのモデル推定プログラムを記憶する。
【発明の効果】
【００１７】
本発明によれば、因子化漸近ベイズ推論に基づいて、隠れ属性モデルに関するモデル選択問題を解決することができる。
【図面の簡単な説明】
【００１８】
【図１】図１は、本発明の隠れ属性モデル推定装置の構成例を示すブロック図である。
【図２】図２は、本発明の処理経過の例を示すフローチャートである。
【図３】図３は、本発明の概要を示すブロック図である。
【発明を実施するための形態】
【００１９】
本発明の貢献を明確化するために、まず隠れ属性モデルについての説明と、因子化漸近ベイズ推論がなぜ直接、隠れ属性モデルに適用できないかの問題点を、具体的に説明する。
【００２０】
以下の説明では、観測データをXとする。ここで、Xは、サンプル数をNとし、次元数をDとしたときに、N行D列の行列として表される。そして、行列の第n行第d列の要素を下付けの添え字ndで表す。例えば、Xの第n行第d列はXndである。
【００２１】
隠れ属性モデルでは、Xを２つの行列（AとZとする。）の積として表現されることを仮定する。すなわち、加法ノイズ項をEとすれば、X = ZA + Eである。ここで、A（サイズはKかけるD）は連続値をとる重みパラメータである。Zはバイナリ値をとる隠れ変数（サイズはNかけるK）である。Kは隠れ状態数を表す。なお、以下の説明では、Eの分布が正規分布に従うと仮定して説明をする。ただし、指数分布族などより広い分布クラスにおいても同様の議論が成立する。
【００２２】
XとZに関する同時確率分布を考えると、同時分布は、以下の式１に示すように分解される。
【００２３】
p(X, Z|j) = p(X | Z, jx) p(Z|jz)  ・・・  （式１）
【００２４】
ただし、jは同時分布のパラメータであり、jxとjzは各分布のパラメータである。加法ノイズ項Eが独立な正規分布に従うと仮定した場合、jxはA及び共分散行列
【数１】
0006398991.tif 000002
であり、p(X | Z, jx)は、平均がZAで共分散行列
【数２】
0006398991.tif 000003
の正規分布となる。Iは単位行列である。ここで、Xndは、平均が
【数３】
0006398991.tif 000004
で分散
【数４】
0006398991.tif 000005
の正規分布にしたがう。ここで重要な点は、パラメータAが隠れ変数のインデックスkに対して相互に依存している点である。
【００２５】
対比のために、混合分布の例を説明する。混合分布は、Xnの分布を
【数５】
0006398991.tif 000006

と表す。ただし、a_kは混合比である。pkはk番目の隠れ変数に対応する分布であり、jkはそのパラメータである。隠れ属性モデルとは異なり、混合分布では隠れ変数のインデックスkに対して、パラメータjkが相互に独立である事がわかる。
【００２６】
このパラメータの依存性の問題を、非特許文献２を例に説明する。非特許文献２では、観測変数と隠れ変数の同時分布をラプラス近似し、対数同時尤度関数を近似する。非特許文献２に記載されている（５）式がその近似式である。ここで重要な点は、隠れ変数が与えられたときに、対数尤度関数の二回微分行列（以下、単にヘッセ行列と呼ぶ。）が、ブロック対角性をもつ点である。換言すれば、重要な点は、各隠れ変数に対応するパラメータが、同一の隠れ変数に対しては依存性を持つが、異なる隠れ変数に対しては独立である場合には、ヘッセ行列の非対角ブロックが全て０になるという点である。この性質によって、pk(Xn|jk)をkに対して別々にラプラス近似し、各因子化情報量基準（非特許文献２に記載されている（１０）式）が導出され、その下界を最大化するアルゴリズムである因子化漸近ベイズ推論アルゴリズムが導出される（非特許文献２のSection 4参照）。しかし、上述のように、隠れ属性モデルでは、隠れ変数に対するパラメータが依存しているため、ヘッセ行列がブロック対角性を持たない。そのため、因子化漸近ベイズ推論の手順を隠れ属性モデルに直接適用することができない点に問題がある。本発明は、公知技術とは異なるヘッセ行列（の行列式）を近似する処理手順を導入することによってこの問題を解決している点で、前述の先行技術と本質的に相違する。
【００２７】
以下、本発明の実施形態を図面を参照して説明する。
【００２８】
図１は、本発明の隠れ属性モデル推定装置の構成例を示すブロック図である。本発明の隠れ属性モデル推定装置１００は、データ入力装置１０１と、隠れ状態数設定部１０２と、初期化処理部１０３と、隠れ変数変分確率計算処理部１０４と、情報量基準近似処理部１０５と、隠れ状態選択処理部１０６と、パラメータ最適化処理部１０７と、最適性判定処理部１０８と、モデル推定結果出力装置１０９とを備えている。隠れ属性モデル推定装置１００には、入力データ１１１が入力される。隠れ属性モデル推定装置１００は、入力データ１１１に対して隠れ属性モデル最適化し、モデル推定結果１１２として出力する。
【００２９】
データ入力装置１０１は、入力データ１１１を入力するための装置である。観測確率の種類や、隠れ状態数の候補値など、モデルの推定に必要なパラメータが入力データ１１１として同時にデータ入力装置１０１に入力される。
【００３０】
隠れ状態数設定部１０２は、モデルの隠れ状態数Kを、入力データ１１１として入力された最大値Kmaxに設定する。すなわち、隠れ状態数設定部１０２によって、K=Kmaxと設定される。
【００３１】
初期化処理部１０３は、推定のための初期化処理を実施する。なお、初期化は任意の方法によって実施することが可能である。例としては、各観測確率のパラメータjをランダムに設定する方法や、隠れ変数の変分確率をランダムに設定する方法が挙げられる。
【００３２】
隠れ変数変分確率計算処理部１０４は、隠れ変数の変分確率を計算する。ここで、パラメータjは初期化処理部１０３あるいはパラメータ最適化処理部１０７で算出されているため、隠れ変数変分確率計算処理部１０４は、その値を利用する。隠れ変数変分確率計算処理部１０４は、以下のように定義する最適化基準Aを最大化することによって変分確率を計算する。最適化基準Aとは、周辺化対数尤度関数を、完全変数に対する推定量（例えば最尤推定量や最大事後確率推定量）に関してラプラス近似した近似量の下界として定義される。
【００３３】
情報量基準近似処理部１０５では、隠れ変数変分確率計算処理部１０４やパラメータ最適か処理部１０７で必要な、ヘッセ行列の行列式の近似処理を行う。情報量基準近似処理部１０５の具体的な処理については以下に説明する。
【００３４】
隠れ変数変分確率計算処理部１０４と情報量基準近似処理部１０５の処理について、以下で具体的に説明する。
【００３５】
本発明では、ベイズ推論にしたがって、周辺化対数尤度を最大化することで、モデル及びパラメータを最適化する。ただし、周辺化対数尤度を直接最適化することは難しいため、まず周辺化対数尤度を、式２のように変形する。
【００３６】
【数６】
0006398991.tif 000007
【００３７】
ただし、Mはモデルであり、q(Z)はZに関する変分分布である。また、max_qは、qに関する最大値であることを表す。ここで、同時周辺化尤度p(X, Z|M)は、パラメータに対する積分形式で、式３のように変形可能である。
【００３８】
【数７】
0006398991.tif 000008
【００３９】
まず、混合分布モデルの同時分布
【数８】
0006398991.tif 000009

を考える。ここで、
【数９】
0006398991.tif 000010

である点に注意すべきである。log p(X, Z|j)に関するヘッセ行列は、jz、jk(k=1, ..., K)に関してブロック対角となる。そのため、log p(X, Z|j)をp(X, Z|j)の最尤推定量周りでテイラー展開して、３次以上の項を無視することによって、log p(X, Z|j)は、以下の式４に示すように近似される。
【００４０】
【数１０】
0006398991.tif 000011
【００４１】
この式は、非特許文献２に記載された（５）式に対応する。ここで、Fz及びFkは、p(Z|jz)及びpk(Xn|jk)のヘッセ行列をそれぞれN及び
【数１１】
0006398991.tif 000012
で割った行列であり、p(X, Z|j)のヘッセ行列のブロック対角項に対応する。式４の近似を（式３）に代入することで、log p(X, Z|M)の近似式として、以下に示す式５が得られる。
【００４２】
【数１２】
0006398991.tif 000013
【００４３】
この式は、非特許文献２に記載された（９）式に対応する。ただし、detは引数の行列式を表す。また、Dz及びDkは、jz及びjkの次元を表すものとする。ここで、Nの極限を考えると、
【数１３】
0006398991.tif 000014
log det(Fz)、log det(Fk)は相対的に小さくなるため、無視することが可能である。式５からこれらに関する項を無視して、式１へ代入すると、因子化情報量基準として、以下の式６が得られる。
【００４４】
【数１４】
0006398991.tif 000015
【００４５】
log p(X, Z|j’)はデータへのフィッティングを表し、
【数１５】
0006398991.tif 000016
はモデルの複雑性を表している。
【００４６】
非特許文献２で提案される因子化漸近ベイズ推論では、j’を任意のjに置き換え、
【数１６】
0006398991.tif 000017
を
【数１７】
0006398991.tif 000018
なる下界に置き換え、以下の式７のようにモデルを推定する。
【００４７】
【数１８】
0006398991.tif 000019
【００４８】
次に、上記の手順を隠れ属性モデルに適用した例を説明する。隠れ属性モデルに関する同時分布
【数１９】
0006398991.tif 000020
について、log p(X, Z|j)を最尤推定量周りでテイラー展開して３次以上の項を無視すると、以下の式８に示す近似式が得られる。
【００４９】
【数２０】
0006398991.tif 000021
【００５０】
ただし、
【数２１】
0006398991.tif 000022
であり、Fdは
【数２２】
0006398991.tif 000023
のjdに関するヘッセ行列である。
【００５１】
ここで、上記で説明した既存技術の処理手順に従うと、以下の式９が得られる。すなわち、式８を式３に代入し、
【数２３】
0006398991.tif 000024
log det(Fz)、log det(Fd)を相対的に小さいとして無視すると、p(X,Z|M)の近似として、以下の式９が得られる。
【００５２】
【数２４】
0006398991.tif 000025
【００５３】
なお、Dd=K+1は、jdの次元数である。また、情報量基準は、式１０で表される。
【００５４】
【数２５】
0006398991.tif 000026
【００５５】
式６に対するモデル推定処理と式１０のモデル推定処理とで本質的に異なるのは、式６の
【数２６】
0006398991.tif 000027
という項が、式１０ではDd log Nになり、モデル複雑性が隠れ変数に依存しなくなっている点である。より具体的に説明する。非特許文献２で提案されている因子化漸近ベイズ推論では、モデルの複雑性が隠れ変数に依存していることによって、不要な隠れ状態の削除やモデルの同定性など、理論的に優れた性質を持つ。なお、不要な隠れ状態の削除については、非特許文献２の" Section 4.4 Shrinkage Mechanism"で説明されている。また、モデルの同定性については、非特許文献２の" Section 4.5 Identifiability"で説明されている。しかし、隠れ属性モデルに関して上述のように得られる式１０では、そのような性質が失われてしまう。
【００５６】
そこで、本発明で提案する隠れ変数変分確率計算処理部１０４と情報量基準近似処理部１０５では、以下のような手順によって情報量基準を算出する。
【００５７】
非特許文献２による手順では、式９のlog det(Fd)を、漸近的に小さいとして、以下に示すように近似している。
【００５８】
【数２７】
0006398991.tif 000028
【００５９】
これに対し、情報量基準近似処理部１０５は、log det(Fd)を以下の式１１のように近似する。
【００６０】
【数２８】
0006398991.tif 000029
【００６１】
式１１を式９へ代入し、
【数２９】
0006398991.tif 000030
log det(Fz)を漸近的に小さいとして無視すると、情報量基準として式１０に代えて、式１２が得られる。
【００６２】
【数３０】
0006398991.tif 000031
【００６３】
式１２は式６と同じ形式である。そして、式１２によれば、モデルの複雑性が隠れ変数に依存していることによって、不要な隠れ状態の削除やモデルの同定性など、理論的に優れた性質を持つ基準となっている。重要な点は、隠れ属性モデルに対して、式１２の基準を得るためには、情報量基準近似処理部１０５による処理（すなわち、式１１の近似処理）が不可欠であり、これは公知技術にはない、本発明の特徴といえる。
【００６４】
隠れ状態選択処理部１０６は、隠れ状態のうち、小さな状態をモデルから除去する。具体的には、k番目の隠れ状態について、
【数３１】
0006398991.tif 000032
が、入力データ１１１として設定された閾値を下回った場合には、その状態をモデルから削除する。
【００６５】
パラメータ最適化処理部１０７は、隠れ変数の変分確率を固定した上で、最適化基準Aに対してjを最適化する。なお、最適化基準Aのjに関する項は、隠れ状態の変分分布によって重み付けされた対数同時尤度関数であり、任意の最適化アルゴリズムによって最適化することが可能である。例えば、上述の例における正規分布では、パラメータ最適化処理部１０７は、平均場近似法によって、パラメータを最適化することができる。また、パラメータ最適化処理部１０７は、最適化されたパラメータに対して最適化A基準を同時に計算する。この際に、パラメータ最適化処理部１０７は、上述の情報量基準近似処理部１０５の近似計算を利用する。すなわち、パラメータ最適化処理部１０７は、式１１によるヘッセ行列の行列式の近似結果を利用する。
【００６６】
最適性判定処理部１０８は、最適化基準Aの収束を判定する。収束の判定は、最適化基準Aの絶対変化量や相対変化量などに閾値を設け、閾値を用いることで実現できる。
【００６７】
モデル推定結果出力装置１０９は、最適な隠れ状態数、観測確率のパラメータ、変分分布などをモデル推定結果出力結果１１２として出力する。
【００６８】
隠れ状態数設定部１０２、初期化処理部１０３、隠れ変数変分確率計算処理部１０４、情報量基準近似処理部１０５、隠れ状態選択処理部１０６、パラメータ最適化処理部１０７および最適性判定処理部１０８は、例えば、隠れ属性モデル推定プログラムに従って動作するコンピュータのＣＰＵによって実現される。この場合、ＣＰＵが隠れ属性モデル推定プログラムを読み込み、そのプログラムに従って、隠れ状態数設定部１０２、初期化処理部１０３、隠れ変数変分確率計算処理部１０４、情報量基準近似処理部１０５、隠れ状態選択処理部１０６、パラメータ最適化処理部１０７および最適性判定処理部１０８として動作すればよい。隠れ属性モデル推定プログラムは、コンピュータが読み取り可能な記録媒体に記憶されていてもよい。また、上記の各要素１０２〜１０８がそれぞれ別々のハードウェアで実現されていてもよい。
【００６９】
図２は、本発明の処理経過の例を示すフローチャートである。データ入力装置１０１を介して入力データ１１１が入力される（ステップＳ１００）。
【００７０】
次に、隠れ状態数設定部１０２は、隠れ状態数の初期値として、入力データ１１１として入力された隠れ状態数の最大値を設定する（ステップＳ１０１）。すなわち、隠れ状態数設定部１０２は、モデルの隠れ状態数Kを、入力された最大値Kmaxに設定する。
【００７１】
次に、初期化処理部１０３は、指定された隠れ状態数に対して、推定のためパラメータ（例えば、各観測確率のパラメータj）や隠れ変数変分確率の初期化処理を実施する（ステップＳ１０２）。
【００７２】
次に、情報量基準近似処理部１０５は、ヘッセ行列の行列式の近似処理を行う（ステップＳ１０３）。情報量基準近似処理部１０５は、式１１の計算によりヘッセ行列の行列式の近似値を計算する。
【００７３】
次に、隠れ変数変分確率計算処理部１０４は、計算されたヘッセ行列の行列式の近似値を利用して、隠れ変数の変分確率を計算する（ステップＳ１０４）。
【００７４】
次に、隠れ状態選択処理部１０６は、前述の閾値判定により、不要な隠れ状態をモデルから除去する（ステップＳ１０５）。すなわち、隠れ状態選択処理部１０６は、k番目の隠れ状態について、
【数３２】
0006398991.tif 000033
が、入力データ１１１として設定された閾値を下回った場合には、その状態をモデルから削除する。
【００７５】
次に、パラメータ最適化処理部１０７は、最適化基準Aを最適とするパラメータを算出する（ステップＳ１０６）。なお、パラメータ最適化処理部１０７が最初にステップＳ１０６を実行する場合に用いる最適化基準Aは、例えば、初期化処理部１０３によってランダムに定められていてもよい。あるいは、初期化処理部１０３が隠れ変数の変分確率をランダムに設定し、１回目のステップＳ１０３〜Ｓ１０９ａのループ処理（図２参照）では、ステップＳ１０６の処理を実行しなくてもよい。
【００７６】
次に、情報量基準近似処理部１０５は、ヘッセ行列の行列式の近似処理を行う（ステップＳ１０７）。情報量基準近似処理部１０５は、式１１の計算によりヘッセ行列の行列式の近似値を計算する。
【００７７】
次に、パラメータ最適化処理部１０７は、ステップＳ１０６で最適化されたパラメータを利用し、最適化基準Aの値を算出する（ステップＳ１０８）。
【００７８】
次に、最適性判定処理部１０８は、最適化基準A が収束したか否かを判定する（ステップＳ１０９）。例えば、最適性判定処理部１０８は、直近のステップＳ１０３〜Ｓ１０９ａのループ処理で得られた最適化基準Aと、その１つ前のステップＳ１０３〜Ｓ１０９ａのループ処理で得られた最適化基準Aとの差を計算し、その差の絶対値が予め定められた閾値以下になっていれば、最適化基準Aが収束したと判定してよい。また、最適性判定処理部１０８は、その差の絶対値が閾値より大きければ、最適化基準Aが収束していないと判定してよい。
【００７９】
最適化基準Aが収束していないと判定した場合（ステップＳ１０９ａのＮｏ）、隠れ属性モデル推定装置１００は、ステップＳ１０３以降の処理を繰り返す。最適化基準Aが収束したと判定した場合（ステップＳ１０９ａのＹｅｓ）、モデル推定結果出力装置１０９は、モデル推定結果を出力し、処理を完了する（ステップＳ１１０）。ステップＳ１１０では、モデル推定結果出力装置１０９は、最適化基準Aが収束したと判定された時点における隠れ状態数、およびその時点で得られているパラメータ、変分分布を出力する。
【００８０】
以下では、本発明で提案する隠れ属性モデル推定装置の応用例について、健康診断データの因子分析を例に説明する。本例では、Xとして、健康診断の被験者を行方向（サンプル）とし、血圧や血糖値、BMIなど、健康診断の検査項目の値を列方向（属性）とする行列を考える。各検査項目値の分布は、例えば年齢や性別のような観測が容易な因子だけではなく、生活習慣のような観測が難しい因子が複雑に絡まって構成される。また、因子の数を事前に決定することは難しい。そして、恣意的な分析を排除するために、因子の数をデータから自動的に決定できることが望ましい。
【００８１】
このようなデータに対して、本発明で提案する隠れ属性モデル推定装置を適用すると、各項目の多変量の依存性を考慮し、各サンプルに対する隠れ属性の変分分布を推定することが可能である。例えば、あるサンプルに対する因子を分析するために、そのサンプルの変分分布の期待値が０．５より大きい因子は「影響あり」とし、期待値が０．５より小さい因子は「影響なし」として、影響度の高い因子を分析することができる。さらに、本発明によれば、因子化漸近ベイズ推論の枠組みにしたがって、周辺化尤度最大化という意味で、隠れ属性の数を適切に決定することが可能である。例えば、主成分分析による因子の分析では、観測変数の最も特徴的な変数を因子として扱うが、本発明に従うと、観測されない因子をデータから自動的に見つけることができるという大きな効果がある。
【００８２】
次に、本発明の概要について説明する。図３は、本発明の概要を示すブロック図である。本発明の隠れ属性モデル推定装置１００は、近似値計算部７１と、変分確率計算部７２と、隠れ状態削除部７３と、パラメータ最適化部７４と、収束判定部７５とを備える。
【００８３】
近似値計算部７１（例えば、情報量基準近似処理部１０５）は、行列として表される観測データに関するヘッセ行列の行列式の近似値を計算する（例えば、式１１の近似計算を行う）。
【００８４】
変分確率計算部７２（例えば、隠れ変数変分確率計算処理部１０４）は、その行列式の近似値を利用して隠れ変数の変分確率を計算する。
【００８５】
隠れ状態削除部７３（例えば、隠れ状態選択処理部１０６）は、変分分布に基づいて隠れ状態を削除する。
【００８６】
パラメータ最適化部７４（例えば、パラメータ最適化処理部１０７）は、周辺化対数尤度関数を完全変数に対する推定量に関してラプラス近似した近似量の下界として定義される基準値（例えば、最適化基準A）に対してパラメータを最適化するとともに、基準値を計算する。
【００８７】
収束判定部７５（例えば、最適性判定処理部１０８）は、基準値が収束したか否かを判定する。
【００８８】
また、近似値計算部７１がヘッセ行列の行列式の近似値を計算し、変分確率計算部７２が隠れ変数の変分確率を計算し、隠れ状態削除部７３が隠れ状態を削除し、パラメータ最適化部７４がパラメータを最適化し、近似値計算部７１がヘッセ行列の行列式の近似値を計算し、パラメータ最適化部７４が基準値を計算し、収束判定部７５が、基準値が収束したか否かを判定するループ処理を、収束判定部７５が、基準値が収束したと判定するまで繰り返す構成であることが好ましい。
【００８９】
この出願は、２０１３年５月２０日に出願された米国特許出願１３／８９８１１８を基礎とする優先権を主張し、その開示の全てをここに取り込む。
【００９０】
以上、実施形態を参照して本願発明を説明したが、本願発明は上記の実施形態に限定されるものではない。本願発明の構成や詳細には、本願発明のスコープ内で当業者が理解し得る様々な変更をすることができる。
【符号の説明】
【００９１】
１０１  データ入力装置
１０２  隠れ状態数設定部
１０３  初期化処理部
１０４  隠れ変数変分確率計算処理部
１０５  情報量基準近似処理部
１０６  隠れ状態選択処理部
１０７  パラメータ最適化処理部
１０８  最適性判定処理部
１０９  モデル推定結果出力装置
【図１】
0006398991.tif 000034
【図２】
0006398991.tif 000035
【図３】
0006398991.tif 000036