(57)【要約】【課題】道や場所の潜在コストが動的に変化する場合でも、その変化に追随して動的な潜在コストを推定する。【解決手段】動的潜在コスト推定装置９０は、ユーザ毎に遷移履歴データが入力される毎に、遷移履歴データのユーザｔの潜在的な特徴量を表す特徴ベクトル及び各状態の潜在的な特徴量を表す状態特徴行列を用いて、ユーザｔの各状態の価値関数が表現されると仮定して、ユーザｔの特徴ベクトル及び状態特徴行列の各々の分布を表す分布パラメタを、各ユーザの特徴ベクトルの分布パラメタ、状態特徴行列の分布パラメタ、状態間の基本遷移確率、及び状態間の遷移回数を用いて表される目的関数を最適化して推定するパラメタ推定部４０と、ユーザｔの特徴ベクトル及び状態特徴行列の各々の分布パラメタを用いて、ユーザｔの価値関数を推定し、推定した価値関数、及び状態間の基本遷移確率に基づいて、ユーザｔの報酬関数を計算する報酬関数計算部５０と、を含む。【選択図】図２
2018197653.tif 000002
【特許請求の範囲】
【請求項１】
ユーザ毎又は時間帯毎にユーザが状態間を遷移した履歴を表す遷移履歴データが入力される毎に、前記入力された前記遷移履歴データのユーザｔ又は時間帯ｔの潜在的な特徴量を表す特徴ベクトル、及び各状態の潜在的な特徴量を表す特徴ベクトルを持つ状態特徴行列を用いて、ユーザｔ又は時間帯ｔの各状態の価値関数が表現されると仮定して、ユーザｔ又は時間帯ｔの潜在的な特徴量を表す特徴ベクトルの分布を表す分布パラメタ、及び前記状態特徴行列の分布を表す分布パラメタを、各ユーザ又は各時間帯の潜在的な特徴量を表す特徴ベクトルの分布パラメタ、前記状態特徴行列の分布パラメタ、前記状態間の基本遷移確率、及び前記遷移履歴データから得られる状態間の遷移回数を用いて表される目的関数を最適化するように推定するパラメタ推定部と、
前記パラメタ推定部によって推定された、前記ユーザｔ又は時間帯ｔの潜在的な特徴量を表す特徴ベクトルの分布パラメタ、及び前記状態特徴行列の分布パラメタを用いて、前記ユーザｔ又は時間帯ｔについての各状態の目的地への到達し易さを表す価値関数を推定し、推定した価値関数、及び前記状態間の基本遷移確率に基づいて、前記ユーザｔ又は時間帯ｔについての各状態の目的地らしさを表す報酬関数を計算する報酬関数計算部と、
を含む動的潜在コスト推定装置。
【請求項２】
前記パラメタ推定部は、前記状態間の基本遷移確率の負の対数の集合を表すグラフ重みデータ、前記状態特徴行列の分布パラメタ、及び前記入力された遷移履歴データに基づいて、前記ユーザｔ又は時間帯ｔの潜在的な特徴量を表す特徴ベクトルの分布パラメタを推定し、前記グラフ重みデータ、前記ユーザｔ又は時間帯ｔの潜在的な特徴量を表す特徴ベクトルの分布パラメタ、及び前記状態特徴行列の分布パラメタに基づいて、前記目的関数に含まれる補助変数を推定することを、予め定められた反復終了条件を満たすまで繰り返し、
前記推定した、前記ユーザｔ又は時間帯ｔの潜在的な特徴量を表す特徴ベクトルの分布パラメタ及び前記補助変数に基づいて、前記状態特徴行列の分布パラメタを推定する請求項１に記載の動的潜在コスト推定装置。
【請求項３】
前記目的関数は、前記状態間の基本遷移確率を用いて表される、最適方策により行動を選択する場合の状態間の遷移確率、前記遷移履歴データから得られる状態間の遷移回数、前記ユーザｔ又は時間帯ｔの潜在的な特徴量を表す特徴ベクトルの分布パラメタ、及び前記状態特徴行列の分布パラメタを用いて表され、
前記報酬関数計算部は、前記推定した価値関数、及び前記状態間の基本遷移確率を用いて表される、最適方策により行動を選択する場合の状態間の遷移確率に基づいて、前記ユーザｔ又は時間帯ｔについての各状態の報酬関数を計算する請求項１又は２に記載の動的潜在コスト推定装置。
【請求項４】
パラメタ推定部が、ユーザ毎又は時間帯毎にユーザが状態間を遷移した履歴を表す遷移履歴データが入力される毎に、前記入力された前記遷移履歴データのユーザｔ又は時間帯ｔの潜在的な特徴量を表す特徴ベクトル、及び各状態の潜在的な特徴量を表す特徴ベクトルを持つ状態特徴行列を用いて、ユーザｔ又は時間帯ｔの各状態の価値関数が表現されると仮定して、ユーザｔ又は時間帯ｔの潜在的な特徴量を表す特徴ベクトルの分布を表す分布パラメタ、及び前記状態特徴行列の分布を表す分布パラメタを、各ユーザ又は各時間帯の潜在的な特徴量を表す特徴ベクトルの分布パラメタ、前記状態特徴行列の分布パラメタ、前記状態間の基本遷移確率、及び前記遷移履歴データから得られる状態間の遷移回数を用いて表される目的関数を最適化するように推定するステップと、
報酬関数計算部が、前記パラメタ推定部によって推定された、前記ユーザｔ又は時間帯ｔの潜在的な特徴量を表す特徴ベクトルの分布パラメタ、及び前記状態特徴行列の分布パラメタを用いて、前記ユーザｔ又は時間帯ｔについての各状態の目的地への到達し易さを表す価値関数を推定し、推定した価値関数、及び前記状態間の基本遷移確率に基づいて、前記ユーザｔ又は時間帯ｔについての各状態の目的地らしさを表す報酬関数を計算するステップと、
を含む動的潜在コスト推定方法。
【請求項５】
前記パラメタ推定部が推定するステップは、前記状態間の基本遷移確率の負の対数の集合を表すグラフ重みデータ、前記状態特徴行列の分布パラメタ、及び前記入力された遷移履歴データに基づいて、前記ユーザｔ又は時間帯ｔの潜在的な特徴量を表す特徴ベクトルの分布パラメタを推定し、前記グラフ重みデータ、前記ユーザｔ又は時間帯ｔの潜在的な特徴量を表す特徴ベクトルの分布パラメタ、及び前記状態特徴行列の分布パラメタに基づいて、前記目的関数に含まれる補助変数を推定することを、予め定められた反復終了条件を満たすまで繰り返し、
前記推定した、前記ユーザｔ又は時間帯ｔの潜在的な特徴量を表す特徴ベクトルの分布パラメタ及び前記補助変数に基づいて、前記状態特徴行列の分布パラメタを推定する請求項４に記載の動的潜在コスト推定方法。
【請求項６】
コンピュータを、請求項１〜３のいずれか１項に記載の動的潜在コスト推定装置の各部として機能させるためのプログラム。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、動的潜在コスト推定装置、方法、及びプログラムに係り、特に、人の移動経路における動的な潜在コストを推定する動的潜在コスト推定装置、方法、及びプログラムに関する。
【背景技術】
【０００２】
近年、人の移動経路の選択をモデリングする手法として、逆強化学習（ＩＲＬ：Inverse Reinforcement Learning）と呼ばれる手法が利用されている（例えば、非特許文献１、２を参照）。ＩＲＬでは、人の移動の際に通る道や場所に潜在的なコストが存在し、人は経路全体で支払うコストの和が小さい経路をより高い確率で選択しやすくなると仮定してモデル化される。
【０００３】
上記の仮定に基づいて、ＩＲＬでは、その人の過去の移動履歴から、その人にとっての各道や場所における潜在的なコストを推定する。例えば、非特許文献２に記載されているように、あるタクシー運転手の移動履歴を入力とすることで、その運転手における各道のコストを推定し、運転手の好む、コストの低い道路特性、例えば、整備されている道路を好むか、カーブのゆるやかな道路を好むか等を知ることができる。また、この運転手がある目的地に至るまでに選択しやすい経路がどれであるかを見つけることができる。ほかにもロボット制御、自然言語処理における特筆すべき成功例も存在する（例えば、非特許文献３、４を参照）。
【先行技術文献】
【非特許文献】
【０００４】
【非特許文献１】AndrewY Ng and Stuart Russell. Algorithms for inverse reinforcement learning. In in Proc. 17th International Conf. on Machine Learning, 2000.
【非特許文献２】BrianD Ziebart, AndrewL Maas, JAndrew Bagnell, and AnindK Dey. Maximum entropy inverse reinforcement learning. In AAAI, pages 1433-1438, 2008.
【非特許文献３】Pieter Abbeel and AndrewY Ng. Apprenticeship learning via inverse reinforcement learning. In Proceedings of the twenty-first international conference on Machine learning, page1. ACM, 2004.
【非特許文献４】Gergely Neu and Csaba Szepesv_ari. Training parsers by inverse reinforcement learning. Machine learning, 77(2-3):303-337, 2009.
【非特許文献５】Emanuel Todorov. Linearly-solvable markov decision problems. In Advances in neural information processing systems, pages 1369-1376, 2006.
【非特許文献６】Krishnamurthy Dvijotham and Emanuel Todorov. Inverse optimal control with linearly-solvable mdps. In Proceedings of the 27th International Conference on Machine Learning (ICML-10), pages 335-342, 2010.
【非特許文献７】Masa-Aki Sato. Online model selection based on the variational bayes. Neural Computation, 13(7):1649-1681, 2001.
【非特許文献８】Matthew Hoffman, FrancisR Bach, and DavidM Blei. Online learning for latent dirichlet allocation. In advances in neural information processing systems, pages 856-864, 2010.
【発明の概要】
【発明が解決しようとする課題】
【０００５】
既存のＩＲＬ手法では、各道や場所における潜在的なコストは固定された静的なものとして定義されている。しかしながら、このことは現実に即さない場合がある。例えば、突発的な事故によってある特定の道が通行不能になるという状況や、混雑によってある特定の道の通過に非常に長い時間を要するという状況が存在する。このような状況では、人はその道路を避けた経路を選択する確率が高まると想定できるが、これは事故の時点から一定期間または渋滞発生中の期間、その道路におけるコストが通常時よりも大きい値になることに相当する。
【０００６】
従って、各道や場所における潜在的なコストを固定ではなく、動的に変化すると考えるほうが望ましい場合があるが、既存のＩＲＬ手法では、このような動的な潜在コストを推定することができない。
【０００７】
本発明は、上記の事情に鑑みてなされたもので、道や場所の潜在コストが動的に変化する場合であっても、その変化に追随して動的な潜在コストを推定することができる動的潜在コスト推定装置、方法、及びプログラムを提供することを目的とする。
【課題を解決するための手段】
【０００８】
上記目的を達成するために、第１の発明に係る動的潜在コスト推定装置は、ユーザ毎又は時間帯毎にユーザが状態間を遷移した履歴を表す遷移履歴データが入力される毎に、前記入力された前記遷移履歴データのユーザｔ又は時間帯ｔの潜在的な特徴量を表す特徴ベクトル、及び各状態の潜在的な特徴量を表す特徴ベクトルを持つ状態特徴行列を用いて、ユーザｔ又は時間帯ｔの各状態の価値関数が表現されると仮定して、ユーザｔ又は時間帯ｔの潜在的な特徴量を表す特徴ベクトルの分布を表す分布パラメタ、及び前記状態特徴行列の分布を表す分布パラメタを、各ユーザ又は各時間帯の潜在的な特徴量を表す特徴ベクトルの分布パラメタ、前記状態特徴行列の分布パラメタ、前記状態間の基本遷移確率、及び前記遷移履歴データから得られる状態間の遷移回数を用いて表される目的関数を最適化するように推定するパラメタ推定部と、前記パラメタ推定部によって推定された、前記ユーザｔ又は時間帯ｔの潜在的な特徴量を表す特徴ベクトルの分布パラメタ、及び前記状態特徴行列の分布パラメタを用いて、前記ユーザｔ又は時間帯ｔについての各状態の目的地への到達し易さを表す価値関数を推定し、推定した価値関数、及び前記状態間の基本遷移確率に基づいて、前記ユーザｔ又は時間帯ｔについての各状態の目的地らしさを表す報酬関数を計算する報酬関数計算部と、を含むものである。
【０００９】
第２の発明に係る動的潜在コスト推定装置は、第１の発明において、前記パラメタ推定部が、前記状態間の基本遷移確率の負の対数の集合を表すグラフ重みデータ、前記状態特徴行列の分布パラメタ、及び前記入力された遷移履歴データに基づいて、前記ユーザｔ又は時間帯ｔの潜在的な特徴量を表す特徴ベクトルの分布パラメタを推定し、前記グラフ重みデータ、前記ユーザｔ又は時間帯ｔの潜在的な特徴量を表す特徴ベクトルの分布パラメタ、及び前記状態特徴行列の分布パラメタに基づいて、前記目的関数に含まれる補助変数を推定することを、予め定められた反復終了条件を満たすまで繰り返し、前記推定した、前記ユーザｔ又は時間帯ｔの潜在的な特徴量を表す特徴ベクトルの分布パラメタ及び前記補助変数に基づいて、前記状態特徴行列の分布パラメタを推定するものである。
【００１０】
第３の発明に係る動的潜在コスト推定装置は、第１又は第２の発明において、前記目的関数が、前記状態間の基本遷移確率を用いて表される、最適方策により行動を選択する場合の状態間の遷移確率、前記遷移履歴データから得られる状態間の遷移回数、前記ユーザｔ又は時間帯ｔの潜在的な特徴量を表す特徴ベクトルの分布パラメタ、及び前記状態特徴行列の分布パラメタを用いて表され、前記報酬関数計算部が、前記推定した価値関数、及び前記状態間の基本遷移確率を用いて表される、最適方策により行動を選択する場合の状態間の遷移確率に基づいて、前記ユーザｔ又は時間帯ｔについての各状態の報酬関数を計算するものである。
【００１１】
第４の発明に係る動的潜在コスト推定方法は、パラメタ推定部が、ユーザ毎又は時間帯毎にユーザが状態間を遷移した履歴を表す遷移履歴データが入力される毎に、前記入力された前記遷移履歴データのユーザｔ又は時間帯ｔの潜在的な特徴量を表す特徴ベクトル、及び各状態の潜在的な特徴量を表す特徴ベクトルを持つ状態特徴行列を用いて、ユーザｔ又は時間帯ｔの各状態の価値関数が表現されると仮定して、ユーザｔ又は時間帯ｔの潜在的な特徴量を表す特徴ベクトルの分布を表す分布パラメタ、及び前記状態特徴行列の分布を表す分布パラメタを、各ユーザ又は各時間帯の潜在的な特徴量を表す特徴ベクトルの分布パラメタ、前記状態特徴行列の分布パラメタ、前記状態間の基本遷移確率、及び前記遷移履歴データから得られる状態間の遷移回数を用いて表される目的関数を最適化するように推定するステップと、報酬関数計算部が、前記パラメタ推定部によって推定された、前記ユーザｔ又は時間帯ｔの潜在的な特徴量を表す特徴ベクトルの分布パラメタ、及び前記状態特徴行列の分布パラメタを用いて、前記ユーザｔ又は時間帯ｔについての各状態の目的地への到達し易さを表す価値関数を推定し、推定した価値関数、及び前記状態間の基本遷移確率に基づいて、前記ユーザｔ又は時間帯ｔについての各状態の目的地らしさを表す報酬関数を計算するステップと、を含むものである。
【００１２】
第５の発明に係る動的潜在コスト推定方法は、第４の発明において、前記パラメタ推定部が推定するステップが、前記状態間の基本遷移確率の負の対数の集合を表すグラフ重みデータ、前記状態特徴行列の分布パラメタ、及び前記入力された遷移履歴データに基づいて、前記ユーザｔ又は時間帯ｔの潜在的な特徴量を表す特徴ベクトルの分布パラメタを推定し、前記グラフ重みデータ、前記ユーザｔ又は時間帯ｔの潜在的な特徴量を表す特徴ベクトルの分布パラメタ、及び前記状態特徴行列の分布パラメタに基づいて、前記目的関数に含まれる補助変数を推定することを、予め定められた反復終了条件を満たすまで繰り返し、前記推定した、前記ユーザｔ又は時間帯ｔの潜在的な特徴量を表す特徴ベクトルの分布パラメタ及び前記補助変数に基づいて、前記状態特徴行列の分布パラメタを推定するものである。
【００１３】
第６の発明に係るプログラムは、コンピュータを、第１〜第３のいずれかの発明に係る動的潜在コスト推定装置が備える各部として機能させるためのプログラムである。
【発明の効果】
【００１４】
以上説明したように、本発明に係る動的潜在コスト推定装置、方法、及びプログラムによれば、道や場所の潜在コストが動的に変化する場合であっても、その変化に追随して動的な潜在コストを推定することができる。
これによって、事故や渋滞等が突発的に発生する場合においても、より正確にユーザが選択する移動経路を予測することができるようになる。
【図面の簡単な説明】
【００１５】
【図１】実施形態に係る線形可解マルコフ決定過程（ＬＭＤＰ）における状態集合、基本遷移関数、及び報酬関数の一例を示す模式図である。
【図２】実施形態に係る動的潜在コスト推定装置の構成の一例を示すブロック図である。
【図３】実施形態に係る動的潜在コスト推定プログラムの処理の流れの一例を示すフローチャートである。
【図４】実施形態に係る動的潜在コスト推定プログラムのサブルーチンの処理の流れの一例を示すフローチャートである。
【発明を実施するための形態】
【００１６】
以下、図面を参照して、本発明を実施するための形態の一例について詳細に説明する。
【００１７】
＜本発明の実施形態に係る原理＞
【００１８】
まず、本発明の実施形態における原理を、上記非特許文献５〜８を参考文献として提示しながら説明する。ここでは、本実施形態に係る構成及び作用の説明において用いる記号、アルゴリズムについて記述する。なお、本実施形態においては、ユーザ（人）毎に異なる行動の内部基準（後述する報酬関数に相当する）が存在する例を用いて説明するが、時間帯毎に異なる内部基準が存在する場合であっても同様に適用することができる。
【００１９】
［線形可解マルコフ決定過程（ＬＭＤＰ：Linearly-solvable Markov Decision Process）の定義］
まず、上記非特許文献５を参照して、線形可解マルコフ決定過程（ＬＭＤＰ）について説明する。線形可解マルコフ決定過程（ＬＭＤＰ）では、状態集合
【００２０】
【数１】
2018197653.tif 000003
【００２１】
、基本遷移関数
【００２２】
【数２】
2018197653.tif 000004
【００２３】
、報酬関数
【００２４】
【数３】
2018197653.tif 000005
【００２５】
が定義された環境を考える。本実施形態では、一例として、ユーザ毎に報酬関数が存在する場合について述べるため、上記の定義を、非特許文献５から若干変更する。ユーザ数をＩ、状態数をＪとする。状態集合を
【００２６】
【数４】
2018197653.tif 000006
【００２７】
、基本遷移関数をＪ×Ｊ行列
【００２８】
【数５】
2018197653.tif 000007
【００２９】
、報酬関数をＩ×Ｊ行列
【００３０】
【数６】
2018197653.tif 000008
【００３１】
と定義する。非特許文献５に記載の定義では、上記定義においてユーザ数Ｉ＝１であり、報酬関数がＪ列ベクトルで与えられる場合を考えていることに相当する。上記３つ組
【００３２】
【数７】
2018197653.tif 000009
【００３３】
の具体例を図１に示す。なお、図１は、本実施形態に係る線形可解マルコフ決定過程（ＬＭＤＰ）における状態集合、基本遷移関数、及び報酬関数の一例を示す模式図である。図１中の各グリッドが、状態集合
【００３４】
【数８】
2018197653.tif 000010
【００３５】
中のいずれかの状態ｓｊに対応し、各グリッド間の“基本的な”遷移の確率（基本遷移確率）が基本遷移関数
【００３６】
【数９】
2018197653.tif 000011
【００３７】
で決められている。“基本的な”と注釈をつけたのは、ＬＭＤＰでは行動
【００３８】
【数１０】
2018197653.tif 000012
【００３９】
を実行し、行動コストを支払うことで、この遷移の確率を変更することが許されているからである。具体的には、行動
【００４０】
【数１１】
2018197653.tif 000013
【００４１】
を実行することで、状態ｊから状態ｋへの遷移の確率を
【００４２】
【数１２】
2018197653.tif 000014
【００４３】
から
【００４４】
【数１３】
2018197653.tif 000015
【００４５】
へ変更することができる。
【００４６】
【数１４】
2018197653.tif 000016
【００４７】
ａｊの値が大きければ大きいほど、状態ｊへの遷移の確率を増加できることがわかる。ただし、
【００４８】
【数１５】
2018197653.tif 000017
【００４９】
の値は確率の和が１となることを満たすもの、すなわち
【００５０】
【数１６】
2018197653.tif 000018
【００５１】
を満たすものに限られる。
【００５２】
【数１７】
2018197653.tif 000019
【００５３】
また、上記のように状態ｊで行動
【００５４】
【数１８】
2018197653.tif 000020
【００５５】
を実行することにより発生するコスト
【００５６】
【数１９】
2018197653.tif 000021
【００５７】
は、確率分布の類似度を図る尺度であるＫＬ（カルバックライブラーダイバージェンス）距離を用いて以下（３）式のように定義される。ＫＬ距離は、値が小さいほど類似していることを表し、一致すれば０となる。
【００５８】
【数２０】
2018197653.tif 000022
【００５９】
つまり、上記（３）式においては遷移の確率を
【００６０】
【数２１】
2018197653.tif 000023
【００６１】
から変化させればさせるほど、支払うコストが多くなる。上記の行動コストとは別に、ユーザと状態のみに依存した状態コストを表すのが、報酬関数
【００６２】
【数２２】
2018197653.tif 000024
【００６３】
である。上記図１の例で示すように、この報酬関数の値ｒｉｊは、例えば、“ユーザｉにとっての状態ｊの目的地らしさ”に対応した値などに対応し、値が小さいほど目的地らしいといえる。さらに、各状態ｓｊにおいて実行する行動
【００６４】
【数２３】
2018197653.tif 000025
【００６５】
を定めたものを方策
【００６６】
【数２４】
2018197653.tif 000026
【００６７】
と呼ぶこととする。
【００６８】
方策πにしたがって行動するユーザｉを考えると、上記（１）式にしたがって確率的に状態間を遷移し、遷移する度に行動コストと状態コストとの和
【００６９】
【数２５】
2018197653.tif 000027
【００７０】
を支払う動きがモデリングできる。状態遷移を十分な回数繰り返した際に、遷移の過程で支払ったコスト（行動コストと状態コストとの和）の累積和の期待値が最小になるような方策を最適方策と呼び、
【００７１】
【数２６】
2018197653.tif 000028
【００７２】
と書く。すなわち、最適方策とは、先ほどの例を用いると、できるだけ少ない行動コストで状態コストの低い目的地に早く到達できる方策であり、直感的には、目的地に到達する最適な手順を表しているとみなせる。このように、基本遷移関数
【００７３】
【数２７】
2018197653.tif 000029
【００７４】
と報酬関数
【００７５】
【数２８】
2018197653.tif 000030
【００７６】
が既知の時に、最適方策
【００７７】
【数２９】
2018197653.tif 000031
【００７８】
を求める問題が、ＬＭＤＰにおける順問題である。ユーザｉの最適方策は、以下（４）式の最適性方程式を満たす価値関数
【００７９】
【数３０】
2018197653.tif 000032
2018197653.tif 000033
【００８０】
を求めることで、次の（６）式のように算出できる。なお、上記（４）式から（５）式への式展開は、上記非特許文献５を参照できる。また、価値関数は、各状態の目的地への到達しやすさを表す。価値関数の求め方は後述する。
【００８１】
【数３１】
2018197653.tif 000034
【００８２】
ただし、γは、割引率
【００８３】
【数３２】
2018197653.tif 000035
【００８４】
であり、将来に得る報酬をどれだけ減衰させて評価するかを定めるものである。ここで、上記（６）式により得られた最適方策の値を上記（１）式に代入することで、最適方策により行動を選択する場合の遷移確率は、以下（７）式で与えられることがわかる。
【００８５】
【数３３】
2018197653.tif 000036
【００８６】
したがって、価値関数
【００８７】
【数３４】
2018197653.tif 000037
【００８８】
が計算できれば、ＬＭＤＰの順問題が解けることになる。
【００８９】
［価値関数計算アルゴリズム］
次に、価値関数を計算するアルゴリズムについて説明する。価値関数を計算するアルゴリズムは上記（５）式から導かれる。両辺の符号を変換したのち、指数をとると、以下の（８）式が求まる。
【００９０】
【数３５】
2018197653.tif 000038
【００９１】
ここで
【００９２】
【数３６】
2018197653.tif 000039
【００９３】
と定義すると、上記（８）式は、以下（９）式のようにベクトル表記できる。
【００９４】
【数３７】
2018197653.tif 000040
【００９５】
そして、以下（１０）式に従い、
【００９６】
【数３８】
2018197653.tif 000041
【００９７】
と収束するまで更新することで、上記（９）式を満たす
【００９８】
【数３９】
2018197653.tif 000042
【００９９】
を求めることができる（上記非特許文献５を参照）。
【０１００】
【数４０】
2018197653.tif 000043
【０１０１】
最後に
【０１０２】
【数４１】
2018197653.tif 000044
【０１０３】
の定義より、
【０１０４】
【数４２】
2018197653.tif 000045
【０１０５】
と変換することで、価値関数
【０１０６】
【数４３】
2018197653.tif 000046
【０１０７】
が求まる。各ユーザｉについて以上の手続きを繰り返すことで、各ユーザｉの価値関数ｖの集合を表す価値関数Ｖが求まる。以上がＬＭＤＰの順問題の解法である。なお、ユーザ数Ｉ＝１であれば、上記の手続きは、上記非特許文献５による価値関数計算方法と同等であり、計算方法に本質的な違いは存在しない。
【０１０８】
［逆問題の解法］
次に、上記非特許文献６を参照して、ＬＭＤＰの逆問題を説明する。ＬＭＤＰの逆問題は、人、ロボットの遷移の系列と基本遷移関数とを入力とし、報酬関数を推定する問題、として定式化される。逆問題における解法の鍵は、人、ロボットの遷移が上記（７）式の最適方策による遷移確率にしたがって得られているとみなすことにある。ユーザｉの遷移履歴データ
【０１０９】
【数４４】
2018197653.tif 000047
【０１１０】
を（遷移前の状態、遷移後の状態）という１ステップの遷移の履歴の集合で表現し、
【０１１１】
【数４５】
2018197653.tif 000048
【０１１２】
を状態ｊが遷移前の状態として現れた回数とし、
【０１１３】
【数４６】
2018197653.tif 000049
【０１１４】
を状態ｊが遷移後の状態として現れた回数とする。ユーザｉの価値関数
【０１１５】
【数４７】
2018197653.tif 000050
【０１１６】
を未知のパラメタとみなすと、上記（７）式から、負の対数尤度関数
【０１１７】
【数４８】
2018197653.tif 000051
【０１１８】
を最小化する
【０１１９】
【数４９】
2018197653.tif 000052
【０１２０】
が、入力の遷移の履歴
【０１２１】
【数５０】
2018197653.tif 000053
【０１２２】
を生成する確率の最も高い最尤推定量
【０１２３】
【数５１】
2018197653.tif 000054
【０１２４】
となる。
【０１２５】
【数５２】
2018197653.tif 000055
【０１２６】
上記の尤度関数の最小化には、最急降下法やニュートン法など任意の方法が適用できる。上記（５）式に最尤推定量
【０１２７】
【数５３】
2018197653.tif 000056
【０１２８】
を代入すれば、報酬関数は、次の（１２）式により求めることができる。
【０１２９】
【数５４】
2018197653.tif 000057
【０１３０】
以上がＬＭＤＰの逆問題の解法である。上記のアプローチでは、ユーザｉの推定にはユーザｉの遷移履歴データ
【０１３１】
【数５５】
2018197653.tif 000058
【０１３２】
のみを利用して推定を行っている。また、基本遷移確率
【０１３３】
【数５６】
2018197653.tif 000059
【０１３４】
は、既知であるとした。
【０１３５】
＜本発明の実施形態に係る概要＞
【０１３６】
本実施形態に係る手法は、ユーザ毎に報酬関数が存在する状況において、各ユーザ毎の報酬関数と基本遷移確率とを動的に推定することを可能とする新しい手法である。本実施形態に係るパラメタ推定の定式化の方法は、最尤推定や変分ベイズ推定等の多様な方法が存在する。以下では、一例として、上記非特許文献７、８等に記載されている変分ベイズ推定法に基づくオンラインアルゴリズムについて示す。但し、最尤推定法などに基づくアルゴリズムも、以下に示すアルゴリズム導出の手続きから容易に類推して導出可能である。まず、コストが静的な場合のアルゴリズムを示した後、コストが動的な場合のアルゴリズムを示す。
【０１３７】
変分ベイズ推定法では、推定するパラメタの確率分布を推定する。総ユーザ数をＩ、あるユーザｉの潜在的な特徴量を表す特徴ベクトルを、
【０１３８】
【数５７】
2018197653.tif 000060
【０１３９】
と書く。同様に、総場所数をＪ、ある場所ｊの潜在的な特徴量を表す特徴ベクトルを、
【０１４０】
【数５８】
2018197653.tif 000061
2018197653.tif 000062
【０１４１】
と書く。Ｒは上記特徴ベクトルの次元を表す。Ｒ＜＜Ｉ，Ｊを満たす。
【０１４２】
上記の特徴ベクトルを各行にもつ行列を、
【０１４３】
【数５９】
2018197653.tif 000063
【０１４４】
と書き、それぞれユーザ特徴行列、状態特徴行列と呼ぶ。これらの行列が推定対象となるパラメタＡ、Ｂ（以下、個別にユーザパラメタＡ、場所パラメタＢともいう。）である。
【０１４５】
ユーザｉの場所ｊにおける価値関数の値が、
【０１４６】
【数６０】
2018197653.tif 000064
【０１４７】
として、上記２つの特徴ベクトルの積として表現されると仮定し、あるユーザ特徴行列、状態特徴行列が与えられた場合に、全ユーザの移動履歴を表すデータ
【０１４８】
【数６１】
2018197653.tif 000065
【０１４９】
を生成する確率は、上記（７）式より、
【０１５０】
【数６２】
2018197653.tif 000066
【０１５１】
と書ける。ただし、
【０１５２】
【数６３】
2018197653.tif 000067
【０１５３】
は、ユーザｉが場所ｊから場所ｋへ遷移した回数を表す。
【０１５４】
また、以後の記述で、
【０１５５】
【数６３】
2018197653.tif 000068
【０１５６】
のように、添え字の一部をドット記号で記述されている場合には, その添え字に関しては和をとった値を表すこととする。すなわち、
【０１５７】
【数６４】
2018197653.tif 000069
【０１５８】
である。また、記述の簡便性のため、 基本遷移確率を
【０１５９】
【数６５】
2018197653.tif 000070
【０１６０】
と定義しなおした。また、ユーザパラメタＡ及び場所パラメタＢの各々の事前分布は、正規分布で与えられるとし、以下の式で定義する。
【０１６１】
【数６６】
2018197653.tif 000071
【０１６２】
ただし、αＡ、αＢは、ハイパーパラメタである。これらより、観測データとパラメタの同時確率は、次の式のように書き下せる。
【０１６３】
【数６７】
2018197653.tif 000072
【０１６４】
前述したように、変分ベイズ推定法による学習では、ユーザパラメタＡ及び場所パラメタＢの各々の分布を推定する。以後、ユーザパラメタＡ及び場所パラメタＢの各々の分布を
【０１６５】
【数６８】
2018197653.tif 000073
【０１６６】
とする。このパラメタ分布を
【０１６７】
【数６９】
2018197653.tif 000074
【０１６８】
のように、分解できるという条件のもと、目的関数
【０１６９】
【数７０】
2018197653.tif 000075
【０１７０】
を最小化することで求める方法が変分ベイズ推定法である。なお、関数ｈは、任意の対数尤度関数の下界となる関数を示し、ここでは以下の定義を利用する(この関数が下界であることの証明は省略する)。
【０１７１】
【数７１】
2018197653.tif 000076
【０１７２】
【数７２】
2018197653.tif 000077
【０１７３】
は、補助変数であり、これらの補助変数も最適化の対象とされる。また、
【０１７４】
【数７３】
2018197653.tif 000078
【０１７５】
と定義した。σ（ ）は、シグモイド関数を表す。変分ヘイズ推定法による解析から、
【０１７６】
【数７４】
2018197653.tif 000079
【０１７７】
が最小となるとき、行列Ａ、Ｂの分布は正規分布であり、その平均と分散のパラメタは、以下の式を満たすことが示される。なお、平均がμ、分散がΣとする。
【０１７８】
【数７５】
2018197653.tif 000080
【０１７９】
なお、上記（１６）〜（１９）式において、
【０１８０】
【数７６】
2018197653.tif 000081
【０１８１】
は、
【０１８２】
【数７７】
2018197653.tif 000082
【０１８３】
であり、
【０１８４】
【数７８】
2018197653.tif 000083
【０１８５】
は、
【０１８６】
【数７９】
2018197653.tif 000084
【０１８７】
である。また、上記（１６）、（１８）式において、上付きのＴは転置を示す。
【０１８８】
これらの式に基づいて、「各分布のパラメタを他の分布のパラメタを固定したもとで更新する」ことを何度も繰り返すことで、パラメタ分布を得ることができる。このアルゴリムが静的なコストの場合における変分ベイズアルゴリズムである。なお、補助変数も上記繰り返し手続きの途中に下記の更新を実行することで最適化することができる。
【０１８９】
【数８０】
2018197653.tif 000085
【０１９０】
ただし、
【０１９１】
【数８１】
2018197653.tif 000086
【０１９２】
である。なお、｜Ωｊ｜は、状態ｊから１ステップで移動可能な状態の数を表す。上記のアルゴリズムに基づいて、遷移履歴データが逐次的に到着し、動的にコストが変化する状況のアルゴリズムを導く。
【０１９３】
ここでは特に、動的に場所のコストが変化し、一例として、逐次的に新たなユーザの移動に関する遷移履歴データが到着する状況を考える。すなわち、過去のユーザの移動を表現していた場所パラメタＢ（場所の特徴ベクトル）と、現在の場所における場所パラメタＢは等しくないという状況である。したがって、過去のデータよりもより現在のデータをより重視する、という方針のもとアルゴリズムを設計する。以下にアルゴリズムの手続きを示す。なお、ここでいう新たなユーザとは、総ユーザ数Ｉに含まれない新たなユーザを想定しているが、過去のユーザを対象としても同様である。
【０１９４】
新たなユーザｔの遷移履歴データから、ユーザｔが場所ｊから場所ｋへ遷移した回数を示すデータ
【０１９５】
【数８２】
2018197653.tif 000087
【０１９６】
が得られたとする。そして、これまでに推定されている場所の場所パラメタＢに基づいて、ユーザｔの特徴ベクトルａｔと関連する補助変数
【０１９７】
【数８３】
2018197653.tif 000088
【０１９８】
を推定する。これらの補助変数に基づいて、場所パラメタＢに関する平均及び分散のパラメタの値を、上記式（１８）、（１９）式のうち過去のユーザに関する項を減衰させたうえで、ユーザｔに関する項を加えることで更新する。具体的な平均及び分散のパラメタの更新式を以下に示す。
【０１９９】
【数８４】
2018197653.tif 000089
【０２００】
このような手続きによって、状態特徴行列Ｂに関する平均及び分散のパラメタは、過去のデータを忘却しながら、現在のデータに追随するよう推定されていくようになる。なお、このようにして推定した平均及び分散のパラメタを用いれば、ユーザｔの場所ｊに関する価値関数の推定値は、
【０２０１】
【数８５】
2018197653.tif 000090
【０２０２】
であるから、報酬関数を求める上記（１２）式の記号
【０２０３】
【数８６】
2018197653.tif 000091
【０２０４】
を、記号
【０２０５】
【数８７】
2018197653.tif 000092
【０２０６】
で置き換えて、
【０２０７】
【数８８】
2018197653.tif 000093
【０２０８】
と計算することで、ユーザｔの場所ｊにおける報酬関数の推定値
【０２０９】
【数８９】
2018197653.tif 000094
【０２１０】
が得られる。上記が本実施形態に係るアルゴリズムである。
【０２１１】
＜本発明の実施形態に係る動的潜在コスト推定装置の構成＞
【０２１２】
図２は、本実施形態に係る動的潜在コスト推定装置９０の構成の一例を示すブロック図である。
図２に示すように、本実施形態に係る動的潜在コスト推定装置９０は、ＣＰＵ（Central Processing Unit）、ＲＡＭ(Random Access Memory)、ＲＯＭ(Read Only Memory)、及びＨＤＤ(Hard Disk Drive)等を備えたコンピュータとして構成される。ＲＯＭには、本実施形態に係る動的潜在コスト推定プログラムが記憶されている。なお、動的潜在コスト推定プログラムは、ＨＤＤに記憶されていてもよい。
【０２１３】
上記の動的潜在コスト推定プログラムは、例えば、動的潜在コスト推定装置９０に予めインストールされていてもよい。この動的潜在コスト推定プログラムは、不揮発性の記憶媒体に記憶して、又は、ネットワークを介して配布して、動的潜在コスト推定装置９０に適宜インストールすることで実現してもよい。なお、不揮発性の記憶媒体の例としては、ＣＤ-ＲＯＭ(Compact Disc Read Only Memory)、光磁気ディスク、ＤＶＤ-ＲＯＭ(Digital Versatile Disc Read Only Memory)、フラッシュメモリ、メモリカード等が挙げられる。
【０２１４】
ＣＰＵは、ＲＯＭに記憶されている動的潜在コスト推定プログラムを読み込んで実行することにより、グラフ重み処理部１０、状態特徴行列処理部２０、遷移履歴処理部３０、パラメタ推定部４０、報酬関数計算部５０、及び報酬関数処理部６０として機能する。また、ＨＤＤは、記録部７０として機能する。記録部７０は、グラフ重み記録部７１、状態特徴行列記録部７２、遷移履歴記録部７３、ユーザ特徴行列記録部７４、及び報酬関数記録部７５を含む。
【０２１５】
動的潜在コスト推定装置９０は、入出力部８０を備える。入出力部８０は、外部装置９２から出力された、状態間の基本遷移確率の負の対数として定義される
【０２１６】
【数９０】
2018197653.tif 000095
【０２１７】
の集合を表すグラフ重みデータ
【０２１８】
【数９１】
2018197653.tif 000096
【０２１９】
及びハイパーパラメタαＡ、αＢの入力を受け付ける。
【０２２０】
グラフ重み処理部１０は、入出力部８０により入力を受け付けたグラフ重みデータ
【０２２１】
【数９２】
2018197653.tif 000097
【０２２２】
及びハイパーパラメタαＡ、αＢをグラフ重み記録部７１に格納する。
【０２２３】
また、入出力部８０は、外部装置９２から出力された、状態特徴行列Ｂの分布パラメタ
【０２２４】
【数９３】
2018197653.tif 000098
【０２２５】
の初期値の入力を受け付ける。
【０２２６】
状態特徴行列処理部２０は、入出力部８０により入力を受け付けた状態特徴行列Ｂの分布パラメタ
【０２２７】
【数９４】
2018197653.tif 000099
【０２２８】
の初期値を状態特徴行列記録部７２に格納する。
【０２２９】
また、入出力部８０は、外部装置９２から出力された、推定対象とされる、新たなユーザｔの遷移履歴データ
【０２３０】
【数９５】
2018197653.tif 000100
【０２３１】
の入力を受け付ける。
【０２３２】
遷移履歴処理部３０は、入出力部８０により入力を受け付けたユーザｔの遷移履歴データ
【０２３３】
【数９６】
2018197653.tif 000101
【０２３４】
を遷移履歴記録部７３に格納する。本実施の形態では、新たなユーザｔについての遷移履歴データが到着する毎に、到着した遷移履歴データが遷移履歴記録部７３に格納される。
【０２３５】
パラメタ推定部４０は、グラフ重み記録部７１に格納されたグラフ重みデータ
【０２３６】
【数９７】
2018197653.tif 000102
【０２３７】
及びハイパーパラメタαＡ、αＢ、状態特徴行列記録部７２に格納された状態特徴行列Ｂの分布パラメタ
【０２３８】
【数９８】
2018197653.tif 000103
【０２３９】
の初期値、及び遷移履歴記録部７３に格納されたユーザｔの遷移履歴データ
【０２４０】
【数９９】
2018197653.tif 000104
【０２４１】
を入力とし、上記（１４）式に示す目的関数を最適化するように、ユーザｔの潜在的な特徴量を表す特徴ベクトルａｔ（以下、ユーザ特徴ベクトルａｔ）の分布を表す分布パラメタ
【０２４２】
【数１００】
2018197653.tif 000105
【０２４３】
と、ユーザ特徴ベクトルａｔと関連する補助変数
【０２４４】
【数１０１】
2018197653.tif 000106
【０２４５】
と、状態特徴行列Ｂの分布を表す分布パラメタ
【０２４６】
【数１０２】
2018197653.tif 000107
【０２４７】
と、を推定する。本実施形態では、目的関数を最小化することで、最適化しているが、利用する目的関数の種類によっては、目的関数を最大化することで、最適化してもよい。なお、上記では、ユーザｔのユーザ特徴ベクトルａｔ及び状態特徴行列Ｂを用いて、ユーザｔの各場所の価値関数が表現されると仮定する。また、上記目的関数は、状態間の基本遷移確率を用いて表される、最適方策により行動を選択する場合の状態間の遷移確率（上記（７）式）、上記遷移履歴データから得られる状態間の遷移回数、ユーザｔのユーザ特徴ベクトルａｔの分布パラメタ、及び状態特徴行列Ｂの分布パラメタを用いて表される。本実施形態に係るパラメタ推定部４０は、上記目的関数を最小化することにより、ユーザｔのユーザ特徴ベクトルａｔの分布パラメタ
【０２４８】
【数１０３】
2018197653.tif 000108
【０２４９】
及びユーザ特徴ベクトルａｔと関連する補助変数
【０２５０】
【数１０４】
2018197653.tif 000109
【０２５１】
を推定し、推定した分布パラメタ
【０２５２】
【数１０５】
2018197653.tif 000110
【０２５３】
及び補助変数
【０２５４】
【数１０６】
2018197653.tif 000111
【０２５５】
をユーザ特徴行列記録部７４に格納する。
【０２５６】
具体的には、パラメタ推定部４０は、予め定められた反復終了条件を満たすまで、上記（１４）式で示される目的関数を最小化するように、ユーザ特徴ベクトルａｔの分布パラメタ、及びユーザ特徴ベクトルａｔと関連する補助変数を更新することを繰り返すことにより、ユーザ特徴ベクトルａｔの分布パラメタ、及びユーザ特徴ベクトルａｔと関連する補助変数を推定する。一例として、変分ベイズ推定法により上記（１４）式で示される目的関数が最小となるとき、ユーザ特徴行列Ａ及び状態特徴行列Ｂの各々は正規分布となり、パラメタ推定部４０は、ユーザ特徴ベクトルａｔの分布パラメタ
【０２５７】
【数１０７】
2018197653.tif 000112
【０２５８】
を、上記（１６）、（１７）式に従って更新する。そして、パラメタ推定部４０は、更新前と更新後の平均パラメタ
【０２５９】
【数１０８】
2018197653.tif 000113
【０２６０】
の差の絶対値の最大値
【０２６１】
【数１０９】
2018197653.tif 000114
【０２６２】
が、平均パラメタの更新量の最大変化幅を示す変数δ１より大きい場合に、変数δ１を上記最大値に更新する。但し、更新前の平均パラメタの値を
【０２６３】
【数１１０】
2018197653.tif 000115
【０２６４】
、更新後の値を
【０２６５】
【数１１１】
2018197653.tif 000116
【０２６６】
と記述した。
【０２６７】
また、パラメタ推定部４０は、ユーザ特徴ベクトルａｔと関連する補助変数
【０２６８】
【数１１２】
2018197653.tif 000117
【０２６９】
を、上記（２０）、（２１）式に従って更新する。そして、パラメタ推定部４０は、更新前と更新後の値の差の絶対値の最大値
【０２７０】
【数１１３】
2018197653.tif 000118
【０２７１】
が、値の更新量の最大変化幅を示す変数δ２より大きい場合に、変数δ２を上記最大値に更新する。この変数δ２は、上記変数δ１と同じでもよいし、異なっていてもよい。以下では、説明を簡単にするため、δ１＝δ２として、単に変数δと記述する。但し、更新前の補助変数
【０２７２】
【数１１４】
2018197653.tif 000119
【０２７３】
の要素を
【０２７４】
【数１１５】
2018197653.tif 000120
【０２７５】
、更新後の要素を
【０２７６】
【数１１６】
2018197653.tif 000121
【０２７７】
と記述した。同様に、更新前の補助変数
【０２７８】
【数１１７】
2018197653.tif 000122
【０２７９】
の要素を
【０２８０】
【数１１８】
2018197653.tif 000123
【０２８１】
、更新後の要素を
【０２８２】
【数１１９】
2018197653.tif 000124
【０２８３】
と記述した。
【０２８４】
パラメタ推定部４０は、変数δが予め定められた閾値εより小さい場合に、ユーザ特徴ベクトルａｔの分布パラメタ及び補助変数の更新を終了する。また、パラメタ推定部４０は、ユーザ特徴ベクトルａｔの分布パラメタ及び補助変数の更新を繰り返し、繰り返し回数をカウントする。この場合、パラメタ推定部４０は、カウントした繰り返し回数が予め定められた最大繰り返し回数を超える場合に、ユーザ特徴ベクトルａｔの分布パラメタ及び補助変数の更新を終了してもよい。
【０２８５】
そして、パラメタ推定部４０は、ユーザ特徴行列記録部７４に格納された、更新後のユーザｔについてのユーザ特徴ベクトルａｔの分布パラメタ、及びユーザ特徴ベクトルａｔと関連する補助変数に基づいて、上記（２２）、（２３）式に従って、状態特徴行列Ｂの分布パラメタ
【０２８６】
【数１２０】
2018197653.tif 000125
【０２８７】
を更新することにより、状態特徴行列Ｂの分布パラメタを推定し、推定した状態特徴行列Ｂの分布パラメタ
【０２８８】
【数１２１】
2018197653.tif 000126
【０２８９】
を状態特徴行列記録部７２に格納する。
【０２９０】
報酬関数計算部５０は、状態特徴行列記録部７２に格納された状態特徴行列Ｂの分布パラメタ
【０２９１】
【数１２２】
2018197653.tif 000127
【０２９２】
、及び、ユーザ特徴行列記録部７４に格納されたユーザｔについてのユーザ特徴ベクトルａｔの分布パラメタ
【０２９３】
【数１２３】
2018197653.tif 000128
【０２９４】
に基づいて、ユーザｔについての場所ｊの目的地への到達し易さを表す価値関数
【０２９５】
【数１２４】
2018197653.tif 000129
【０２９６】
を推定し、推定した価値関数、及び状態間の基本遷移確率（グラフ重みデータＷ）を用いて表される、最適方策により行動を選択する場合の状態間の遷移確率に基づいて、上記（２４）式を用いて、ユーザｔについての場所ｊの目的地らしさを表す報酬関数
【０２９７】
【数１２５】
2018197653.tif 000130
【０２９８】
を計算し、報酬関数記録部７５に格納する。
【０２９９】
報酬関数処理部６０は、報酬関数記録部７５に格納された報酬関数
【０３００】
【数１２６】
2018197653.tif 000131
【０３０１】
を、入出力部８０を介して外部装置９２に出力する。
【０３０２】
なお、本実施形態では、上記遷移履歴データがユーザ毎に入力される場合について説明したが、上記遷移履歴データが時間帯毎に入力される場合についても同様に適用することができる。この場合、ユーザｔを、時間帯ｔに置き換えればよい。つまり、パラメタ推定部４０は、時間帯ｔについての特徴ベクトルａｔの分布パラメタ、及び状態特徴行列Ｂの分布パラメタを推定する。報酬関数計算部５０は、パラメタ推定部４０により推定された、時間帯ｔについての特徴ベクトルａｔの分布パラメタ、及び状態特徴行列Ｂの分布パラメタに基づいて、時間帯ｔについての各状態の報酬関数を計算する。
【０３０３】
＜本発明の実施形態に係る動的潜在コスト推定装置の作用＞
【０３０４】
次に、図３を参照して、本実施形態に係る動的潜在コスト推定装置９０の作用について説明する。なお、図３は、本実施形態に係る動的潜在コスト推定プログラムの処理の流れの一例を示すフローチャートである。
【０３０５】
まず、図３のステップ１００では、グラフ重み処理部１０が、外部装置９２から入出力部８０を介してグラフ重みデータＷ及びハイパーパラメタαＡ、αＢの入力を受け付け、グラフ重み記録部７１に格納する。一方、状態特徴行列処理部２０は、外部装置９２から入出力部８０を介して状態特徴行列Ｂの分布パラメタの初期値の入力を受け付け、状態特徴行列記録部７２に格納する。
【０３０６】
ステップ１０２では、遷移履歴処理部３０が、外部装置９２から入出力部８０を介してユーザｔの遷移履歴データＤｔの到着を受け付けたか否かを判定する。ユーザｔの遷移履歴データＤｔの到着を受け付けたと判定した場合（肯定判定の場合）、ユーザｔの遷移履歴データＤｔを遷移履歴記録部７３に格納し、ステップ１０４に移行する。一方、ユーザｔの遷移履歴データＤｔの到着を受け付けていないと判定した場合（否定判定の場合）、ステップ１０２で待機となる。
【０３０７】
ステップ１０４では、パラメタ推定部４０が、グラフ重み記録部７１に格納されたグラフ重みデータＷ、ハイパーパラメタαＡ、αＢ、状態特徴行列記録部７２に格納された状態特徴行列Ｂの分布パラメタの初期値、及び遷移履歴記録部７３に格納されたユーザｔの遷移履歴データＤｔに基づいて、上記（１４）式に示す目的関数を最小化することで、ユーザｔのユーザ特徴ベクトルａｔの分布パラメタ、及び状態特徴行列Ｂの分布パラメタを推定する。具体的には、ユーザｔのユーザ特徴ベクトルａｔの分布パラメタ、及び状態特徴行列Ｂの分布パラメタの更新処理として、図４に示すサブルーチンが実行される。
【０３０８】
図４は、本実施形態に係る動的潜在コスト推定プログラムのサブルーチンの処理の流れの一例を示すフローチャートである。
【０３０９】
まず、図４のステップ２００では、パラメタ推定部４０が、ユーザｔのユーザ特徴ベクトルａｔ、及びユーザ特徴ベクトルａｔと関連する補助変数ηｔ、ξｔを初期化する。
【０３１０】
ステップ２０２では、パラメタ推定部４０が、終了条件に用いる変数δを初期化する。そして、パラメタ推定部４０は、終了条件となる閾値ε、及び最大繰り返し回数を設定する。
【０３１１】
ステップ２０４では、パラメタ推定部４０が、ユーザｔのユーザ特徴ベクトルａｔの分布パラメタを、上記（１６）、（１７）式に従って更新する。このとき、更新前と更新後の平均パラメタの差の絶対値の最大値
【０３１２】
【数１２７】
2018197653.tif 000132
【０３１３】
が変数δより大きい場合、変数δを上記最大値に更新する。
【０３１４】
ステップ２０６では、パラメタ推定部４０が、ユーザｔのユーザ特徴ベクトルａｔと関連する補助変数ηｔ、ξｔを、上記（２０）、（２１）式に従って更新する。このとき、更新前と更新後の値の差の絶対値の最大値
【０３１５】
【数１２８】
2018197653.tif 000133
【０３１６】
が変数δより大きい場合、変数δを上記最大値に更新する。
【０３１７】
ステップ２０８では、パラメタ推定部４０が、ユーザｔのユーザ特徴ベクトルａｔの分布パラメタ、及びユーザ特徴ベクトルａｔと関連する補助変数ηｔ、ξｔの繰り返し回数をカウントアップする。
【０３１８】
ステップ２１０では、パラメタ推定部４０が、終了条件を満たしているか否かを判定する。終了条件を満たしていると判定した場合（肯定判定の場合）、更新したユーザｔのユーザ特徴ベクトルａｔの分布パラメタ、及びユーザ特徴ベクトルａｔと関連する補助変数ηｔ、ξｔをユーザ特徴行列記録部７４に格納し、ステップ２１２に移行する。一方、終了条件を満たしていないと判定した場合（否定判定の場合）、ステップ２１４に移行し、変数δ←０として、ステップ２０４に戻り処理を繰り返す。なお、ここでいう終了条件としては、上述したように、繰り返し回数が最大繰り返し回数を超えた場合、又は、変数δが閾値εより小さい場合に、終了と判定される。
【０３１９】
ステップ２１２では、パラメタ推定部４０が、状態特徴行列Ｂの分布パラメタを、ユーザｔのユーザ特徴ベクトルａｔの分布パラメタ、及びユーザ特徴ベクトルａｔと関連する補助変数ηｔ、ξｔに基づいて、上記（２２）、（２３）式に従って更新し、図３のステップ１０６に移行する。
【０３２０】
次に、図３に戻り、ステップ１０６では、報酬関数計算部５０が、グラフ重みデータＷ、ハイパーパラメタαＡ、αＢ、ユーザｔのユーザ特徴ベクトルａｔの分布パラメタ、及び状態特徴行列Ｂの分布パラメタに基づいて、上記（２４）式を用いて、ユーザｔについての報酬関数を計算し、計算した報酬関数を報酬関数記録部７５に格納する。そして、報酬関数処理部６０は、報酬関数記録部７５に格納された報酬関数を、入出力部８０を介して外部装置９２に出力する。
【０３２１】
なお、上記では、各因子行列Ａ、Ｂの更新式として、一例として変分ベイズ推定法に基づいて導出した式を利用しているが、上記（１３）式で表されるモデルのパラメタを推定する他のいかなる方法を用いても良い。
【０３２２】
本実施形態によれば、道や場所の潜在コストが動的に変化する場合であっても、その変化に追随して動的な潜在コストを推定することができる。
これによって、事故や渋滞等が突発的に発生する場合においても、より正確にユーザが選択する移動経路を予測することができるようになる。
【０３２３】
以上、実施形態として動的潜在コスト推定装置を例示して説明した。実施形態は、コンピュータを、動的潜在コスト推定装置が備える各部として機能させるためのプログラムの形態としてもよい。実施形態は、このプログラムを記憶したコンピュータが読み取り可能な記憶媒体の形態としてもよい。
【０３２４】
その他、上記実施形態で説明した動的潜在コスト推定装置の構成は、一例であり、主旨を逸脱しない範囲内において状況に応じて変更してもよい。
【０３２５】
また、上記実施形態で説明したプログラムの処理の流れも、一例であり、主旨を逸脱しない範囲内において不要なステップを削除したり、新たなステップを追加したり、処理順序を入れ替えたりしてもよい。
【０３２６】
また、上記実施形態では、プログラムを実行することにより、実施形態に係る処理がコンピュータを利用してソフトウェア構成により実現される場合について説明したが、これに限らない。実施形態は、例えば、ハードウェア構成や、ハードウェア構成とソフトウェア構成との組み合わせによって実現してもよい。
【符号の説明】
【０３２７】
１０   グラフ重み処理部
２０   状態特徴行列処理部
３０   遷移履歴処理部
４０   パラメタ推定部
５０   報酬関数計算部
６０   報酬関数処理部
７０   記録部
７１   グラフ重み記録部
７２   状態特徴行列記録部
７３   遷移履歴記録部
７４   ユーザ特徴行列記録部
７５   報酬関数記録部
８０   入出力部
９０   動的潜在コスト推定装置
９２   外部装置
【図１】
2018197653.tif 000134
【図２】
2018197653.tif 000135
【図３】
2018197653.tif 000136
【図４】
2018197653.tif 000137