JP6935765B2 - 動的分布推定装置、方法、及びプログラム - Google Patents

動的分布推定装置、方法、及びプログラム Download PDF

Info

Publication number
JP6935765B2
JP6935765B2 JP2018023593A JP2018023593A JP6935765B2 JP 6935765 B2 JP6935765 B2 JP 6935765B2 JP 2018023593 A JP2018023593 A JP 2018023593A JP 2018023593 A JP2018023593 A JP 2018023593A JP 6935765 B2 JP6935765 B2 JP 6935765B2
Authority
JP
Japan
Prior art keywords
component
data
statistic
sample
update
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Active
Application number
JP2018023593A
Other languages
English (en)
Other versions
JP2019139597A (ja
Inventor
匡宏 幸島
匡宏 幸島
寛 清武
寛 清武
達史 松林
達史 松林
塩原　寿子
寿子 塩原
浩之 戸田
浩之 戸田
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Nippon Telegraph and Telephone Corp
Original Assignee
Nippon Telegraph and Telephone Corp
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Nippon Telegraph and Telephone Corp filed Critical Nippon Telegraph and Telephone Corp
Priority to JP2018023593A priority Critical patent/JP6935765B2/ja
Priority to PCT/JP2019/004677 priority patent/WO2019159845A1/ja
Priority to US16/969,052 priority patent/US20210035000A1/en
Publication of JP2019139597A publication Critical patent/JP2019139597A/ja
Application granted granted Critical
Publication of JP6935765B2 publication Critical patent/JP6935765B2/ja
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06N—COMPUTING ARRANGEMENTS BASED ON SPECIFIC COMPUTATIONAL MODELS
    • G06N7/00—Computing arrangements based on specific mathematical models
    • G06N7/01—Probabilistic graphical models, e.g. probabilistic networks
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
    • G06F17/10—Complex mathematical operations
    • G06F17/18—Complex mathematical operations for evaluating statistical data, e.g. average values, frequency distributions, probability functions, regression analysis

Landscapes

  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Theoretical Computer Science (AREA)
  • Data Mining & Analysis (AREA)
  • Mathematical Analysis (AREA)
  • Computational Mathematics (AREA)
  • Mathematical Physics (AREA)
  • Mathematical Optimization (AREA)
  • Pure & Applied Mathematics (AREA)
  • Software Systems (AREA)
  • Algebra (AREA)
  • Probability & Statistics with Applications (AREA)
  • General Engineering & Computer Science (AREA)
  • Computing Systems (AREA)
  • Evolutionary Computation (AREA)
  • Artificial Intelligence (AREA)
  • Life Sciences & Earth Sciences (AREA)
  • Bioinformatics & Cheminformatics (AREA)
  • Bioinformatics & Computational Biology (AREA)
  • Evolutionary Biology (AREA)
  • Operations Research (AREA)
  • Databases & Information Systems (AREA)
  • Complex Calculations (AREA)

Description

本発明は、動的分布推定装置、方法、及びプログラムに関する。
打ち切りデータとは観測値がある閾値以上（またはある値以下）であるサンプルについては、値が観測されず、閾値以上である、という情報しか得られないデータのことを指す。病気の発症や人の死亡などを記述する臨床データや、インターネット回線利用者の契約履歴データ、Ｅコマースサイトのサービス利用履歴データなど多くのデータが打ち切りデータとして表現される。上記の例と同様に、有名アーティストの音楽ライブや人気スポーツの国際試合などのイベントの当日に収集される観客のイベント周辺への到着時間に関するデータも打ち切りデータとして表現される。図７に具体例を示す。チケット総販売数から分かる、総来場者予定者数をＮ人と書き、ライブ当日の現在時刻までに観測された来場者数をＭ人と書く。到着済みのＭ人に関しては到着時間のデータが得られているが、残りのＮ−Ｍ人については、現在時刻までには到着していない、ということしか分からない。これは典型的な打ち切りデータである。
打ち切りデータから混合モデルのパラメタを（バッチ的に）推定するという技術は、非特許文献２及び非特許文献３で提案されている。ここでは一例として、代表的な混合モデルの１つである混合正規分布の既存技術について述べる。
＜モデル＞
入力データが右側打ち切りされている状況を考える。右側打切りとは、サンプルの中で値がある既知の閾値
Figure 0006935765
以上となるサンプルについては値が分からない、という状況のことを指す。得られた全データを
Figure 0006935765
と書く。ｄｉがｉ番目データを表し、ｄｉ＝（ｗｉ，Ｘｉ）とｉ番目サンプルの値が観測されたか否かを表す変数ｗｉ∈｛０，１｝と観測された値
Figure 0006935765
の２つからなる。ｗｉ＝１が値を観測されたこと、ｗｉ＝０が値が観測されなかったことを表す。(値の観測されないものも含めた)全サンプル数をＮ、そのうち値の観測されたサンプルの個数を
Figure 0006935765
と書く。本研究で考える設定では閾値Ｃは既知であり、Ｘ，Ｗの２つが観測変数である。 一般に混合モデルの確率密度関数は次の式で定義される。
Figure 0006935765
Ｋはコンポーネント数、
Figure 0006935765
がモデルのパラメタを表す。
Figure 0006935765
はそれぞれｋ番目のコンポーネントの混合比とコンポーネントのパラメタを表す。本稿では特にコンポーネントとして正規分布を採用した場合を考える（以下の議論は指数分布など任意の指数型分布族に属する分布の混合モデルを考える場合でも同様に成り立つ。）。正規分布の確率密度関数は平均μｋと標準偏差σｋの2種類コンポーネントのパラメタを用いて、次の式で与えられる。
Figure 0006935765
また、以後正規分布の累積密度関数を関数Ｆで表す。
Figure 0006935765
打ち切りデータの生成過程は次の４ステップから成る。まず初めに、各データｉについて、ｉ番目データが所属するコンポーネントを表す潜在変数
Figure 0006935765
が、下記の多項分布に従い生成される。なお、ｉ番目のデータが第ｋ番目コンポーネントに属するならばｚｉｋ＝１、それ以外のｋ’≠ｋについてはｚｉｋ＝０である。
Figure 0006935765
次に、値が観測されるか否かを表す観測変数ｗｉが下記の所属コンポーネントの累積密度関数をパラメタに持つベルヌーイ分布に従って生成される。
Figure 0006935765
なお、累積密度
Figure 0006935765
は、確率変数が閾値Ｃ以下となる確率を表す。
さらに、ｗｉ＝１、すなわち観測可能となったデータｉは、観測変数
Figure 0006935765
が切断正規分布に従い生成される。
Figure 0006935765
なお、切断正規分布
Figure 0006935765
は範囲［ａ，ｂ］以外には値のとらない以下の確率密度関数で定義される。
Figure 0006935765
最後に、ｗｉ＝０、すなわち観測不可能となったデータｉは、潜在変数ｙｉが切断正規分布に従い生成される。
Figure 0006935765
以上を全てのデータｉに関して繰り返すことで、観測変数Ｘ，Ｗと潜在変数Ｚ，Ｙが生成される。
以後表記の簡便さのため、生成されたデータは
Figure 0006935765
でｗｉ＝１、
Figure 0006935765
ではｗｊ＝０となるように並び替えてあるとする。このとき、式（４）（５）（６）（８）を用いて完全データの尤度関数は次の式で与えられる。
Figure 0006935765
＜バッチ型ＥＭアルゴリズム＞
Expectation-Maximization（ＥＭ）アルゴリズムは、潜在変数を含むモデルの推定に広く利用される手法である。潜在変数の事後確率の算出とそれを用いた期待値の計算からなるＥステップと、Ｑ関数と呼ばれる、対数尤度関数を潜在変数の事後確率に関して平均した関数を最大化するＭステップの２ステップからなる。
本モデルのＥステップにおいては、観測値が得られた場合の事後確率Ｐ（ｚｉ｜ｘｉ，ｗｉ＝１，θ）と得られなかった場合のＰ（ｚｉ｜ｘｉ，ｗｉ＝０，θ）の２つが必要となり、これらはそれぞれ以下の式で与えられる。
Figure 0006935765
上記の事後確率を用いて、下記の式でｚｉ，ｚｊの負担率γ，ηとｙｊのモーメント｛νｋ，ξｋ｝を計算できる。
Figure 0006935765
ただし、
Figure 0006935765
は事後確率
Figure 0006935765
の出方に関する平均を表す。この平均操作に関しては切断正規分布の１次モーメントと２次モーメントの結果を利用している。また、上記式（１２）（１３）から明らかなように
Figure 0006935765
は添え字ｊに依存しないため以後
Figure 0006935765
と書く。これらを用いるとＭステップで最大化するＱ関数は以下の式で表現される。
Figure 0006935765
ただし、
Figure 0006935765
偏微分をゼロと置いて解くとＱ関数を最大化するパラメタは
Figure 0006935765
で与えられる。これにより打ち切りデータに対する混合モデルのバッチ型ＥＭアルゴリズムが求められた。以下に示すAlgorithm1に手続きをまとめる。Ｅステップ、Ｍステップによってパラメタの更新を繰り返し、各反復において、対数尤度関数は単調増加し、（局所）最適解への収束が保証される。
Figure 0006935765
Didier Chauveau. , "A stochastic em algorithm for mixtures with censored data." , Journal of statistical planning and inference, 46(1):p.1-25, 1995. Gyemin Lee and Clayton Scott. "Em algorithms for multivar-iate gaussian mixture models with truncated and censored data.", ComputationalStatistics & Data Analysis, 56(9):p.2816-2829, 2012.
既存技術は、打ち切りデータに対してバッチ型の推定を行うことしかできなかった。
本発明は、上記の事情を鑑みてなされたものであり、高速、かつ省メモリ、かつパラメタが時間連続性を有する状態で、打ち切りデータを含むモデルのパラメタを推定することができる動的分布推定装置、方法、及びプログラムを提供することを目的とする。
上記の目的を達成するために本発明に係る動的分布推定装置は、観測されるデータの分布を表す、指数型分布族に属する任意の分布を混合した、混合モデルのパラメタをオンラインで推定する動的分布推定装置であって、新たに観測されたサンプルのデータに基づいて、前記観測されていない各サンプルのデータが各コンポーネントに所属すると仮定した場合の、前記観測されていないサンプルのデータの十分統計量の、切断されたコンポーネントの分布による期待値を更新する期待値更新部と、前記新たに観測されたサンプルのデータと、前記期待値更新部によって更新された前記期待値とに基づいて、各コンポーネントに関する統計量を更新する統計量更新部と、前記統計量更新部によって更新された前記統計量に基づいて、各コンポーネントについて、前記コンポーネントに関するパラメタを更新するパラメタ更新部と、を含み、予め定められたパラメタ更新タイミングが到来する毎に、前記期待値更新部による更新、前記統計量更新部による更新、及び前記パラメタ更新部による更新を繰り返す。
また、本発明に係る動的分布推定装置は、観測されるデータの分布を表す、複数のコンポーネントを混合した混合ガウスモデルのパラメタをオンラインで推定する動的分布推定装置であって、新たに観測されたサンプルのデータに基づいて、前記新たに観測されたサンプルのデータが各コンポーネントに所属する度合いを表す負担率、及びまだ観測されていない各サンプルのデータが各コンポーネントに所属する度合いを表す負担率を更新する負担率更新部と、前記観測されていない各サンプルのデータが各コンポーネントに所属すると仮定した場合の、前記観測されていないサンプルのデータのモーメントを更新するモーメント更新部と、前記新たに観測されたサンプルのデータが各コンポーネントに所属する度合いを表す負担率に基づいて、観測されたサンプルのうち、各コンポーネントに所属するサンプル数の統計量を更新し、前記観測されていない各サンプルのデータが各コンポーネントに所属する度合いが表す負担率に基づいて、全サンプルのうち、各コンポーネントに所属するサンプル数の統計量を更新し、前記新たに観測されたサンプルのデータが各コンポーネントに所属する度合いが表す負担率に基づいて、各コンポーネントについて、前記コンポーネントに所属する、観測されたサンプルのデータの統計量を更新し、各コンポーネントについて、前記新たに観測された各サンプルのデータが前記コンポーネントに所属すると仮定した場合の、前記観測されていないサンプルのデータのモーメント、観測されたサンプルのうち、前記コンポーネントに所属するサンプル数の統計量、及び全サンプルのうち、前記コンポーネントに所属するサンプル数の統計量に基づいて、前記コンポーネントに所属する、前記観測されていないサンプルのデータの統計量を更新する統計量更新部と、各コンポーネントについて、全サンプルのうち、前記コンポーネントに所属するサンプル数の統計量、前記コンポーネントに所属する、前記観測されたサンプルのデータの統計量、及び前記コンポーネントに所属する、前記観測されていないサンプルのデータの統計量に基づいて、前記コンポーネントに関するパラメタを更新するパラメタ更新部と、を含み、予め定められたパラメタ更新タイミングが到来する毎に、前記負担率更新部による更新、前記モーメント更新部による更新、前記統計量更新部による更新、及び前記パラメタ更新部による更新を繰り返す。
本発明に係る動的分布推定装置は、観測されるデータの分布を表す、複数のコンポーネントを混合した混合ガウスモデルのパラメタをオンラインで推定する動的分布推定装置であって、新たに観測されたサンプルのデータに基づいて、前記新たに観測されたサンプルのデータについての、各コンポーネントの潜在変数に関する変分分布のパラメタ、及び新たに観測されたサンプルを含む、既に観測されたサンプル集合のデータについての、各コンポーネントの潜在変数に関する変分分布のパラメタを更新する潜在変数パラメタ更新部と、前記新たに観測されたサンプルのデータについての、各コンポーネントに対する変分分布のパラメタに基づいて、観測されたサンプルのうち、各コンポーネントに所属するサンプル数の統計量を更新し、既に観測されたサンプル集合のデータについての、各コンポーネントに対する変分分布のパラメタに基づいて、まだ観測されていないサンプルのうち、各コンポーネントに所属するサンプル数の統計量を更新し、前記新たに観測されたサンプルのデータについての、各コンポーネントに対する変分分布のパラメタに基づいて、各コンポーネントについて、前記コンポーネントに所属する、観測されたサンプルのデータの統計量を更新し、前記既に観測されたサンプル集合のデータと、前記観測されていないサンプルのうち、各コンポーネントに所属するサンプル数の統計量とに基づいて、前記コンポーネントに所属する、前記観測されていないサンプルのデータの統計量を更新する統計量更新部と、各コンポーネントについて、全サンプルのうち、前記コンポーネントに所属するサンプル数、前記コンポーネントに所属する、前記観測されたサンプルのデータの統計量、及び前記コンポーネントに所属する、前記観測されていないサンプルのデータの統計量に基づいて、前記コンポーネントのパラメタに関する変分分布のパラメタを更新するパラメタ更新部と、を含み、予め定められたパラメタ更新タイミングが到来する毎に、前記潜在変数パラメタ更新部による更新、前記統計量更新部による更新、及び前記パラメタ更新部による更新を繰り返す。
本発明の前記パラメタ更新タイミングは、前記新たに観測されたサンプルのデータが得られたタイミング、前記新たに観測されたサンプルのデータが予め定められた個数だけ得られたタイミング、及び予め定められた更新時期が到来したタイミングの何れかであるようにすることができる。
本発明の動的分布推定方法は、観測されるデータの分布を表す、指数型分布族に属する任意の分布を混合した、混合モデルのパラメタをオンラインで推定する動的分布推定装置であって、期待値更新部が、新たに観測されたサンプルのデータに基づいて、前記観測されていない各サンプルのデータが各コンポーネントに所属すると仮定した場合の、前記観測されていないサンプルのデータの十分統計量の、切断されたコンポーネントの分布による期待値を更新するステップと、統計量更新部が、前記新たに観測されたサンプルのデータと、前記期待値更新部によって更新された前記期待値とに基づいて、各コンポーネントに関する統計量を更新するステップと、パラメタ更新部が、前記統計量更新部によって更新された前記統計量に基づいて、各コンポーネントについて、前記コンポーネントに関するパラメタを更新するステップと、を含み、予め定められたパラメタ更新タイミングが到来する毎に、前記期待値更新部による更新、前記統計量更新部による更新、及び前記パラメタ更新部による更新を繰り返す。
本発明の動的分布推定方法は、観測されるデータの分布を表す、複数のコンポーネントを混合した混合ガウスモデルのパラメタをオンラインで推定する動的分布推定装置における動的分布推定方法であって、負担率更新部が、新たに観測されたサンプルのデータに基づいて、前記新たに観測されたサンプルのデータが各コンポーネントに所属する度合いを表す負担率、及びまだ観測されていない各サンプルのデータが各コンポーネントに所属する度合いを表す負担率を更新するステップと、モーメント更新部が、前記観測されていない各サンプルのデータが各コンポーネントに所属すると仮定した場合の、前記観測されていないサンプルのデータのモーメントを更新するステップと、統計量更新部が、前記新たに観測されたサンプルのデータが各コンポーネントに所属する度合いを表す負担率に基づいて、観測されたサンプルのうち、各コンポーネントに所属するサンプル数の統計量を更新し、前記観測されていない各サンプルのデータが各コンポーネントに所属する度合いが表す負担率に基づいて、全サンプルのうち、各コンポーネントに所属するサンプル数の統計量を更新し、前記新たに観測されたサンプルのデータが各コンポーネントに所属する度合いが表す負担率に基づいて、各コンポーネントについて、前記コンポーネントに所属する、観測されたサンプルのデータの統計量を更新し、各コンポーネントについて、前記新たに観測された各サンプルのデータが前記コンポーネントに所属すると仮定した場合の、前記観測されていないサンプルのデータのモーメント、観測されたサンプルのうち、前記コンポーネントに所属するサンプル数の統計量、及び全サンプルのうち、前記コンポーネントに所属するサンプル数の統計量に基づいて、前記コンポーネントに所属する、前記観測されていないサンプルのデータの統計量を更新するステップと、パラメタ更新部が、各コンポーネントについて、全サンプルのうち、前記コンポーネントに所属するサンプル数の統計量、前記コンポーネントに所属する、前記観測されたサンプルのデータの統計量、及び前記コンポーネントに所属する、前記観測されていないサンプルのデータの統計量に基づいて、前記コンポーネントに関するパラメタを更新するステップと、を含み、予め定められたパラメタ更新タイミングが到来する毎に、前記負担率更新部による更新、前記モーメント更新部による更新、前記統計量更新部による更新、及び前記パラメタ更新部による更新を繰り返す。
本発明の動的分布推定方法は、観測されるデータの分布を表す、複数のコンポーネントを混合した混合ガウスモデルのパラメタをオンラインで推定する動的分布推定装置における動的分布推定方法であって、潜在変数パラメタ更新部が、新たに観測されたサンプルのデータに基づいて、前記新たに観測されたサンプルのデータについての、各コンポーネントの潜在変数に関する変分分布のパラメタ、及び新たに観測されたサンプルを含む、既に観測されたサンプル集合のデータについての、各コンポーネントの潜在変数に関する変分分布のパラメタを更新するステップと、統計量更新部が、前記新たに観測されたサンプルのデータについての、各コンポーネントに対する変分分布のパラメタに基づいて、観測されたサンプルのうち、各コンポーネントに所属するサンプル数の統計量を更新し、既に観測されたサンプル集合のデータについての、各コンポーネントに対する変分分布のパラメタに基づいて、まだ観測されていないサンプルのうち、各コンポーネントに所属するサンプル数の統計量を更新し、前記新たに観測されたサンプルのデータについての、各コンポーネントに対する変分分布のパラメタに基づいて、各コンポーネントについて、前記コンポーネントに所属する、観測されたサンプルのデータの統計量を更新し、前記既に観測されたサンプル集合のデータと、前記観測されていないサンプルのうち、各コンポーネントに所属するサンプル数の統計量とに基づいて、前記コンポーネントに所属する、前記観測されていないサンプルのデータの統計量を更新するステップと、パラメタ更新部が、各コンポーネントについて、全サンプルのうち、前記コンポーネントに所属するサンプル数、前記コンポーネントに所属する、前記観測されたサンプルのデータの統計量、及び前記コンポーネントに所属する、前記観測されていないサンプルのデータの統計量に基づいて、前記コンポーネントのパラメタに関する変分分布のパラメタを更新するステップと、を含み、予め定められたパラメタ更新タイミングが到来する毎に、前記潜在変数パラメタ更新部による更新、前記統計量更新部による更新、及び前記パラメタ更新部による更新を繰り返す。
本発明に係るプログラムは、本発明の動的分布推定装置の各部として機能させるためのプログラムである。
以上説明したように、本発明の動的分布推定装置、方法、及びプログラムによれば、複数のコンポーネントを混合した、指数型分布族に属する任意の分布のパラメタをオンラインで推定することにより、高速、かつ省メモリ、かつパラメタが時間連続性を有する状態で、打ち切りデータを含むモデルのパラメタを推定することができる、という効果が得られる。
逐次更新型オンラインアルゴリズムを説明するための説明図である。 更新のタイミングを説明するための説明図である。 第１の実施の形態に係る動的分布推定装置の構成例を示す概略図である。 第１の実施の形態に係る動的分布推定装置における動的分布推定処理ルーチンの内容を示すフローチャートである。 第２の実施の形態に係る動的分布推定装置の構成例を示す概略図である。 第２の実施の形態に係る動的分布推定装置における動的分布推定処理ルーチンの内容を示すフローチャートである。 打ち切りデータを説明するための説明図である。
以下、図面を参照して本発明の実施の形態を詳細に説明する。
＜本発明の実施の形態の概要＞
イベント当日にデータを収集している状況においては、時間経過につれ、新たに到着した観客のデータが観測でき、データが時事刻々と更新されていく。このような状況で到着時間分布のモデルパラメタを推定するうえでは、新たに到着したデータを反映して逐次パラメタを更新する、オンラインアルゴリズム（例えば、参考文献（Olivier Cappe and Eric Moulines. "On-line expectation-maximization algorithm for latent data models.", Journal of the Royal Statistical Society: Series B (Statistical Methodology), 71(3)p.593-613, 2009.）を参照。）が有用である。
そこで、本発明の実施形態では、時事刻々と更新されていく打ち切りデータから到着時間のモデルパラメタをオンラインに推定するアルゴリズムであるオンラインＥＭＣＭアルゴリズム(online EM algorithm for Censored Mixture models)を構築した。この技術は、次の３つの点でバッチ型の手法に対する優位性を持つ。
（１）省メモリである点。
本実施形態の提案アルゴリズムは、十分統計量のみの保持でパラメタ更新が可能であり、到着済み観客の到着時間全てを保持する必要がない。これはプライバシー保護の観点からも優れている。
（２）高速である点。
本実施形態の提案アルゴリズムは、前述した統計量と新たに観測されたデータから計算される量を用いて更新される。全てのデータを利用して計算するバッチ型手法と比べて短い処理でパラメタ更新処理を行うことができる。特に前述のようなライブやスポーツイベントでは来場者数は数万人規模に及び、各時刻での全データを利用したバッチ処理は避けられることが好ましい。
（３）推定パラメタの時間連続性を有する点。
本実施形態の提案アルゴリズムによって各時刻毎に出力されるパラメタは、前時刻におけるパラメタから連続的に変化した値となる。各時刻でバッチ型手法を適用し直す処理を行うと、目的関数の異なる局所最適解に到達することで前時刻とは全く異なるパラメタが出力される可能性があり、これは実用上好ましくないが、本実施形態にはそのような問題がない。
本実施形態では、さらに多様な実システムの実装形態に合わせられるようパラメタ更新のタイミングの異なる、（ａ）逐次更新型、（ｂ）ミニバッチ型、（ｃ）スケジュール型の３種類のアルゴリズムを示す。これら３つは全て前述の３つの優位性を持つアルゴリズムである。これによって、新データ入手時に即パラメタ更新を行う場合、いくつかのまとまったデータが集まってからパラメタ更新を行う場合、決まったタイミングでパラメタ更新を行う場合といういずれの場合であっても本技術を適用できるようになる。また、上記では到着時間分布の推定の例として説明したが、本実施形態は広く打ち切りデータのパラメタ推定に利用可能である。
なお、本実施形態では、到着時間のモデルには混合モデルを採用した。なぜなら、前述のようなイベントにおいては、イベント開始前にアーティストグッズやユニフォームなどの物販購入をするか、イベント開始ちょうどに間に合うようにするか、などに応じて観客の到着時間分布は多峰性を持つことが想像されるからである。
なお、本実施形態は、指数分布や対数正規分布など正規分布以外の他の確率分布の混合モデルを考える場合であっても、ほぼ同様に適用することができる。
上記Algorithm1のバッチ型ＥＭアルゴリズムは、Ｅステップでメモリの全データに対して負担率を計算し、それらを用いてＭステップで統計量
Figure 0006935765
を計算することを繰り返している。これはすなわち、データＸの値全てをメモリに保持し、各反復でこのメモリ全体を読みにいくことを必要としていることになる。それに対して我々の提案するアルゴリズムであるオンラインＥＭＣＭアルゴリズム(online Expectation-Maximization algorithm for Censored Mixture models)は、データ全てをメモリに保持する必要がなく、新たに観測されたデータのみを利用して、負担率や統計量を計算してパラメタ更新を行う。
＜逐次更新型オンラインＥＭアルゴリズム＞
まず、新たにデータｘｔが観測されるたびにパラメタを更新する、逐次更新型のアルゴリズムについて説明する。このアルゴリズムでは、図１に示すようにデータ観測と更新のタイミングが一致する。提案アルゴリズムはバッチ型のアルゴリズムの統計量が逐次的な形で書けることを利用して導出する。データｘｔ−１が観測された時点での統計量を上付き添え字（ｔ−１）で表すと、具体的な逐次計算式は以下で与えられる。
Figure 0006935765
まず、Ｅステップで、新しい観測データｘｔを用いて計算される、負担率を計算する。このデータの観測によって、その時点で未観測なデータもｘｔ以上ということが分かるので、それに合わせてモーメントの計算を行う。その後にＭステップで上記の逐次式を用いて統計量を計算し、パラメタを更新する。これらを新しいデータが到着するたびに行うことでパラメタを推定していく。提案アルゴリズムの手続きをAlgorithm2にまとめる。
Figure 0006935765
＜ミニバッチ型とスケジュール型のアルゴリズム＞
前節では、データが更新するたびのパラメタ更新を行っていた。しかし、毎回のデータは必ずしも必須でなく、多様な実システムの実装形態に合わせられるようパラメタ更新のタイミングの異なるアルゴリズムを導出できる。したがってこの節では、前節の（ａ）逐次更新型に加えて、図２に示す（ｂ）ミニバッチ型、（ｃ）スケジュール型の２種類のアルゴリズムを示す。
まず、ミニバッチ型について説明する。この方法では、あらかじめパラメタ更新までに蓄えるデータの数Ｂ（これをミニバッチサイズと称する。）を定めておき、この数のデータが蓄えられた時点でパラメタ更新を行う。Ｅステップの計算は、Ｍステップは下記のようにミニバッチから計算される負担率、モーメントのデータを求めて下記のように統計量を更新する。逐次更新型に比べてＭステップの実行回数が減るため、処理時間をより少なくすることができる。
Figure 0006935765
提案アルゴリズムの手続きをAlgorithm3にまとめる。なお、記号
Figure 0006935765
は入力の値を越えない整数を返す床関数を表す。次に（ｃ）スケジュール型について説明する。アルゴリズムは（ｂ）ミニバッチ型とほぼ同様である。
Figure 0006935765
ただし、統計量
Figure 0006935765
が異なっている。
Figure 0006935765
提案アルゴリズムの手続きをAlgorithm4にまとめる。なお、上記３種類のアルゴリズムの更新方法をミックスさせた方法、たとえばミニバッチと更新スケジュールの両方を利用する方法も同様に構築可能であるが、割愛する。
Figure 0006935765
＜逐次更新型オンライン変分ベイズアルゴリズム＞
上記までの記述ではバッチ型のＥＭアルゴリズムを発展させたオンライン型のＥＭアルゴリズムを示した。混合モデルの推定アルゴリズムにはＥＭアルゴリズムのほかにも変分ベイズ（Variational Bayes，ＶＢ）アルゴリズムと呼ばれるアルゴリズムも存在し、本発明と同様のアプローチによって、ＶＢアルゴリズムをオンラインアルゴリズムとすることもできる。したがって本発明の範囲はオンラインＥＭアルゴリズムに限定されず、打ち切りデータに対する混合モデルのオンライン推定アルゴリズム全般を含む。以下にＶＢアルゴリズムのバッチ型アルゴリズムを基に逐次更新型のアルゴリズムを導く例を以下に記す。
＜バッチ型変分ベイズ（ＶＢ）アルゴリズム＞
変分ベイズアルゴリズムではモデルのパラメタ
Figure 0006935765
に事前分布
Figure 0006935765
が設定されていることを考える。ただし、
Figure 0006935765
は精度パラメタであり、この章では標準偏差
Figure 0006935765
の代わりに精度
Figure 0006935765
を用いて正規分布の確率密度関数が
Figure 0006935765
と表現されているとする。
Figure 0006935765
はそれぞれ（対称）ディリクレ分布と正規-ガンマ分布であり、以下の式で定義される。
Figure 0006935765
上記の式と上記式（９）とを組み合わせると、パラメタと完全データの生成確率は以下の式で表せる。
Figure 0006935765
（バッチ型の）ＶＢアルゴリズムは、パラメタと潜在変数の事後分布を近似する、変分分布を推定する方法である。打ち切りデータに対するＶＢアルゴリズムでは、変分分布が
Figure 0006935765
と分解されるという条件のもと汎関数
Figure 0006935765
を最小化することで変分分布を推定することを考える。
Figure 0006935765
変分法による解析から所望の変分分布は次の最適性条件を満たさなければならないことが示される。
Figure 0006935765
上記を計算すると、
Figure 0006935765
の(最適)変分分布はそれぞれ以下のディリクレ分布、正規-ガンマ分布、多項-切断正規分布で与えられることが示される。
Figure 0006935765
上記の式中の統計量等は次の通りである。
Figure 0006935765
ただし、
Figure 0006935765
はディガンマ関数を表す。これによってAlgorithm5に示すようにバッチ型のVBアルゴリズムを構築できる。
Figure 0006935765
＜逐次更新型ＶＢアルゴリズム＞
上記のバッチ型ＶＢアルゴリズムを基にオンラインアルゴリズムを導く。統計量
Figure 0006935765
はオンラインＥＭアルゴリズムの時と同様に逐次式の形で書くことができる。データｘｔ−１が観測された時点での統計量を上付き添え字（ｔ−１）で表すと、具体的な逐次計算式は以下で与えられる。
Figure 0006935765
したがって、Algorithm6に示すように逐次更新型のアルゴリズムを構築できる。同様にＶＢアルゴリズムに基づくミニバッチ型とスケジュール型のアルゴリズムを導出することが可能であるが割愛する。
Figure 0006935765
＜第１の実施形態の動的分布推定装置１００の構成＞
第１の実施形態の動的分布推定装置１００は、逐次更新型オンラインＥＭアルゴリズムを用いてパラメタの推定を行う。
図３に示すように、第１の実施の形態に係る動的分布推定装置１００は、ＣＰＵ（Central Processing Unit）と、ＲＡＭ（Random Access Memory）と、後述する動的分布推定ルーチンを実行するためのプログラムを記憶したＲＯＭ（Read Only Memory）とを備えたコンピュータ１０と外部装置３０とを含んで構成されている。コンピュータ１０は、機能的には、記憶部１２と、初期化処理部１７と、更新処理部１８と、パラメタ処理部２３と、入出力部２４とを備えている。
記憶部１２は、パラメタ記録部１３と、観測データ数記録部１４と、閾値記録部１５と、統計量記録部１６とを含む。
パラメタ記録部１３には、モデルのパラメタ
Figure 0006935765
が格納される。
観測データ数記録部１４には、観測されたデータの数Ｍが格納される。
閾値記録部１５には、観測されたデータの閾値Ｃが格納される。
統計量記録部１６には、各統計量
Figure 0006935765
Figure 0006935765
が格納される。
初期化処理部１７は、パラメタ記録部１３に格納された変分パラメタと、統計量記録部１６に格納された各統計量とを初期化する。
更新処理部１８は、観測されるデータの分布を表す、複数のコンポーネントを混合した混合ガウスモデルのパラメタをオンラインで推定する。更新処理部１８は、負担率更新部１９と、モーメント更新部２０と、統計量更新部２１と、パラメタ更新部２２とを備えている。モーメント更新部２０は、期待値更新部の一例である。
負担率更新部１９は、新たに観測されたサンプルのデータｘｔに基づいて、新たに観測されたサンプルのデータｘｔが各コンポーネントに所属する度合いを表す負担率
Figure 0006935765
、及びまだ観測されていない各サンプルのデータが各コンポーネントに所属する度合いを表す負担率
Figure 0006935765
を、上記式（２３）及び式（２４）に従って更新する。
モーメント更新部２０は、観測されていない各サンプルのデータが各コンポーネントに所属すると仮定した場合の、観測されていないサンプルのデータのモーメント
Figure 0006935765
を、上記式（１２）及び式（１３）に従って更新する。
統計量更新部２１は、新たに観測されたサンプルのデータｘｔが各コンポーネントに所属する度合いを表す負担率
Figure 0006935765
に基づいて、観測されたサンプルのうち、各コンポーネントに所属するサンプル数の統計量
Figure 0006935765
を、上記式（２３）に従って更新する。
統計量更新部２１は、まだ観測されていない各サンプルのデータが各コンポーネントに所属する度合いを表す負担率
Figure 0006935765
に基づいて、全サンプルのうち、各コンポーネントに所属するサンプル数の統計量
Figure 0006935765
を、上記式（２４）に従って更新する。
統計量更新部２１は、新たに観測されたサンプルのデータが各コンポーネントに所属する度合いを表す負担率
Figure 0006935765
に基づいて、各コンポーネントについて、当該コンポーネントに所属する、観測されたサンプルのデータの統計量
Figure 0006935765
Figure 0006935765
を、上記式（２５）及び上記式（２６）に従って更新する。
統計量更新部２１は、各コンポーネントについて、観測されていない各サンプルのデータが当該コンポーネントに所属すると仮定した場合の、観測されていないサンプルのデータのモーメント
Figure 0006935765
、観測されたサンプルのうち、当該コンポーネントに所属するサンプル数の統計量
Figure 0006935765
、及び全サンプルのうち、当該コンポーネントに所属するサンプル数の統計量
Figure 0006935765
に基づいて、当該コンポーネントに所属する、まだ観測されていないサンプルのデータの統計量
Figure 0006935765
を、上記式（２７）及び式（２８）に従って更新する。
パラメタ更新部２２は、各コンポーネントについて、全サンプルのうち、統計量更新部２１によって更新された、コンポーネントに所属するサンプル数の統計量
Figure 0006935765
、コンポーネントに所属する、観測されたサンプルのデータの統計量
Figure 0006935765
Figure 0006935765
、及びコンポーネントに所属する、まだ観測されていないサンプルのデータの統計量
Figure 0006935765
に基づいて、コンポーネントに関するパラメタ
Figure 0006935765
を、上記式（２０）〜（２２）に従って更新する。
入出力部２４は、パラメタ更新部２２によって更新されたパラメタπｎｅｗ ｋ，μｎｅｗ ｋ，（σｎｅｗ ｋ）２を、外部装置３０へ出力する。
外部装置３０は、入出力部２４から出力されたパラメタを結果として出力する。
＜動的分布推定装置１００の作用＞
次に、本実施の形態に係る動的分布推定装置１００の作用について説明する。まず、動的分布推定装置１００の初期化処理部１７は、パラメタ記録部１３に格納されたパラメタと、統計量記録部１６に格納された各統計量とを初期化する。そして、動的分布推定装置１００は、既に観測されたデータに基づいてＥＭアルゴリズムを用いて、モデルのパラメタと観測データ数と閾値と各統計量とを推定し、パラメタ記録部１３、観測データ数記録部１４、閾値記録部１５、及び統計量記録部１６へ格納する。そして、動的分布推定装置１００は、新たに観測されたデータｘｔが入力されると、図４に示す動的分布推定処理ルーチンを実行する。
まず、ステップＳ１００において、新たに観測されたデータｘｔを取得する。
ステップＳ１０２において、観測データ数Ｍと閾値Ｃとを更新する。
ステップＳ１０４において、負担率更新部１９は、上記ステップＳ１００で取得されたデータｘｔに基づいて、上記式（１０）に従って、負担率γ（ｚｔ）を更新する。また、負担率更新部１９は、上記ステップＳ１００で取得されたデータｘｔに基づいて、上記式（１１）に従って、負担率η（ｚｊ）を更新する。
ステップＳ１０６において、モーメント更新部２０は、上記ステップＳ１００で取得されたデータｘｔに基づいて、モーメントνｋ（ｙｊ;ｘｔ），ξｋ（ｙｊ;ｘｔ）を、上記式（１２）及び上記式（１３）に従って更新する。
ステップＳ１０８において、統計量更新部２１は、上記ステップＳ１０４で更新された負担率γ（ｚｔ）に基づいて、観測されたサンプルのうち、各コンポーネントに所属するサンプル数の統計量Ｍ（ｔ） ｋを、上記式（２３）に従って更新する。また、統計量更新部２１は、上記ステップＳ１０４で更新された負担率η（ｚｊ）に基づいて、全サンプルのうち、各コンポーネントに所属するサンプル数の統計量Ｎ（ｔ） ｋを、上記式（２４）に従って更新する。また、統計量更新部２１は、上記ステップＳ１０４で更新された負担率γ（ｚｔ）に基づいて、各コンポーネントについて、当該コンポーネントに所属する、観測されたサンプルのデータの統計量Ｓ（ｔ） ｋ１，Ｓ（ｔ） ｋ２を、上記式（２５）及び上記式（２６）に従って更新する。また、統計量更新部２１は、各コンポーネントについて、上記ステップＳ１０６で更新されたモーメントνｋ（ｙｊ;ｘｔ），ξｋ（ｙｊ;ｘｔ）、更新された統計量Ｍ（ｔ） ｋ及び統計量Ｎ（ｔ） ｋに基づいて、当該コンポーネントに所属する、まだ観測されていないサンプルのデータの統計量Ｕ（ｔ） ｋ１，Ｕ（ｔ） ｋ２を、上記式（２７）及び上記式（２８）に従って更新する。
ステップＳ１１０において、入出力部２４は、上記ステップＳ１１０で更新されたパラメタπｎｅｗ ｋ，μｎｅｗ ｋ，（σｎｅｗ ｋ）２を、外部装置３０へ出力して処理を終了する。
以上説明したように、第１の実施の形態に係る動的分布推定装置によれば、観測されるデータの分布を表す、複数のコンポーネントを混合した混合ガウスモデルのパラメタをオンラインで推定する。具体的には、第１の実施の形態に係る動的分布推定装置は、新たに観測されたサンプルのデータｘｔに基づいて負担率を更新し、観測されていないサンプルのデータのモーメントを更新し、負担率及びモーメントの少なくとも一方に基づいて各統計量を更新し、各統計量に基づいて、コンポーネントに関するパラメタを更新する。これにより、打ち切りデータに対してオンライン型のアルゴリズムを用いて、高速であり、かつ省メモリであり、かつパラメタが時間連続性を有する状態で、モデルのパラメタを推定することができる。
このように、本実施形態は、バッチ型の推定アルゴリズムと比較し、高速であって、かつ省メモリであり、かつパラメタの時間連続性を有するという３つの優位性を持つ。
なお、本発明は、上述した実施形態に限定されるものではなく、この発明の要旨を逸脱しない範囲内で様々な変形や応用が可能である。
例えば、上記の第１の実施の形態では、パラメタ更新タイミングが新たに観測されたサンプルのデータが得られたタイミングとなる、逐次更新型オンラインＥＭアルゴリズムを用いた場合を例に説明したがこれに限定されるものではない。例えば、パラメタ更新タイミングが新たに観測されたサンプルのデータが予め定められた個数だけ得られたタイミングとなる、ミニバッチ型オンラインＥＭアルゴリズムを用いてもよい。この場合には、上記Algorithm3に従って、予め定められた個数だけ新たに観測されたサンプルのデータｘｔに基づいて、負担率とモーメントとを更新し、負担率及びモーメントの少なくとも一方に基づいて各統計量を更新し、各統計量に基づいて、コンポーネントに関するパラメタを更新するようにすればよい。
また、パラメタ更新タイミングが予め定められた更新時期が到来したタイミングとなる、ミニバッチ型オンラインＥＭアルゴリズムを用いてもよい。この場合には、上記Algorithm4に従って、更新時期が到来するまでの間に新たに観測されたサンプルのデータｘｔに基づいて、負担率とモーメントとを更新し、負担率及びモーメントの少なくとも一方に基づいて各統計量を更新し、各統計量に基づいて、コンポーネントに関するパラメタを更新するようにすればよい。
＜第２の実施形態の動的分布推定装置の構成＞
第２の実施形態の動的分布推定装置は、逐次更新型オンラインＶＢアルゴリズムを用いてパラメタの推定を行う。
図５に示すように、第２の実施の形態に係る動的分布推定装置２００は、ＣＰＵ（Central Processing Unit）と、ＲＡＭ（Random Access Memory）と、後述する動的分布推定ルーチンを実行するためのプログラムを記憶したＲＯＭ（Read Only Memory）とを備えたコンピュータ２１０と外部装置３０とを含んで構成されている。コンピュータ２１０は、機能的には、記憶部２１２と、初期化処理部２１７と、更新処理部２１８と、パラメタ処理部２３と、入出力部２４とを備えている。
記憶部１２は、パラメタ記録部２１３と、観測データ数記録部１４と、閾値記録部１５と、統計量記録部２１６とを含む。
パラメタ記録部２１３には、変分パラメタ
Figure 0006935765
が格納される。
統計量記録部２１６には、各統計量
Figure 0006935765
が格納される。
初期化処理部２１７は、パラメタ記録部２１３に格納された変分パラメタと、統計量記録部１６に格納された各統計量とを初期化する。
更新処理部２１８は、観測されるデータの分布を表す、複数のコンポーネントを混合した混合ガウスモデルのパラメタをオンラインで推定する。更新処理部２１８は、潜在変数パラメタ更新部２１９と、統計量更新部２２１と、パラメタ更新部２２２とを備えている。
潜在変数パラメタ更新部２１９は、新たに観測されたサンプルのデータｘｔに基づいて、新たに観測されたサンプルのデータｘｔについての、各コンポーネントの潜在変数に関する変分分布のパラメタ
Figure 0006935765
、及び新たに観測されたサンプルを含む、既に観測されたサンプル集合のデータについての、各コンポーネントの潜在変数に関する変分分布のパラメタ
Figure 0006935765
を、上記式（５５）に従って更新する。
統計量更新部２２１は、潜在変数パラメタ更新部２１９によって更新された変分分布のパラメタ
Figure 0006935765
に基づいて、観測されたサンプルのうち、各コンポーネントに所属するサンプル数の統計量
Figure 0006935765
を、上記式（５７）に従って更新する。
また、統計量更新部２２１は、潜在変数パラメタ更新部２１９によって更新された変分分布のパラメタ
Figure 0006935765
に基づいて、まだ観測されていないサンプルのうち、各コンポーネントに所属するサンプル数の統計量
Figure 0006935765
を、上記式（５７）に従って更新する。
また、統計量更新部２２１は、潜在変数パラメタ更新部２１９によって更新された変分分布のパラメタ
Figure 0006935765
に基づいて、各コンポーネントについて、当該コンポーネントに所属する、観測されたサンプルのデータの統計量
Figure 0006935765
を、上記式（５９）及び上記式（６０）に従って更新する。
また、統計量更新部２２１は、潜在変数パラメタ更新部２１９によって更新された
Figure 0006935765
に基づいて、当該コンポーネントに所属する、まだ観測されていないサンプルのデータの統計量
Figure 0006935765
を、上記式（５９）及び上記式（６０）に従って更新する。
パラメタ更新部２２２は、各コンポーネントについて、全サンプルのうち、統計量更新部２２１によって更新されたサンプル数の統計量
Figure 0006935765
Figure 0006935765
、統計量更新部２２１によって更新された、当該コンポーネントに所属する、観測されたサンプルのデータの統計量
Figure 0006935765
、及び当該コンポーネントに所属する、まだ観測されていないサンプルのデータの統計量
Figure 0006935765
に基づいて、当該コンポーネントのパラメタに関する変分分布のパラメタ
Figure 0006935765
を、上記式（４４）〜式（５０）に従って更新する。
パラメタ処理部２２３は、予め定められたパラメタ更新タイミングが到来する毎に、潜在変数パラメタ更新部２１９による更新、統計量更新部２２１による更新、及びパラメタ更新部２２２による更新を繰り返すように、各処理部を制御する。例えば、パラメタ処理部２２３は、予め定められたパラメタ更新タイミングとして、新たに観測された各サンプルのデータが取得されたときに、潜在変数パラメタ更新部２１９による更新、統計量更新部２２１による更新、及びパラメタ更新部２２２による更新を繰り返すように、各処理部を制御する。
＜動的分布推定装置２００の作用＞
次に、本実施の形態に係る動的分布推定装置２００の作用について説明する。まず、動的分布推定装置２００の初期化処理部１７は、パラメタ記録部２１３に格納されたパラメタと、統計量記録部２１６に格納された各統計量とを初期化する。そして、動的分布推定装置２００は、既に観測されたデータに基づいてＶＢアルゴリズムを用いて、モデルのパラメタと観測データ数と閾値と各統計量とを推定し、パラメタ記録部２１３、観測データ数記録部１４、閾値記録部１５、及び統計量記録部２１６へ格納する。そして、動的分布推定装置２００は、新たに観測されたデータｘｔが入力されると、図６に示す動的分布推定処理ルーチンを実行する。
まず、ステップＳ１００において、新たに観測されたデータｘｔを取得する。
ステップＳ１０２において、観測データ数Ｍと閾値Ｃとを更新する。
ステップＳ２０４において、潜在変数パラメタ更新部２１９は、新たに観測されたサンプルのデータｘｔに基づいて、潜在変数に関する変分分布のパラメタ
Figure 0006935765
Figure 0006935765
を、上記式（５５）に従って更新する。
ステップＳ２０６において、統計量更新部２２１は、上記ステップＳ２０４で更新された変分分布のパラメタに基づいて、
ステップＳ２０８において、統計量更新部２１は、上記ステップＳ１０４で更新された負担率γ（ｚｔ）に基づいて、各統計量
Figure 0006935765
を、上記式（４４）〜式（５０）に従って更新する。
ステップＳ２１０において、入出力部２４は、上記ステップＳ２０８で更新されたパラメタ
Figure 0006935765
を、外部装置３０へ出力して処理を終了する。
以上説明したように、第２の実施の形態に係る動的分布推定装置によれば、観測されるデータの分布を表す、複数のコンポーネントを混合した混合ガウスモデルのパラメタをオンラインで推定する。具体的には、第２の実施の形態に係る動的分布推定装置は、新たに観測されたサンプルのデータｘｔに基づいて、潜在変数に関する変分分布のパラメタと変分パラメタとを更新し、潜在変数に関する変分分布のパラメタ及び変分パラメタの少なくとも一方に基づいて各統計量を更新し、各統計量に基づいて、コンポーネントに関するパラメタを更新する。これにより、ＶＢアルゴリズムを用いる際に、高速であり、かつ省メモリであり、かつパラメタが時間連続性を有する状態で、モデルのパラメタを推定することができる。
なお、本発明は、上述した実施形態に限定されるものではなく、この発明の要旨を逸脱しない範囲内で様々な変形や応用が可能である。
例えば、上記の第２の実施形態では、パラメタ更新タイミングが新たに観測されたサンプルのデータが得られたタイミングとなる、逐次更新型オンラインＶＢアルゴリズムを用いたである場合を例に説明したがこれに限定されるものではない。例えば、パラメタ更新タイミングが、新たに観測されたサンプルのデータが予め定められた個数だけ得られたタイミングとなる、ミニバッチ型のオンラインＶＢアルゴリズムを用いてもよい。この場合には、予め定められた個数だけ新たに観測されたサンプルのデータｘｔに基づいて、潜在変数に関する変分分布のパラメタと変分パラメタとを更新し、潜在変数に関する変分分布のパラメタ及び変分パラメタの少なくとも一方に基づいて各統計量を更新し、各統計量に基づいて、コンポーネントに関するパラメタを更新すればよい。
また、パラメタ更新タイミングが予め定められた更新時期が到来したタイミングとなる、ミニバッチ型のオンラインＶＢアルゴリズムを用いてもよい。この場合には、更新時期が到来するまでの間に新たに観測されたサンプルのデータｘｔに基づいて、潜在変数に関する変分分布のパラメタと変分パラメタとを更新し、潜在変数に関する変分分布のパラメタ及び変分パラメタの少なくとも一方に基づいて各統計量を更新し、各統計量に基づいて、コンポーネントに関するパラメタを更新すればよい。
また、上記実施形態ではコンポーネントの分布としてガウス分布を利用する場合を例に説明したが、これに限定されず、任意の指数型分布族を利用する場合を含む。指数型分布族とは、密度関数が以下の式（６７）の形式で表されるものである。
Figure 0006935765
Figure 0006935765
は自然パラメタ、Ｔ（ｘ）は十分統計量、ｈ（ｘ）はベース測度、
Figure 0006935765
は対数分配関数と呼ばれる既知の関数であり、式中の記号”・”はベクトルの内積を表す。ガウス分布も指数型分布族に属する確率分布であり、以下の式（６８）のように自然パラメタ、十分統計量、ベース測度、対数分配関数が定義された場合、上記式（２）に示されるガウス分布の密度関数と上記式（６７）は等しくなる。
Figure 0006935765
これを踏まえると、混合ガウスモデルの推定アルゴリズムの中で計算する切断正規分布のモーメント（上記式（１２）及び上記式（１３）、上記式（６１）及び上記式（６２））は、正規分布を指数型分布族の形式で表現した場合に、十分統計量Ｔ（ｘ）の各次元に対応する値（ｘとｘの二乗）のｘが切断正規分布に従う場合の期待値を計算していると見なせる。ガウス分布以外の指数型分布に属する分布をコンポーネントの分布として利用する場合でも、モーメント計算の処理を十分統計量の各次元に対応する値の切断分布による期待値を計算する処理に置き変えることで、混合ガウス分布モデルの場合と同様に推定を行うことが可能である。なお、上記式（６１）及び上記式（６２）を計算する統計量更新部２２１は、期待値更新部の一例である。
また、本発明は、周知のコンピュータに媒体もしくは通信回線を介して、プログラムをインストールすることによっても実現可能である。
また、上述の装置は、内部にコンピュータシステムを有しているが、「コンピュータシステム」は、ＷＷＷシステムを利用している場合であれば、ホームページ提供環境（あるいは表示環境）も含むものとする。
また、本願明細書中において、プログラムが予めインストールされている実施形態として説明したが、当該プログラムを、コンピュータ読み取り可能な記録媒体に格納して提供することも可能である。
１０, ２１０ コンピュータ
１２,２１２ 記憶部
１３, ２１３ パラメタ記録部
１４ 観測データ数記録部
１５ 閾値記録部
１６, ２１６ 統計量記録部
１７, ２１７ 初期化処理部
１８, ２１８ 更新処理部
１９ 負担率更新部
２０ モーメント更新部
２１,２２１ 統計量更新部
２２, ２２２ パラメタ更新部
２３, ２２３ パラメタ処理部
２４ 入出力部
３０ 外部装置
１００,２００ 動的分布推定装置
２１９ 潜在変数パラメタ更新部

Claims (6)

  1. 観測されるデータの分布を表す、複数のコンポーネントを混合した混合ガウスモデルのパラメタをオンラインで推定する動的分布推定装置であって、
    新たに観測されたサンプルのデータに基づいて、前記新たに観測されたサンプルのデータが各コンポーネントに所属する度合いを表す負担率、及びまだ観測されていない各サンプルのデータが各コンポーネントに所属する度合いを表す負担率を更新する負担率更新部と、
    前記観測されていない各サンプルのデータが各コンポーネントに所属すると仮定した場合の、前記観測されていないサンプルのデータのモーメントを更新するモーメント更新部と、
    前記新たに観測されたサンプルのデータが各コンポーネントに所属する度合いを表す負担率に基づいて、観測されたサンプルのうち、各コンポーネントに所属するサンプル数の統計量を更新し、
    前記観測されていない各サンプルのデータが各コンポーネントに所属する度合いが表す負担率に基づいて、全サンプルのうち、各コンポーネントに所属するサンプル数の統計量を更新し、
    前記新たに観測されたサンプルのデータが各コンポーネントに所属する度合いが表す負担率に基づいて、各コンポーネントについて、前記コンポーネントに所属する、観測されたサンプルのデータの統計量を更新し、
    各コンポーネントについて、前記新たに観測された各サンプルのデータが前記コンポーネントに所属すると仮定した場合の、前記観測されていないサンプルのデータのモーメント、観測されたサンプルのうち、前記コンポーネントに所属するサンプル数の統計量、及び全サンプルのうち、前記コンポーネントに所属するサンプル数の統計量に基づいて、前記コンポーネントに所属する、前記観測されていないサンプルのデータの統計量を更新する統計量更新部と、
    各コンポーネントについて、全サンプルのうち、前記コンポーネントに所属するサンプル数の統計量、前記コンポーネントに所属する、前記観測されたサンプルのデータの統計量、及び前記コンポーネントに所属する、前記観測されていないサンプルのデータの統計量に基づいて、前記コンポーネントに関するパラメタを更新するパラメタ更新部と、
    を含み、
    予め定められたパラメタ更新タイミングが到来する毎に、前記負担率更新部による更新、前記モーメント更新部による更新、前記統計量更新部による更新、及び前記パラメタ更新部による更新を繰り返す
    動的分布推定装置。
  2. 観測されるデータの分布を表す、複数のコンポーネントを混合した混合ガウスモデルのパラメタをオンラインで推定する動的分布推定装置であって、
    新たに観測されたサンプルのデータに基づいて、前記新たに観測されたサンプルのデータについての、各コンポーネントの潜在変数に関する変分分布のパラメタ、及び新たに観測されたサンプルを含む、既に観測されたサンプル集合のデータについての、各コンポーネントの潜在変数に関する変分分布のパラメタを更新する潜在変数パラメタ更新部と、
    前記新たに観測されたサンプルのデータについての、各コンポーネントに対する変分分布のパラメタに基づいて、観測されたサンプルのうち、各コンポーネントに所属するサンプル数の統計量を更新し、
    既に観測されたサンプル集合のデータについての、各コンポーネントに対する変分分布のパラメタに基づいて、まだ観測されていないサンプルのうち、各コンポーネントに所属するサンプル数の統計量を更新し、
    前記新たに観測されたサンプルのデータについての、各コンポーネントに対する変分分布のパラメタに基づいて、各コンポーネントについて、前記コンポーネントに所属する、観測されたサンプルのデータの統計量を更新し、
    前記既に観測されたサンプル集合のデータと、前記観測されていないサンプルのうち、各コンポーネントに所属するサンプル数の統計量とに基づいて、前記コンポーネントに所属する、前記観測されていないサンプルのデータの統計量を更新する統計量更新部と、
    各コンポーネントについて、全サンプルのうち、前記コンポーネントに所属するサンプル数、前記コンポーネントに所属する、前記観測されたサンプルのデータの統計量、及び前記コンポーネントに所属する、前記観測されていないサンプルのデータの統計量に基づいて、前記コンポーネントのパラメタに関する変分分布のパラメタを更新するパラメタ更新部と、
    を含み、
    予め定められたパラメタ更新タイミングが到来する毎に、前記潜在変数パラメタ更新部による更新、前記統計量更新部による更新、及び前記パラメタ更新部による更新を繰り返す
    動的分布推定装置。
  3. 前記パラメタ更新タイミングは、前記新たに観測されたサンプルのデータが得られたタイミング、前記新たに観測されたサンプルのデータが予め定められた個数だけ得られたタイミング、及び予め定められた更新時期が到来したタイミングの何れかである請求項１又は請求項２に記載の動的分布推定装置。
  4. 観測されるデータの分布を表す、複数のコンポーネントを混合した混合ガウスモデルのパラメタをオンラインで推定する動的分布推定装置における動的分布推定方法であって、
    負担率更新部が、新たに観測されたサンプルのデータに基づいて、前記新たに観測されたサンプルのデータが各コンポーネントに所属する度合いを表す負担率、及びまだ観測されていない各サンプルのデータが各コンポーネントに所属する度合いを表す負担率を更新するステップと、
    モーメント更新部が、前記観測されていない各サンプルのデータが各コンポーネントに所属すると仮定した場合の、前記観測されていないサンプルのデータのモーメントを更新するステップと、
    統計量更新部が、前記新たに観測されたサンプルのデータが各コンポーネントに所属する度合いを表す負担率に基づいて、観測されたサンプルのうち、各コンポーネントに所属するサンプル数の統計量を更新し、
    前記観測されていない各サンプルのデータが各コンポーネントに所属する度合いが表す負担率に基づいて、全サンプルのうち、各コンポーネントに所属するサンプル数の統計量を更新し、
    前記新たに観測されたサンプルのデータが各コンポーネントに所属する度合いが表す負担率に基づいて、各コンポーネントについて、前記コンポーネントに所属する、観測されたサンプルのデータの統計量を更新し、
    各コンポーネントについて、前記新たに観測された各サンプルのデータが前記コンポーネントに所属すると仮定した場合の、前記観測されていないサンプルのデータのモーメント、観測されたサンプルのうち、前記コンポーネントに所属するサンプル数の統計量、及び全サンプルのうち、前記コンポーネントに所属するサンプル数の統計量に基づいて、前記コンポーネントに所属する、前記観測されていないサンプルのデータの統計量を更新するステップと、
    パラメタ更新部が、各コンポーネントについて、全サンプルのうち、前記コンポーネントに所属するサンプル数の統計量、前記コンポーネントに所属する、前記観測されたサンプルのデータの統計量、及び前記コンポーネントに所属する、前記観測されていないサンプルのデータの統計量に基づいて、前記コンポーネントに関するパラメタを更新するステップと、
    を含み、
    予め定められたパラメタ更新タイミングが到来する毎に、前記負担率更新部による更新、前記モーメント更新部による更新、前記統計量更新部による更新、及び前記パラメタ更新部による更新を繰り返す
    動的分布推定方法。
  5. 観測されるデータの分布を表す、複数のコンポーネントを混合した混合ガウスモデルのパラメタをオンラインで推定する動的分布推定装置における動的分布推定方法であって、
    潜在変数パラメタ更新部が、新たに観測されたサンプルのデータに基づいて、前記新たに観測されたサンプルのデータについての、各コンポーネントの潜在変数に関する変分分布のパラメタ、及び新たに観測されたサンプルを含む、既に観測されたサンプル集合のデータについての、各コンポーネントの潜在変数に関する変分分布のパラメタを更新するステップと、
    統計量更新部が、前記新たに観測されたサンプルのデータについての、各コンポーネントに対する変分分布のパラメタに基づいて、観測されたサンプルのうち、各コンポーネントに所属するサンプル数の統計量を更新し、
    既に観測されたサンプル集合のデータについての、各コンポーネントに対する変分分布のパラメタに基づいて、まだ観測されていないサンプルのうち、各コンポーネントに所属するサンプル数の統計量を更新し、
    前記新たに観測されたサンプルのデータについての、各コンポーネントに対する変分分布のパラメタに基づいて、各コンポーネントについて、前記コンポーネントに所属する、観測されたサンプルのデータの統計量を更新し、
    前記既に観測されたサンプル集合のデータと、前記観測されていないサンプルのうち、各コンポーネントに所属するサンプル数の統計量とに基づいて、前記コンポーネントに所属する、前記観測されていないサンプルのデータの統計量を更新するステップと、
    パラメタ更新部が、各コンポーネントについて、全サンプルのうち、前記コンポーネントに所属するサンプル数、前記コンポーネントに所属する、前記観測されたサンプルのデータの統計量、及び前記コンポーネントに所属する、前記観測されていないサンプルのデータの統計量に基づいて、前記コンポーネントのパラメタに関する変分分布のパラメタを更新するステップと、
    を含み、
    予め定められたパラメタ更新タイミングが到来する毎に、前記潜在変数パラメタ更新部による更新、前記統計量更新部による更新、及び前記パラメタ更新部による更新を繰り返す
    動的分布推定方法。
  6. コンピュータを、請求項１〜請求項３の何れか１項に記載の動的分布推定装置の各部として機能させるためのプログラム。
JP2018023593A 2018-02-13 2018-02-13 動的分布推定装置、方法、及びプログラム Active JP6935765B2 (ja)

Priority Applications (3)

Application Number Priority Date Filing Date Title
JP2018023593A JP6935765B2 (ja) 2018-02-13 2018-02-13 動的分布推定装置、方法、及びプログラム
PCT/JP2019/004677 WO2019159845A1 (ja) 2018-02-13 2019-02-08 動的分布推定装置、方法、及びプログラム
US16/969,052 US20210035000A1 (en) 2018-02-13 2019-02-08 Dynamic distribution estimation device, method, and program

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
JP2018023593A JP6935765B2 (ja) 2018-02-13 2018-02-13 動的分布推定装置、方法、及びプログラム

Publications (2)

Publication Number Publication Date
JP2019139597A JP2019139597A (ja) 2019-08-22
JP6935765B2 true JP6935765B2 (ja) 2021-09-15

Family

ID=67619409

Family Applications (1)

Application Number Title Priority Date Filing Date
JP2018023593A Active JP6935765B2 (ja) 2018-02-13 2018-02-13 動的分布推定装置、方法、及びプログラム

Country Status (3)

Country Link
US (1) US20210035000A1 (ja)
JP (1) JP6935765B2 (ja)
WO (1) WO2019159845A1 (ja)

Families Citing this family (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP7438008B2 (ja) 2020-04-28 2024-02-26 花王株式会社 商品の販売数予測方法、装置及びプログラム、発注量決定方法、装置及びプログラム

Also Published As

Publication number Publication date
JP2019139597A (ja) 2019-08-22
US20210035000A1 (en) 2021-02-04
WO2019159845A1 (ja) 2019-08-22

Similar Documents

Publication Publication Date Title
CN106503623B (zh) 基于卷积神经网络的人脸图像年龄估算方法
Fan et al. Learning continuous-time social network dynamics
TW201237647A (en) Method and system for identifying rare-event failure rates
JP5870189B1 (ja) 個別電気機器稼働状態推定装置、およびその方法
US11593693B2 (en) System and method for tree-based machine learning
Choi et al. Efficient targeting of sensor networks for large-scale systems
Zougab et al. A Bayesian approach to bandwidth selection in univariate associate kernel estimation
Lim Stochastic approximation over multidimensional discrete sets with applications to inventory systems and admission control of queueing networks
JP6935765B2 (ja) 動的分布推定装置、方法、及びプログラム
WO2020013236A1 (ja) データ解析装置、方法、及びプログラム
JP2001312712A (ja) 非線形時系列予測方法及び非線形時系列予測プログラムを記録した記録媒体
Hasler et al. Vine Copulas for Imputation of Monotone Non‐response
JP6099099B2 (ja) 収束判定装置、方法、及びプログラム
JP2017192041A (ja) 品質劣化推定装置、品質劣化推定方法、及びプログラム
Landauskas et al. Modelling of stock prices by the Markov chain Monte Carlo method
Noubary A procedure for prediction of sports records
JP6558860B2 (ja) 推定装置、予測装置、方法、及びプログラム
Shittaka et al. An evolutionary game model and MCMC estimation for analyzing stochastic properties of traffic state on a road network
Suparman et al. Hierarchical Bayesian of ARMA models using simulated annealing algorithm
Xu et al. Modeling and predicting chaotic circuit data
Lee A Bayesian inference for fixed effect panel probit model
Galdino Interval continuous-time Markov chains simulation
Ouzienko et al. Prediction of attributes and links in temporal social networks
WO2022074711A1 (ja) 学習方法、推定方法、学習装置、推定装置、及びプログラム
Oh Stochastic population projections on an uncertainty for the future Korea

Legal Events

Date Code Title Description
A621 Written request for application examination

Free format text: JAPANESE INTERMEDIATE CODE: A621

Effective date: 20200225

A131 Notification of reasons for refusal

Free format text: JAPANESE INTERMEDIATE CODE: A131

Effective date: 20210427

A521 Request for written amendment filed

Free format text: JAPANESE INTERMEDIATE CODE: A523

Effective date: 20210625

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

Free format text: JAPANESE INTERMEDIATE CODE: A01

Effective date: 20210727

A61 First payment of annual fees (during grant procedure)

Free format text: JAPANESE INTERMEDIATE CODE: A61

Effective date: 20210809

R150 Certificate of patent or registration of utility model

Ref document number: 6935765

Country of ref document: JP

Free format text: JAPANESE INTERMEDIATE CODE: R150