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

Description

本発明は、行列形式で与えられるデータから、データが従う分布を推定する分布推定装置、分布推定方法、及び分布推定プログラムに係り、特に、ユーザの特定の場所における滞在時間のデータから、データが従う分布を推定する分布推定装置、分布推定方法、及び分布推定プログラムに関する。

ＰＯＳ（Point of Sales）データに代表される購買履歴等の構造化されたデータ、テキストデータ・画像データ等の構造化されていないデータ等の多くは、前処理によって実数値を要素に持つ行列形式で表現できることが知られている。このように行列表現されたデータ中に存在するクラスタを発見する手法として、非負値行列分解（Non-negative Matrix Factorization, NMF）と呼ばれる手法の有用性がこれまで示されている（例えば、非特許文献１を参照）。

ＮＭＦを適用する際に入力される行列データは、当該行列データより低次の階数の行列の積に分解される。この各々の低次の階数の行列が、それぞれ各行、各列に対応する事物のクラスタへの寄与度を表しており、この寄与度によりクラスタ発見が可能となる。従って、例えば、ユーザの訪問履歴に関するデータを行列形式により表現することにより、おすすめのショップのリストを作成すること等が可能となる。

図１０に、ユーザの訪問履歴Ｘに関するデータをＮＭＦに適用した適応例を示す。図１０に示すように、ユーザの場所毎の訪問回数を表すユーザ訪問場所行列Ｘは、行列中の各行ｉがユーザを表し、各列ｊが訪問場所を表し、各値が訪問回数を表すＩ行Ｊ列の行列である。図１０に示す例では、１行目に対応するユーザがユーザ１であり、１列目に対応する訪問場所が場所１であり、ユーザ１の場所１への訪問回数が４である。

ユーザ訪問場所行列Ｘにおいては、訪問回数の値が大きいほど、その場所が人気のスポットであることを表す。このように、ユーザ訪問場所行列ＸをＮＭＦに適用することで、

となるＩ行Ｒ列のユーザ特徴行列Ａ＝｛ａ_ｉｒ｝と、Ｊ行Ｒ列の訪問場所特徴行列Ｂ＝｛ｂ_ｊｒ｝が求まる。ただし、記号の上付きの記号Ｔは行列の転置を表す。

また、ここで、記号

で表現した類似の尺度について説明する。上記非特許文献１にも記述されているように、行列の類似の尺度には、ユークリッド距離に基づくものや一般化カルバックライブラーダイバージェンス（ＫＬ距離）により定義される距離尺度が用いられ、値が小さいほど両者が類似していることを表す。

図１０において、ユーザ特徴行列Ａのクラスタ１に対応する１列目に着目すると、ユーザ１に対応する１行目、ユーザ２に対応する２行目、及びユーザ３に対応する３行目の値がそれぞれ０より大きい値となっている。これは、ユーザ１、ユーザ２、及びユーザ３がクラスタ１に所属することを示している。

また、図１０において、訪問場所特徴行列Ｂに着目すると、場所１に対応する１列目、場所２に対応する２列目、及び場所３に対応する３列目では、クラスタ１に対応する１行目の値がクラスタ２に対応する２行目の値より大きい値となっている。これは、場所１、場所２、及び場所３は、ユーザ１、ユーザ２、及びユーザ３が訪れやすい場所であるというクラスタ１が持つ特徴を表している。

これらを踏まえ、図１１に示すように、場所１、場所２、及び場所３をまとめてクラスタ１の場所特徴とする。同様に、クラスタ１に所属するユーザ１、ユーザ２、及びユーザ３をクラスタ１のユーザ特徴とする。以下、クラスタ１の場所特徴及びユーザ特徴をまとめてクラスタ１の特徴ともいう。

このように、ＮＭＦの適用によって得られたユーザ特徴行列Ａ、及び訪問場所特徴行列Ｂに基づき、図１１に示すようなクラスタ抽出が可能となる。

なお、クラスタの総数に相当する訪問場所特徴行列Ｂの階数は、解析する前に予め決定しておくものとする。

ＮＭＦでは、行列を分解する際に、行列に対して類似の尺度を定義しており、定義した尺度のもとで値が小さいほど両者が類似していることを表す。従って、ＮＭＦでは、採用した類似の尺度を最小化する行列（例えば、ユーザ特徴行列Ａ、及び訪問場所特徴行列Ｂを求める手法として定式化される。

また、ＮＭＦでは、データが持つ性質を考慮して、利用する距離が決定される。例えば、距離としてＫＬダイバージェンスが採用される場合は、行列の各要素ｘ_ｉｊは平均μ＝ｘ＾_ｉｊのポアソン分布Ｐｏ（ｘ_ｉｊ｜μ_ｉｊ）に従っていると仮定していることに相当する。ポアソン分布Ｐｏ（ｘ_ｉｊ｜μ_ｉｊ）に従っているデータとしては、上述したユーザの場所毎の訪問回数等が例示される。

なお、行列Ｘと行列Ｘ＾のＫＬダイバージェンスＤ_ＫＬは、下記（１）式で定義される。

実数値が従う確率分布としては正規分布が有用であり、頻度を表す離散値が従う確率分布としてはポアソン分布が有用であることは広く認識された事実である。

澤田宏, "非負値行列因子分解NMFの基礎とデータ／信号解析への応用", 電子情報通信学会誌, Vol. 95, No. 9, pp. 829-833, 2012. Kyosuke Nishida, Hiroyuki Toda, Takeshi Kurashima, and Yoshihiko Suhara, "Probabilistic Identification of Visited Point-of-Interest for Personalized Automatic Check-in", pp.631-642, in Proceedings of the 2014 ACM International Joint Conference on Pervasive and Ubiquitous Computing (UbiComp 2014).

従来、ある時間帯にある場所を訪れたユーザがどのくらいの時間その場所に滞在したかを表す入力データから、滞在時間の従う確率分布の推定を行う技術が用いられている。この場合の入力データは、図１２の左図に示すように、各行ｉが時間帯を表し、各列ｊが場所を表し、要素としてベクトルｘ_ｉｊ＝（ｘ_ｉｊ１，ｘ_ｉｊ２，…，ｘ_{ｉｊＫｉｊ}）を持つＩ行Ｊ列のベクトル値行列Ｘとして表現される。ただし、Ｋ_ｉｊは時間帯ｉに場所ｊを訪れたユーザの総数であり、ｘ_ｉｊｋが時間帯ｉに場所ｊにｋ番目に訪れた人の滞在時間を表す。以後、このベクトル値行列Ｘを時間帯場所ベクトル値行列という。

なお、ここでいう時間帯ｉは、例えばｉ＝１が［７：００−８：００］の間の１時間、ｉ＝２が［８：００−９：００］の１時間のように、２４時間をなんらかの基準で分割したものであってもよいし、例えばｉ＝１が休日の午前中、ｉ＝２が休日の午後、ｉ＝３が平日の午前中、ｉ＝４が平日の午後のように、平日及び休日といった日に関する情報と、時間に関する情報とを組み合わせて作成したものであってもよい。

また、場所ｊは、喫茶店、レストラン等の店舗に対応させても良く、地図を何らかの基準でグリッド化したグリッドのＩＤに対応させても良い。

このような入力データから、データが従う分布の推定を行う場合、時間帯場所ベクトル値行列Ｘを、通常の要素がスカラである行列に変換する平均化処理を適用した上で、既存のＮＭＦを適用することが考えられる。より正確に述べると、時間帯場所ベクトル値行列Ｘにおいてユーザ毎に平均を取った時間帯場所平均行列

を下記（２）式に従って作成し、時間帯場所平均行列ＸにＮＭＦを適用する。

しかし、この手法では、各要素の分散の情報を潰しており、時間帯ｉに場所ｊに滞在した滞在時間の平均値のみに基づいて分布の推定を行うことになる。すなわち、この手法では、以下の（Ａ）及び（Ｂ）により、滞在時間が従う確率分布を推定したい場合には不適切だといえる。

（Ａ）上記手法は、正規分布から得られたと仮定されるデータに対して推定を行うためである。上記非特許文献２でも示されるように、ユーザがある場所に滞在した滞在時間の分布は対数正規分布に従うことが知られている。つまり、上記手法により分布の推定を行うと、実際には図１２の右図に示すように、ユーザがある場所に滞在した滞在時間の分布は正規分布とは限らない。しかし、図１３に示すように、ユーザがある場所に滞在した滞在時間の分布を正規分布として推定を行ってしまい、滞在時間が従う分布の推定を精度良く行うことができない。

（Ｂ）上記手法では、平均値を推定しているため、分散値も考慮した分布の形を推定する場合には不適切である。

本発明は、以上のような事情に鑑みてなされたものであり、要素の分散を考慮して、分布を推定することができる分布推定装置、分布推定方法、及び分布推定プログラムを提供することを目的とする。

上記目的を達成するために、本発明の分布推定装置は、個々を識別可能な第１の個体群に含まれる個体ｉ（１≦ｉ≦Ｉ，Ｉは１以上の整数）とオブジェクトｊ（１≦ｊ≦Ｊ、Ｊは１以上の整数）との関連度を表すと共に、平均パラメータμ_ｉｊ及び分散パラメータを持つ対数正規分布に従う要素ｘ_ｉｊｋからなるベクトルｘ_ｉｊを要素として持つＩ×Ｊのベクトル値行列Ｘから、前記個体ｉが、クラスタｒ（１≦ｒ≦Ｒ、Ｒは１以上の整数）に所属することを表す非負値の要素ａ_ｉｒを持つＩ×Ｒの第１の特徴行列Ａと、前記オブジェクトｊが、前記クラスタｒに所属することを表す非負値の要素ｂ_ｊｒを持つＪ×Ｒの第２の特徴行列Ｂと、を抽出する分布推定装置であって、前記ベクトル値行列Ｘの各ベクトルｘ_ｉｊの各要素ｘ_ｉｊｋの対数値、前記第１の特徴行列Ａ、前記第２の特徴行列Ｂ、及び前記分散パラメータを用いて表される目的関数を最適化するように、前記第１の特徴行列Ａ及び前記第２の特徴行列Ｂを推定する特徴行列推定部と、予め定めた反復終了条件を満足するまで、前記特徴行列推定部による推定を繰り返す反復判定部と、を備える。

なお、本発明のｖ装置において、前記ベクトルｘ_ｉｊの要素ｘ_ｉｊｋは、平均パラメータμ_ｉｊ及び個体ｉに依存する分散パラメータσ_ｉを持つ対数正規分布、又は、平均パラメータμ_ｉｊ及びオブジェクトｊに依存する分散パラメータσ_ｊを持つ対数正規分布に従い、前記特徴行列推定部は、前記ベクトル値行列Ｘの各ベクトルｘ_ｉｊの各要素ｘ_ｉｊｋの対数値、前記第１の特徴行列Ａ、前記第２の特徴行列Ｂ、及び前記分散パラメータσ_ｉ又は前記分散パラメータσ_ｊを用いて表される目的関数を最適化するように、前記第１の特徴行列Ａ及び前記第２の特徴行列Ｂを推定するようにしても良い。

また、本発明の分布推定装置において、前記ベクトル値行列Ｘは、時間帯ｉにおける各ユーザの場所ｊの滞在時間を表す要素ｘ_ｉｊｋからなるベクトルｘ_ｉｊを要素として持ち、前記第１の特徴行列Ａは、前記時間帯ｉが、前記クラスタｒに所属することを表す非負値の要素ａ_ｉｒを持ち、前記第２の特徴行列Ｂは、前記場所ｊが、前記クラスタｒに所属することを表す非負値の要素ｂ_ｊｒを持つようにしても良い。

また、本発明の分布推定装置において、前記ベクトルｘ_ｉｊの要素ｘ_ｉｊｋは、平均パラメータμ_ｉｊ及び場所ｊに依存する分散パラメータσ_ｊを持つ対数正規分布に従い、前記特徴行列推定部は、前記ベクトル値行列Ｘの各ベクトルｘ_ｉｊの各要素ｘ_ｉｊｋの対数値、前記第１の特徴行列Ａ、前記第２の特徴行列Ｂ、及び前記分散パラメータσ_ｊを用いて表される目的関数を最小化するように、前記第１の特徴行列Ａ及び前記第２の特徴行列Ｂを推定するようにしても良い。

また、本発明の分布推定装置において、前記ベクトルｘ_ｉｊの要素ｘ_ｉｊｋは、平均パラメータμ_ｉｊ、及び時間帯ｉに依存する分散パラメータσ_ｉを持つ対数正規分布に従い、前記特徴行列推定部は、前記ベクトル値行列Ｘの各ベクトルｘ_ｉｊの各要素ｘ_ｉｊｋの対数値、前記第１の特徴行列Ａ、前記第２の特徴行列Ｂ、及び前記分散パラメータσ_ｊを用いて表される目的関数を最小化するように、前記第１の特徴行列Ａ及び前記第２の特徴行列Ｂを推定するようにしても良い。

上記目的を達成するために、本発明の分布推定方法は、個々を識別可能な第１の個体群に含まれる個体ｉ（１≦ｉ≦Ｉ，Ｉは１以上の整数）とオブジェクトｊ（１≦ｊ≦Ｊ、Ｊは１以上の整数）との関連度を表すと共に、平均パラメータμ_ｉｊ及び分散パラメータを持つ対数正規分布に従う要素ｘ_ｉｊｋからなるベクトルｘ_ｉｊを要素として持つＩ×Ｊのベクトル値行列Ｘから、前記個体ｉが、クラスタｒ（１≦ｒ≦Ｒ、Ｒは１以上の整数）に所属することを表す非負値の要素ａ_ｉｒを持つＩ×Ｒの第１の特徴行列Ａと、前記オブジェクトｊが、前記クラスタｒに所属することを表す非負値の要素ｂ_ｊｒを持つＪ×Ｒの第２の特徴行列Ｂと、を抽出する分布推定装置における分布推定方法であって、特徴行列推定部が、前記ベクトル値行列Ｘの各ベクトルｘ_ｉｊの各要素ｘ_ｉｊｋの対数値、前記第１の特徴行列Ａ、前記第２の特徴行列Ｂ、及び前記分散パラメータを用いて表される目的関数を最適化するように、前記第１の特徴行列Ａ及び前記第２の特徴行列Ｂを推定する特徴行列推定ステップと、反復判定部が、予め定めた反復終了条件を満足するまで、前記特徴行列推定ステップによる推定を繰り返す反復判定ステップと、を行う。

また、本発明の分布推定方法において、前記ベクトル値行列Ｘは、時間帯ｉにおける各ユーザの場所ｊの滞在時間を表す要素ｘ_ｉｊｋからなるベクトルｘ_ｉｊを要素として持ち、前記第１の特徴行列Ａは、前記時間帯ｉが、前記クラスタｒに所属することを表す非負値の要素ａ_ｉｒを持ち、前記第２の特徴行列Ｂは、前記場所ｊが、前記クラスタｒに所属することを表す非負値の要素ｂ_ｊｒを持つようにしても良い。

上記目的を達成するために、本発明の文書分類プログラムは、コンピュータを、上記分布推定装置の各部として機能させるためのプログラムである。

本発明によれば、要素の分散を考慮して、分布を推定することができる。

実施形態に係る分布推定装置における場所及び滞在時間の分布を推定する推定方法の一例を示す機能ブロック図である。 対数正規分布の一例を示す模式図である。 実施形態に係る分布推定装置における時間帯場所ベクトル値行列を時間帯場所ベクトル値対行列に変換する変換方法の一例を示す模式図である。 実施形態に係る分布推定装置における時間帯場所ベクトル値対数行列を時間帯場所対数平均値ベクトルに変換する変換方法の一例を示す模式図である。 実施形態に係る分布推定装置における時間帯場所行列及び場所特徴行列を推定する推定方法の一例を示す模式図である。 実施形態に係る分布推定装置における時間帯場所ベクトル値行列の平均値を推定する推定方法の一例を示す模式図である。 実施形態に係る分布推定装置の構成の一例を示す機能ブロック図である。 実施形態に係る分布推定装置により実行される全体処理の流れを示すフローチャートである。 実施形態に係る分布推定装置により実行される特徴行列推定処理の流れを示すフローチャートである。 ユーザの訪問履歴に関するデータをＮＭＦに適用した行列分解の一例を示す模式図である。 ユーザの訪問履歴に関するデータをＮＭＦに適用した行列分解によりクラスタを抽出方法の一例を示す模式図である。 時間帯場所ベクトル値行列の一例を示す模式図である。 場所及び対外時間の分布を推定する推定方法の一例である。

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

本実施形態に係る分布推定装置は、時間帯場所ベクトル値行列Ｘにおける滞在時間が従う分布の推定を行う非負値行列分解を用いる。特に、本実施形態では、図１に示すように、要素である滞在時間を対数正規分布に従って得られると仮定して、ユーザの滞在時間の分布の推定を行う。

対数正規分布は、ユーザの滞在時間等のモデリングにおいて広く利用される確率分布である（上記非特許文献２を参照）。対数正規分布は、正規分布とは異なり、図２に示すように、右に裾の長いデータを表現することが可能な確率分布である。対数正規分布は、この性質によって、例えばカフェ等における『話し込んでしまい長居をしてしまう』（滞在時間が長くなる）ことを表現することが可能となる。

本実施形態では、時間帯場所ベクトル値行列Ｘにおける各要素ｘ_ｉｊｋが、平均パラメータμ_ｉｊ、分散パラメータσ_ｉｊを持つ対数正規分布ＬＮ（ｘ_ｉｊｋ｜μ_ｉｊ，σ_ｉｊ）に従い、平均パラメータμ_ｉｊは、時間帯特徴行列Ａ、場所特徴行列Ｂの要素を用いて下記（３）式のように表現されるとしてモデリングを行う。

時間帯特徴行列Ａ、場所特徴行列Ｂ、及び分散パラメータσ_ｉｊは、下記（４）式に示す目的関数を用いて最小化を行うことにより推定できる。

このように、本実施形態では、平均パラメータμ_ｉｊが、時間帯特徴行列Ａと場所特徴行列Ｂとの積で表現される。これは一見、入力行列である時間帯場所ベクトル値行列Ｘの行列分解を行う通常のＮＭＦとは全く異なる手法に見える。しかしながら、本実施形態の技術とＮＭＦには強いつながりがある。それは、対数正規分布と正規分布との関係、すなわち「対数正規分布ＬＮ（ｘ｜μ，σ）に従う変数Ｘの対数関数をとった場合、新たな確率変数Ｙ＝ｌｏｇ（Ｘ）が正規分布Ｎ（ｙ｜μ，σ）に従う」という関係に由来するものである。

まず、図３に示すように、時間帯場所ベクトル値行列Ｘを対数変換した時間帯場所ベクトル値対数行列Ｙを算出する。時間帯場所ベクトル値対数行列Ｙのｙ_ｉｊの要素には、訪問したＫ_ｉｊ人のユーザの訪問時間情報（対数値）が、時間帯場所ベクトル値行列Ｘと同様にＫ_ｉｊ次元ベクトル要素として格納されている。

また、図４に示すように、時間帯場所ベクトル値対数行列Ｙを、時間帯場所平均対数値行列

に変換する。時間帯場所平均対数値行列

の

の要素には、訪問したＫ_ｉｊ人のユーザの平均訪問時間情報（対数平均値）が格納されている。なお、この場合の

の要素は、ベクトルではなくスカラである。

本実施形態によるＮＭＦは、図５に示すように、ｙ_ｉｊとの誤差を上記（４）式に示す目的関数に従って最小化するように、ＮＭＦによって因子行列である時間帯特徴行列Ａ及び場所特徴行列Ｂを導出する。これらの因子行列によって再構築されるμ_ｉｊの要素は、上記（４）式に従って、分散パラメータσ_ｉｊも小さくするようにパラメータが抽出される。なお、分散パラメータσ_ｉｊは、ｙ_ｉｊｋが平均パラメータμ_ｉｊからどれだけばらつくかという分散を表したパラメータとなっている。このように、本実施形態では、ＮＭＦの拡張技術であるとみなすことができる。

さらに、上記（４）式では、データ数Ｋ_ｉｊの和を取っていることがわかる。このことから、データ数Ｋ_ｉｊが多ければ多いほど、目的関数の値が大きくなってしまう。そのため、目的関数を小さくするためには分散パラメータσ_ｉｊを小さくすることになる。すなわち、データ数Ｋ_ｉｊが多い場所に関する分布を推定する場合ほど、分散パラメータσ_ｉｊが小さくなり、正確な推定を行えるような目的関数となっている。

このように、本実施形態では、時間帯場所ベクトル値行列Ｘを入力行列とした場合の、各場所及び各時間に対するユーザの滞在時間が従う分布の推定を行うことで、滞在時間の推定が行える。また、対数正規分布という滞在時間のモデリングに適した分布を利用していることから、ユーザの滞在時間の推定精度の向上も期待できる。これは、図６に示すように、時間帯ｉに場所ｊに滞在する時間を表すベクトルｘ_ｉｊを生成する確率が最も高い対数正規分布のパラメータを推定していることに等しい。パラメータ推定結果を用いることで、時間帯及び場所毎の平均滞在時間を対数正規分布の平均値

として推定できる。なお、

の要素の平均値は、下記（５）式で表される。

また、分散パラメータσ_ｉｊを時間帯及び場所によって異なるように取得した場合、例えば、「テーマパークＤでは、午前中早くに来場する客は長時間楽しむことを目的としているが、夕方に来場する客は、より安価に購入できるナイトパスを利用して、短時間で少数のアトラクションを楽しむ傾向にある」という分析及びパターンの抽出が可能になる。

本実施形態では、さらに、特殊な場合として、時間帯場所ベクトル値対数行列Ｙのｊ列目の任意のベクトルｘ_ｉｊに対して分散パラメータσ_ｉｊが共通である場合についても考えることもできる。この場合には分散パラメータσ_ｉｊはｊにしか依存しないので、分散パラメータσ_ｊとして、上述した方法を適用することができる。この場合には、最小化の目的関数は、下記（６）式のように表される。

分散パラメータσ_ｉｊをｊ成分について共通化するということは、ｉ成分に時間帯の情報を持ち、ｊ成分に場所の情報を持つデータにおいて、分散パラメータσ_ｉｊは時間帯に依存しないと考えていることとなる。データの密度によって偏りがあるデータ、つまり場所及び時間帯によって偏りがあるデータ等では、このようにパラメータを共通化することは非常に有効である。

これにより、分散パラメータσ_ｉｊを、時間帯依存がなく、場所のみに特徴があると仮定して抽出した場合、例えば「サービスエリアＡでは、どの時間帯でも平均的に１時間滞在する傾向があるので、データ数は少ないが、同じような利用をされている場所が、道の駅、国道沿いのコンビニエンスストア等でも抽出できた。』という分析及びパターンの抽出が可能になる。

また、同様に、ｉ行の任意ベクトル値に対して分散パラメータσ_ｉｊが共通である場合を考えた場合、分散パラメータσ_ｉｊを分散パラメータσ_ｉと書き直し、最小化する目的関数を下記（７）式のように設定することで、上述した方法を適応することができる。

これにより、分散パラメータσ_ｉｊを、場所依存がなく、時間帯のみに特徴があると仮定して抽出した場合、例えば「ある因子に属する場所では、昼の時間帯は３０分ほどなのに、夜の利用客は２時間以上の利用がされていて、これは多くのレストランでは、昼は客の回転数を上げる努力がなされている一方、夜は客の単価を上げるように努力している」という分析及びパターンの抽出が可能になる。

次に、本実施形態に係る分布推定装置１の機能について説明する。

図７に示すように、本実施形態に係る分布推定装置１は、時間帯場所情報処理部１０、特徴行列推定部２０、反復判定部２０ａ、特徴行列処理部３０、記憶部４０、及び、入出力部５０を有する。入出力部５０には、入力装置、表示装置等の外部装置２に接続されており、入出力部５０は、外部装置２に対して情報の入出力を行う。

記憶部４０は、時間帯場所情報テーブル４１、時間帯特徴テーブル４２、場所特徴テーブル４３、及び、場所分散テーブル４４を有する。

以下に各テーブルについて説明する。なお、テーブル形式のデータは行列形式にて表現できることから、以下では、各テーブルと各特徴行列を同一視し、区別せずに用いる。また、時間帯場所行列の各成分はベクトル値であり、各特徴行列の各成分は実数値であるが、同じ「行列」という表現を用い、文脈によって使い分ける。

＜時間帯場所情報テーブル４１＞

時間帯場所情報テーブル４１は、時間帯フィールド、場所ＩＤフィールド、及び、滞在時間ベクトルフィールドを有する。時間帯フィールドは、時間帯場所情報処理部１０により追加された時間帯を特定する識別子が設定されるフィールドである。場所ＩＤフィールドには、時間帯場所情報処理部１０により追加された場所を特定する識別子が設定されるフィールドである。滞在時間ベクトルフィールドは、時間帯場所情報処理部１０により各々の時間帯に各々の場所に滞在した複数人のユーザの滞在時間の値（以下、滞在時間ベクトルという。）が設定されるフィールドである。なお、各ユーザの滞在時間の値には、０又は正の実数を設定できるが、負の実数又は虚数を設定することはできない。

＜時間帯特徴テーブル４２＞

時間帯特徴テーブル４２は、時間帯フィールド、平均値基底ベクトルＩＤフィールド、及び、時間帯特徴値フィールドを有する。時間帯フィールドは、特徴行列推定部２０により時間帯を特定する識別子が設定されるフィールドである。平均値基底ベクトルＩＤフィールドは、特徴行列推定部２０により各データが従う分布の平均値を表現する基底を特定する識別子が設定されるフィールドである。時間帯特徴値フィールドは、特徴行列推定部２０により算出された当該時間帯の平均値を表現する基底ベクトルの特徴値が設定されるフィールドである。

＜場所特徴テーブル４３＞

場所特徴テーブル４３は、場所ＩＤフィールド、係数ＩＤフィールド、及び、場所特徴値フィールドを有する。場所ＩＤフィールドは、特徴行列推定部２０により場所を特定する識別子が設定されるフィールドである。係数ＩＤフィールドは、特徴行列推定部２０によりどの平均値基底ベクトルに対する係数かを特定する識別子が設定されるフィールドである。場所特徴値フィールドは、特徴行列推定部２０により算出された当該場所の平均値基底ベクトルに対する重みを特定する特徴値が設定されるフィールドである。

＜場所分散テーブル４４＞

場所分散テーブル４４は、時間帯フィールド、場所ＩＤフィールド、及び、場所分散値フィールドを有する。時間帯フィールドは、特徴行列推定部２０により時間帯を特定する識別子が設定されるフィールドである。場所ＩＤフィールドは、特徴行列推定部２０により場所を特定する識別子が設定されるフィールドである。場所分散値フィールドは、特徴行列推定部２０により算出された当該場所の場所分散パラメータの値が設定されるフィールドである。

時間帯場所情報処理部１０は、入出力部５０から、入力行列として、時間帯場所行列を入力する。また、時間帯場所情報処理部１０は、入力した時間帯場所行列に基づき、時間帯場所情報テーブル４１に、追加された場所、滞在時間に応じて、時間帯フィールド、場所ＩＤフィールド、滞在時間の値を設定した行を挿入する。

特徴行列推定部２０は、時間帯場所情報テーブル４１に格納された情報を取得し、後述する特徴行列推定処理を行うことにより特徴行列を推定し、得られた時間帯の特徴値を時間帯特徴テーブル４２に格納し、場所の特徴値を場所特徴テーブル４３に格納し、場所分散パラメータの値を場所分散テーブル４４に格納する。特徴行列推定部２０が特徴行列を推定する際には、時間帯（個体）ｉと場所（オブジェクト）ｊとの関連度を表すと共に、平均パラメータμ_ｉｊ及び分散パラメータσ_ｉｊを持つ対数正規分布に従う要素ｘ_ｉｊｋからなるベクトルｘ_ｉｊを要素として持つＩ×Ｊの時間帯場所ベクトル値行列Ｘの各ベクトルｘ_ｉｊの各要素ｘ_ｉｊｋの対数値、時間帯ｉがクラスタｒに所属することを表す非負値の要素ａ_ｉｒを持つＩ×Ｒの時間帯特徴行列（第１の特徴行列）Ａ、場所ｊがクラスタｒに所属することを表す非負値の要素ｂ_ｊｋを持つＪ×Ｒの場所特徴行列（第２の特徴行列）Ｂ、及び分散パラメータσ_ｉｊを用いて表される上記（４）式の目的関数を最小化するように、時間帯特徴行列Ａ及び場所特徴行列Ｂを推定する。

反復判定部２０ａは、予め定めた反復終了条件を満足するまで、特徴行列推定部２０による推定を繰り返すように制御する。

特徴行列処理部３０は、時間帯特徴テーブル４２、及び場所特徴テーブル４３を参照し、リクエストの引数に対応する特徴行列を出力する。本実施形態では、例えば、外部装置２から特徴出力のリクエストが入力された場合に、特徴行列の出力する実行するが、特徴行列の出力を実行するタイミングはこれに限らず、予め定めた時間が経過する毎に、特徴行列を出力しても良い。また、出力する特徴は全ての特徴であっても一部の特徴であっても良く、全ての特徴を出力する場合には、時間帯特徴テーブル４２、及び場所特徴テーブル４３の全ての行を出力すればよい。

入出力部５０は、外部装置２から、入力行列として時間帯場所行列を入力し、入力した時間帯場所行列を時間帯場所情報処理部１０に引き渡す。また、入出力部５０は、時間帯特徴テーブル４２に格納された特徴値、及び場所特徴テーブル４３に格納された特徴値を外部装置２に対して出力する。

なお、本実施形態に係る分布推定装置１は、例えば、ＣＰＵ（Central Processing Unit）、ＲＡＭ（Random Access Memory）、各種プログラムを記憶するＲＯＭ（Read Only Memory）を備えたコンピュータ装置で構成される。また分布推定装置１を構成するコンピュータは、ハードディスクドライブ、不揮発性メモリ等の記憶部を備えていても良い。本実施形態では、ＣＰＵがＲＯＭ、ハードディスク等の記憶部に記憶されているプログラムを読み出して実行することにより、上記のハードウェア資源とプログラムとが協働し、上述した機能が実現される。

本実施形態に係る分布推定装置１は、時間帯場所行列を入力として特徴行列を推定し、特徴行列を出力する全体処理を行う。

まず、本実施形態に係る分布推定装置１による全体処理の流れを、図８に示すフローチャートを用いて説明する。

本実施形態では、時間帯場所情報処理部１０が時間帯場所情報を更新するタイミングは、例えば、システム管理者が外部装置２から供給されるデータをもとに手動で管理できるようにしてもよいし、あるユーザが新たな場所に滞在した場合に外部装置２が自動的に処理を起動するようにしてもよい。

ステップＳ１０１では、時間帯場所情報処理部１０が、入出力部５０から、入力行列として、時間帯場所行列を入力する。また、時間帯場所情報処理部１０が、入力した時間帯場所行列に基づき、時間帯場所情報テーブル４１に、追加された場所、滞在時間に応じて、時間帯フィールド、場所ＩＤフィールド、滞在時間の値を設定した行を挿入する。

ステップＳ１０３では、特徴行列推定部２０が、後述する特徴行列推定処理を行うことにより、各特徴行列、及び場所分散を推定する。

ステップＳ１０５では、特徴行列処理部３０が、時間帯特徴テーブル４２、及び場所特徴テーブル４３を参照し、リクエストの引数に対応する特徴行列を出力する。本実施形態では、例えば、外部装置２から特徴出力のリクエストが入力された場合に、特徴行列の出力する実行するが、特徴行列の出力を実行するタイミングはこれに限らず、予め定めた時間が経過する毎に、特徴行列を出力しても良い。また、出力する特徴は全ての特徴であっても一部の特徴であっても良く、全ての特徴を出力する場合には、時間帯特徴テーブル４２、及び場所特徴テーブル４３の全ての行を出力すればよい。

ここで、特徴行列推定部２０は、上述した特徴行列推定処理として、以下の方法で各特徴行列、場所分散を推定し、記憶部４０の時間帯特徴テーブル４２、場所特徴テーブル４３、及び場所分散テーブル４４に格納する処理を行う。

本実施形態に係る分布推定装置１による特徴行列推定処理の流れを、図９に示すフローチャートを用いて説明する。

なお、本実施形態では、時間帯場所情報テーブル４１中に存在する全データを下記（８）式のように表す。また、時間帯特徴行列Ａ、及び場所特徴行列Ｂをそれぞれ下記（９）式及び（１０）式のように表す。また、場所分散パラメータσを下記（１１）式のように表す。各式におけるＩは全時間帯数であり、Ｊは全場所数である。また、各式におけるｉは時間帯を特定する識別子であり、ｊは場所を特定する識別子であり、ｒは平均値を表現する基底ベクトルを特定する識別子である。

ステップＳ２０１では、特徴行列推定部２０が、時間帯特徴行列Ａ、及び場所特徴行列Ｂをそれぞれ初期化する。同様に終了条件の閾値ε、最大繰り返し回数を設定する。

ステップＳ２０３では、特徴行列推定部２０が、終了条件に用いる変数として特徴更新の最大変化幅を示す変数δをδ＝０として初期化する。

ステップＳ２０５では、特徴行列推定部２０が、下記（１２）式に従い、時間帯特徴行列Ａを更新する。その後、更新前の時間帯特徴行列Ａの要素の値と更新後の時間帯特徴行列Ａの要素の値の差の絶対値の最大値である

がδより大きければ、

と更新する。

ただし、

と書いた、μ_ｉｊは、時間帯特徴行列Ａ、及び場所特徴行列Ｂによるｌｏｇ（ｘ_ｉｊｋ）の推定値であると見なせる。

なお、上記（１２）式における記号「←」は、右辺の計算結果を左辺の変数に代入する処理を意味する。また、各式において、代入処理前のユーザ特徴行列Ａの要素の値をａ^ｏｌｄ _ｉｒとし、代入処理後の値をａ^ｎｅｗ _ｉｒとして記述している。また、上記（１２）式におけるμ_ｉｊは、時間帯特徴行列Ａ及び場所特徴行列Ｂによるｌｏｇ（ｘ_ｉｊｋ）の推定値であると見なせる。

ステップＳ２０７では、特徴行列推定部２０が、下記（１３）式に従い、場所特徴行列Ｂを更新する。その後、更新前の場所特徴行列Ｂの要素の値と更新後の場所特徴行列Ｂの要素の値の差の絶対値の最大値

がδより大きければ、

と更新する。なお、各式において、代入処理前の場所特徴行列Ｂの要素の値をｂ^ｏｌｄ _ｉｒとし、代入処理後の値をｂ^ｎｅｗ _ｉｒとして記述している。

ステップＳ２０９では、特徴行列推定部２０が、下記（１４）式に従い、場所分散パラメータσ_ｉｊを更新する。

ステップＳ２１１では、特徴行列推定部２０が、計算繰り返し回数を更新する。

ステップＳ２１３では、反復判定部２０ａが、計算繰り返し回数が予め定めた最大繰り返し数を超えたか否か、又は、特徴行列を更新することにより、上述した変数δが予め定めた閾値εより小さいか否かを判定する。

ステップＳ２１３で計算繰り返し回数が予め定めた最大繰り返し数を超えていないと判定した場合、又は、上述した変数δが予め定めた閾値ε以上であると判定した場合（Ｓ２１３，Ｎ）は、ステップＳ２０５に戻る。また、ステップＳ２１３で計算繰り返し回数が予め定めた最大繰り返し数を超えたと判定した場合、又は、上述した変数δが予め定めた閾値εより小さいと判定した場合（Ｓ２１３，Ｙ）は、本特徴行列推定処理のプログラムの実行を終了する。

このようにして、分布推定装置１により、滞在時間が従う対数正規分布が推定される。

ここで、本実施形態の特殊な例として、データＸの場所分散パラメータσを、ｉに関して共通化したものを考える場合、すなわち目的関数が上記（６）式に従う場合には、特徴行列を更新するための式である上記（１２）乃至（１４）式の代わりに、下記（１５）乃至（１７）式を用いて、後述する手法で特徴行列を更新する。

また、目的関数が上記（７）式に従う場合は、特徴行列を更新するための式である更新式である上記（１２）乃至（１４）式の代わりに、下記（１８）乃至（２０）式を用いて下記のように更新する。

上記（１２）乃至（２０）式の各式において、全ての時間帯ｉ、場所ｊ、ユーザｋについてμ_ｉｊ＝ｌｏｇ（ｘ_ｉｊｋ）が成立する場合、各式の左辺と右辺が一致し、更新の最大変化幅を示す変数δの値が閾値ε以下となるため、更新が停止する。

なお、上記の実施の形態では、時間帯場所行列を表現した行列から分布を推定する例を示しているが、この例に限定されることはない。例えば、地域と年齢とに応じてユーザの年収を表現する行列など、ユーザ、年齢、地域等のように１つ１つにＩＤ番号を付与して識別可能であり行列形式としてデータを表現することが可能な事物であるならば、あらゆるものが本装置による分布推定が可能である。また、本実施形態の図１１に示す分布推定装置の各構成要素の動作をプログラムとして構築し、分布推定装置として利用されるコンピュータにインストールして実行させる、または、ネットワークを介して流通させることが可能である。

また、本実施形態では、図１に示す機能の構成要素の動作をプログラムとして構築し、分布推定装置１として利用されるコンピュータにインストールして実行させるが、これに限らず、ネットワークを介して流通させても良い。

また、構築されたプログラムをハードディスクやフレキシブルディスク、ＣＤ−ＲＯＭ等の可搬記憶媒体に格納し、コンピュータにインストールしたり、配布したりしても良い。

１ 分布推定装置
２ 外部装置
１０ 時間帯場所情報処理部
２０ 特徴行列推定部
２０ａ 反復判定部
３０ 特徴行列処理部
４０ 記憶部
４１ 時間帯場所情報テーブル
４２ 時間帯特徴テーブル
４３ 場所特徴テーブル
４４ 場所分散テーブル
５０ 入出力部

Claims (8)

1. 個々を識別可能な第１の個体群に含まれる個体ｉ（１≦ｉ≦Ｉ，Ｉは１以上の整数）とオブジェクトｊ（１≦ｊ≦Ｊ、Ｊは１以上の整数）との関連度を表すと共に、平均パラメータμ_ｉｊ及び分散パラメータを持つ対数正規分布に従う要素ｘ_ｉｊｋからなるベクトルｘ_ｉｊを要素として持つＩ×Ｊのベクトル値行列Ｘから、
   前記個体ｉが、クラスタｒ（１≦ｒ≦Ｒ、Ｒは１以上の整数）に所属することを表す非負値の要素ａ_ｉｒを持つＩ×Ｒの第１の特徴行列Ａと、前記オブジェクトｊが、前記クラスタｒに所属することを表す非負値の要素ｂ_ｊｒを持つＪ×Ｒの第２の特徴行列Ｂと、を抽出する分布推定装置であって、
   前記ベクトル値行列Ｘの各ベクトルｘ_ｉｊの各要素ｘ_ｉｊｋの対数値、前記第１の特徴行列Ａ、前記第２の特徴行列Ｂ、及び前記分散パラメータを用いて表される目的関数を最適化するように、前記第１の特徴行列Ａ及び前記第２の特徴行列Ｂを推定する特徴行列推定部と、
   予め定めた反復終了条件を満足するまで、前記特徴行列推定部による推定を繰り返す反復判定部と、
   を備えた分布推定装置。
2. 前記ベクトルｘ_ｉｊの要素ｘ_ｉｊｋは、平均パラメータμ_ｉｊ及び個体ｉに依存する分散パラメータσ_ｉを持つ対数正規分布、又は、平均パラメータμ_ｉｊ及びオブジェクトｊに依存する分散パラメータσ_ｊを持つ対数正規分布に従い、
   前記特徴行列推定部は、前記ベクトル値行列Ｘの各ベクトルｘ_ｉｊの各要素ｘ_ｉｊｋの対数値、前記第１の特徴行列Ａ、前記第２の特徴行列Ｂ、及び前記分散パラメータσ_ｉ又は前記分散パラメータσ_ｊを用いて表される目的関数を最適化するように、前記第１の特徴行列Ａ及び前記第２の特徴行列Ｂを推定する
   請求項１記載の分布推定装置。
3. 前記ベクトル値行列Ｘは、時間帯ｉにおける各ユーザの場所ｊの滞在時間を表す要素ｘ_ｉｊｋからなるベクトルｘ_ｉｊを要素として持ち、
   前記第１の特徴行列Ａは、前記時間帯ｉが、前記クラスタｒに所属することを表す非負値の要素ａ_ｉｒを持ち、
   前記第２の特徴行列Ｂは、前記場所ｊが、前記クラスタｒに所属することを表す非負値の要素ｂ_ｊｒを持つ
   請求項１記載の分布推定装置。
4. 前記ベクトルｘ_ｉｊの要素ｘ_ｉｊｋは、平均パラメータμ_ｉｊ及び場所ｊに依存する分散パラメータσ_ｊを持つ対数正規分布に従い、
   前記特徴行列推定部は、下記の式で表される前記目的関数を最小化するように、前記第１の特徴行列Ａ及び前記第２の特徴行列Ｂを推定する
   

   

   請求項３記載の分布推定装置。
5. 前記ベクトルｘ_ｉｊの要素ｘ_ｉｊｋは、平均パラメータμ_ｉｊ、及び時間帯ｉに依存する分散パラメータσ_ｉを持つ対数正規分布に従い、
   前記特徴行列推定部は、下記の式で表される前記目的関数を最小化するように、前記第１の特徴行列Ａ及び前記第２の特徴行列Ｂを推定する
   

   

   請求項３記載の分布推定装置。
6. 個々を識別可能な第１の個体群に含まれる個体ｉ（１≦ｉ≦Ｉ，Ｉは１以上の整数）とオブジェクトｊ（１≦ｊ≦Ｊ、Ｊは１以上の整数）との関連度を表すと共に、平均パラメータμ_ｉｊ及び分散パラメータを持つ対数正規分布に従う要素ｘ_ｉｊｋからなるベクトルｘ_ｉｊを要素として持つＩ×Ｊのベクトル値行列Ｘから、
   前記個体ｉが、クラスタｒ（１≦ｒ≦Ｒ、Ｒは１以上の整数）に所属することを表す非負値の要素ａ_ｉｒを持つＩ×Ｒの第１の特徴行列Ａと、前記オブジェクトｊが、前記クラスタｒに所属することを表す非負値の要素ｂ_ｊｒを持つＪ×Ｒの第２の特徴行列Ｂと、を抽出する分布推定装置における分布推定方法であって、
   特徴行列推定部が、前記ベクトル値行列Ｘの各ベクトルｘ_ｉｊの各要素ｘ_ｉｊｋの対数値、前記第１の特徴行列Ａ、前記第２の特徴行列Ｂ、及び前記分散パラメータを用いて表される目的関数を最適化するように、前記第１の特徴行列Ａ及び前記第２の特徴行列Ｂを推定する特徴行列推定ステップと、
   反復判定部が、予め定めた反復終了条件を満足するまで、前記特徴行列推定ステップによる推定を繰り返す反復判定ステップと、
   を行う分布推定方法。
7. 前記ベクトル値行列Ｘは、時間帯ｉにおける各ユーザの場所ｊの滞在時間を表す要素ｘ_ｉｊｋからなるベクトルｘ_ｉｊを要素として持ち、
   前記第１の特徴行列Ａは、前記時間帯ｉが、前記クラスタｒに所属することを表す非負値の要素ａ_ｉｒを持ち、
   前記第２の特徴行列Ｂは、前記場所ｊが、前記クラスタｒに所属することを表す非負値の要素ｂ_ｊｒを持つ
   請求項６記載の分布推定方法。
8. コンピュータを、請求項１〜５の何れか１項記載の分布推定装置の各部として機能させるためのプログラム。

Images