特許 5963493 コスト推定システム、方法及びプログラム

(19)【発行国】日本国特許庁(JP)

(12)【公報種別】特許公報(B2)

(11)【特許番号】5963493

(24)【登録日】2016年7月8日

(45)【発行日】2016年8月3日

(54)【発明の名称】コスト推定システム、方法及びプログラム

(51)【国際特許分類】

   G06F 17/18 20060101AFI20160721BHJP

【ＦＩ】

   G06F17/18 Z

【請求項の数】18

【全頁数】16

(21)【出願番号】特願2012-70503(P2012-70503)

(22)【出願日】2012年3月27日

(65)【公開番号】特開2013-205848(P2013-205848A)

(43)【公開日】2013年10月7日

【審査請求日】2014年9月3日

【国等の委託研究の成果に係る記載事項】（出願人による申告）平成２１年度、総務省、「自動車二酸化炭素排出量削減のための大規模モビリティ社会シミュレータの研究開発」委託事業、産業技術力強化法第１９条の適用を受ける特許出願

(73)【特許権者】

【識別番号】390009531

【氏名又は名称】インターナショナル・ビジネス・マシーンズ・コーポレーション

【氏名又は名称原語表記】ＩＮＴＥＲＮＡＴＩＯＮＡＬ ＢＵＳＩＮＥＳＳ ＭＡＣＨＩＮＥＳ ＣＯＲＰＯＲＡＴＩＯＮ

(74)【代理人】

【識別番号】100108501

【弁理士】

【氏名又は名称】上野 剛史

(74)【代理人】

【識別番号】100112690

【弁理士】

【氏名又は名称】太佐 種一

(74)【代理人】

【識別番号】100091568

【弁理士】

【氏名又は名称】市位 嘉宏

(72)【発明者】

【氏名】高橋 力矢

【審査官】 田中 幸雄

(56)【参考文献】

【文献】 国際公開第２０１１／０７４３６９（ＷＯ，Ａ１）

(58)【調査した分野】（Int.Cl.，ＤＢ名）

Ｇ０６Ｆ １７／１８

(57)【特許請求の範囲】

【請求項1】

コンピュータの処理によって、グラフ・データのリンクに関連付けられた値の分布から、該リンクにおける該値の確率密度関数を推定する方法であって、
推定する手段が、リンク別に分けたデータから、前記リンクに対応して、前記リンクの値の分布をあらわす基底関数を推定するステップと、
与える手段が、前記グラフ・データのリンクの類似するリンク間で前記基底関数を補間するように、前記基底関数をその信頼性を示すよう決定された重要度スカラーで混合する式に従って、前記リンクに対応する前記確率密度関数を与えるステップとを含み、
ここで前記リンク間の類似は、前記グラフ・データのグラフの負のラプラシアン行列とパラメータとを用いて疎行列のべき乗として与えられるカーネル行列の各成分として求められる、
方法。

【請求項2】

前記推定するステップは、
グラフ・データのリンクに関連付けられた値を集めた後、該集めたデータを分割することにより、複数の分割されたデータの集合を生成するステップと、
前記分割されたデータの集合毎に確率密度関数を推定するステップと、
前記基底関数を前記分割されたデータの集合毎の前記確率密度関数の線形結合としてあらわすための係数を凸クラスタリングの技法により決定するステップと、
前記決定された係数で以って、前記基底関数を前記分割されたデータの集合毎の前記確率密度関数の線形結合としてあらわすステップを含む、請求項１に記載の方法。

【請求項3】

前記分割されたデータの集合毎の確率密度関数が、ガンマ分布に従う確率変数の確率密度関数、または対数正規分布に従う確率変数の確率密度関数であり、前記基底関数が混合ガンマ分布に従う変数の関数または混合対数正規分布に従う変数の関数である、請求項２に記載の方法。

【請求項4】

前記重要度スカラーは、対応する前記リンクに対して観測されたサンプル数が多いほど高い値をもつように決定される、請求項１に記載の方法。

【請求項5】

前記重要度スカラーは、Kullback-Leibler Importance Estimation Procedureで定義される目的関数を前記式に適用することにより決定される、請求項１に記載の方法。

【請求項6】

前記グラフ・データが道路網のデータであり、前記リンクに関連付けられた値が旅行時間である、請求項１に記載の方法。

【請求項7】

グラフ・データのリンクに関連付けられた値の分布から、該リンクにおける該値の確率密度関数を推定するシステムにおいて、
コンピュータを、
リンク別に分けたデータから、前記リンクに対応して、前記リンクの値の分布をあらわす基底関数を推定する手段と、
前記グラフ・データのグラフの負のラプラシアン行列とパラメータとを用いて疎行列のべき乗として与えられるカーネル行列の各成分としてリンク間の類似を求める手段と、
前記グラフ・データのリンクの類似するリンク間で前記基底関数を補間するように、前記基底関数をその信頼性を示すよう決定された重要度スカラーで混合する式に従って、前記リンクに対応する前記確率密度関数を与える手段として機能させるための
プログラム。

【請求項8】

前記推定する手段は、
グラフ・データのリンクに関連付けられた値を集めた後、該集めたデータを分割することにより、複数の分割されたデータの集合を生成し、
前記分割されたデータの集合毎に確率密度関数を推定し、
前記基底関数を前記分割されたデータの集合毎の前記確率密度関数の線形結合としてあらわすための係数を凸クラスタリングの技法により決定し、及び、
前記決定された係数で以って、前記基底関数を前記分割されたデータの集合毎の前記確率密度関数の線形結合としてあらわす、請求項７に記載のプログラム。

【請求項9】

前記分割されたデータの集合毎の確率密度関数が、ガンマ分布に従う確率変数の確率密度関数、または対数正規分布に従う確率変数の確率密度関数であり、前記基底関数が混合ガンマ分布に従う変数の関数または混合対数正規分布に従う変数の関数である、請求項８に記載のプログラム。

【請求項10】

前記重要度スカラーは、対応する前記リンクに対して観測されたサンプル数が多いほど高い値をもつように決定される、請求項７に記載のプログラム。

【請求項11】

前記重要度スカラーは、Kullback-Leibler Importance Estimation Procedureで定義される目的関数を前記式に適用することにより決定される、請求項７に記載のプログラム。

【請求項12】

前記グラフ・データが道路網のデータであり、前記リンクに関連付けられた値が旅行時間である、請求項７に記載のプログラム。

【請求項13】

コンピュータの処理によって、グラフ・データのリンクに関連付けられた値の分布から、該リンクにおける該値の確率密度関数を推定するシステムであって、
リンク別に分けたデータから、前記リンクに対応して、前記リンクの値の分布をあらわす基底関数を推定する手段と、
前記グラフ・データのリンクの類似するリンク間で前記基底関数を補間するように、前記基底関数をその信頼性を示すよう決定された重要度スカラーで混合する式に従って、前記リンクに対応する前記確率密度関数を与える手段を有し、
ここで前記リンク間の類似は、前記グラフ・データのグラフの負のラプラシアン行列とパラメータとを用いて疎行列のべき乗として与えられるカーネル行列の各成分として求められる、
システム。

【請求項14】

前記推定する手段は、
グラフ・データのリンクに関連付けられた値を集めた後、該集めたデータを分割することにより、複数の分割されたデータの集合を生成する手段と、
前記分割されたデータの集合毎に確率密度関数を推定する手段と、
前記基底関数を前記分割されたデータの集合毎の前記確率密度関数の線形結合としてあらわすための係数を凸クラスタリングの技法により決定する手段と、
前記決定された係数で以って、前記基底関数を前記分割されたデータの集合毎の前記確率密度関数の線形結合としてあらわす手段とを含む、請求項１３に記載のシステム。

【請求項15】

前記分割されたデータの集合毎の確率密度関数が、ガンマ分布に従う確率変数の確率密度関数、または対数正規分布に従う確率変数の確率密度関数であり、前記基底関数が混合ガンマ分布に従う変数の関数または混合対数正規分布に従う変数の関数である、請求項１４に記載のシステム。

【請求項16】

前記重要度スカラーは、対応する前記リンクに対して観測されたサンプル数が多いほど高い値をもつように決定される、請求項１３に記載のシステム。

【請求項17】

前記重要度スカラーは、Kullback-Leibler Importance Estimation Procedureで定義される目的関数を前記式に適用することにより決定される、請求項１３に記載のシステム。

【請求項18】

前記グラフ・データが道路網のデータであり、前記リンクに関連付けられた値が旅行時間である、請求項１３に記載のシステム。

【発明の詳細な説明】

【技術分野】

【0001】

この発明は、交通路などのグラフ状の経路において、ある経路に沿ったコストを評価または予測するシステム、方法及びプログラムに関するものである。ここでいうコストとは、典型的には所要時間であるが、消費電力、燃費、二酸化炭素排出量などでもよい。

【背景技術】

【0002】

都市計画において、交通量を予測することは重要な課題である。その課題のうち典型的なこととして、保安、利便性などの目的で、ある経路に沿って所要時間を予測することがある。それに付随して、消費電力、燃費、二酸化炭素排出量などを予測する必要が出てくることがある。このような課題に対して、以下のような従来技術が知られている。

【0003】

特開２００８−７７６３６号公報は、出発地から目的地までの取り得る経路を、辺によって接続されるノードからなる確率的グラフとして表し、各辺は、その辺のコストに関する独立確率分布を有するものとし、目的地に到達するための制約を定義し、そのグラフを、比較的小さな１組の決定論的最小コスト問題に変形し、それを解いて、その制約の中で目的地に到達する確率を最大にする最適経路を求めることを開示する。

【0004】

特開２０１１−１２３００３号公報は、走行する車両の燃費がよい走り方の指標を算出する燃費向上度指標算出装置において、微小時間における車速の変化量を計算する車速演算部と、微小時間における燃料消費量の変化量を計算する燃料消費量演算部と、前記車速の変化量の分散と前記燃料消費量の変化量の分散とを計算する分散分析部とを備え、前記分散に基づいて、燃費がよい走り方の指標を算出することを開示する。

【0005】

しかし、これらの従来技術は、予測するために、辺のコストに関する十分な確率分布の情報を必要とした。ＷＯ２０１１／０７４３６９号公報は、過去の経路の情報が不十分であっても、ある始点と終点の間のコストを予測可能とする技法であって、この技法によれば、訓練データとして、（経路、その経路の費用）の集合が与えられたときに、それを元に、任意リンクeに沿った費用c_eを計算するサブルーチンが用意される。c_eは、f_eと表記するある変数から一意的に計算でき、すなわち、c_eが与えられると、f_eも一意的に決まると仮定する。最初のステップでは、データＤに含まれるすべての始点・終点ペアについて、現在の {f_e}から、最小費用経路が求められる。その結果、データＤは、（経路、その経路の費用）の集合に変換される。そこで、変換されたＤをＤ'と表す。ステップでは、コンピュータの処理により、Ｄ'から、上記サブルーチンを用いて、{ f_e}が再計算される。今回計算された{f_e}が、前回計算された{f_e}と比較され、その変化がある閾値以上であれば、最小費用経路を求めるステップに戻る。そうでなければ、{f_e}が確定する。

【0006】

しかし、ＷＯ２０１１／０７４３６９号公報の技法では、最小二乗法を前提とした回帰を用いているため、リンクに対する旅行時間の分布が正規分布に限られるという制限があった。すると、極めて渋滞した場合の分布の裾野の重さや位相が考慮できない、という問題があった。

【0007】

旅行時間の推定に使用されるその他の技法として、Risi Imre Kondor, John Lafferty, "Diffusion Kernels on Graphs and Other Discrete Spaces", International Conference on Machine Learning, 2002に記述されている、拡散カーネルの計算技法が知られているが、道路ネットワークは一般に巨大であるため、そのような拡散カーネルの計算技法は、計算時間が許容範囲を超えてしまう。

【先行技術文献】

【特許文献】

【0008】

【特許文献1】特開２００８−７７６３６号公報

【特許文献2】特開２０１１−１２３００３号公報

【特許文献3】ＷＯ２０１１／０７４３６９号公報

【非特許文献】

【0009】

【非特許文献1】Risi Imre Kondor, John Lafferty, "Diffusion Kernels on Graphs and Other Discrete Spaces", International Conference on Machine Learning, 2002

【非特許文献2】Rikiya Takahashi, "Sequential minimal optimization in convex clustering repetitions", Statistical Analysis and Data Mining, volume 5, issue 1 (Special Issue: Best Papers of SDM'11), pages 70-89, 2012

【非特許文献3】Sugiyama, M., Suzuki, T., Nakajima, S., Kashima, H., von Bunau, P. & Kawanabe, M.Direct importance estimation for covariate shift adaptation. Annals of the Institute of Statistical Mathematics, vol.60, no.4, pp.699-746, 2008

【発明の概要】

【発明が解決しようとする課題】

【0010】

この発明の目的は、コストの観測値が不足したエッジに対しても、可能な限り高い精度で道路リンク別旅行時間を予測できる技法を提供することにある。

【0011】

この発明の他の目的は、リンクに対する旅行時間の分布として正規分布を前提としなくてもよい、道路リンク別旅行時間予測技法を提供することにある。

【0012】

この発明のさらに他の目的は、多数の観測値が存在しても、現実的な計算時間でリンク別旅行時間予測を行うことができる技法を提供することにある。

【課題を解決するための手段】

【0013】

この発明は、道路などをあらわすグラフの一部のエッジに対して、ある確率変数の観測値が与えられたときに、任意のエッジにおける該当変数の確率分布を、任意の分布に対応可能な確率密度関数によって推定する技法を提供する。

【0014】

その確率密度関数とは、エッジ別に分けたデータから推定した基底関数を、エッジ間の類似度スカラーと各エッジの重要度スカラーをもって混合する、下記のような式である。

【数1】

この式で、eは分布を推定したいエッジ、Y_eはエッジeにおける該当確率変数、e_π[1],...,e_π[m]は観測値の存在するエッジ、φ₁,...,φ_mは各エッジに割り当てられる基底関数、λ₁,...,λ_mは各エッジの非負重要度スカラー、K(e,e_π[i])はエッジeとe_π[i]の非負類似度スカラー、φ₀はエッジに非依存な大域的基底関数、λ₀は大域的基底関数の重要度スカラーである。

【0015】

見てとれるように、この確率密度関数は、類似しているエッジ同士で基底関数を補間する形である。

【0016】

基底関数φ₀,φ₁,...,φ_mは、全区間で積分すると1になるという確率密度関数の条件をみたせば正規分布確率密度関数でなくても任意の関数でよく、その引数は多次元でもよい。

【0017】

この発明に係るシステムは、この式を適用する前に、データから関数φ₀,φ₁,...,φ_mを決定する。次にこの発明に係るシステムは、関数φ₀,φ₁,...,φ_mから、例えば、Rikiya Takahashi, "Sequential minimal optimization in convex clustering repetitions", Statistical Analysis and Data Mining, volume 5, issue 1 (Special Issue: Best Papers of SDM'11), pages 70-89, 2012に記述されている高速凸クラスタリング技法により、λ₀,λ₁,...,λ_mを決定する。

【0018】

このようにして決定された重要度スカラーλ_iは、エッジ毎の推定の信頼性の違いを反映するために必要であり、例えば多くのサンプルが観測されているエッジは、高い重要度をもつ。

【0019】

エッジに非依存な基底関数φ₀と重要度スカラーλ₀は、観測値の存在するどのエッジとも類似性がないエッジeに対して分布を計算するために必要である。

【0020】

本発明の一側面によれば、確率密度関数内のパラメータが有限のデータから高速且つ安定的に最適化される。すなわち、ガンマ分布密度関数などの既知の関数形であるL個の確率密度関数を用意し、上記の各基底関数をそのL個の確率密度関数の線形結合としてあらわす。この線形結合の結合定数も、好適には上記高速凸クラスタリング技法により決定する。

【0021】

類似度スカラーは、好適には、エッジ間隣接グラフのラプラシアン行列Hと所定のパラメータにより与えられる。

【0022】

こうして基底関数と類似度スカラーが得られると、リンク重要度λ₀,λ₁,...,λ_mを決定するために、好適には、Sugiyama, M., Suzuki, T., Nakajima, S., Kashima, H., von Bunau, P. & Kawanabe, M.Direct importance estimation for covariate shift adaptation. Annals of the Institute of Statistical Mathematics, vol.60, no.4, pp.699-746, 2008に記述されているKullback-Leibler Importance Estimation Procedureで定義される目的関数を最大化するように最適化する。この場合にも、好適には、上記高速凸クラスタリング技法が使用される。

【0023】

こうして確率密度関数が構成されると、その確率密度関数を用いて、ある経路に沿ったリンク毎の確率を計算することにより、リンクにおける旅行時間分布などの統計量を計算することができる。

【発明の効果】

【0024】

この発明によれば、観測値が十分あるエッジに対する予測精度を損なうことなく、観測値がないか不足するエッジに対しても十分な分布精度を与えることが可能となる。

【0025】

また、従来のノンパラメトリック密度推定と比較して、大幅に低い計算負荷で確率密度関数が最適化できるので、大量データを効率的に処理することが可能となる。

【図面の簡単な説明】

【0026】

【図1】この発明を実施するためのハードウェア構成の一例のブロック図である。

【図2】この発明を実施するための処理ルーチンを示す機能ブロック図を示す図である。

【図3】グラフ・データとリンクの旅行時間の関係を示す図である。

【図4】基本推定部ルーチンの処理のフローチャートを示す図である。

【図5】交叉検定によるパラメータ選択部ルーチンの処理のフローチャートを示す図である。

【発明を実施するための形態】

【0027】

以下、図面に従って、本発明の実施例を説明する。これらの実施例は、本発明の好適な態様を説明するためのものであり、発明の範囲をここで示すものに限定する意図はないことを理解されたい。また、以下の図を通して、特に断わらない限り、同一符号は、同一の対象を指すものとする。

【0028】

図１を参照すると、本発明の一実施例に係るシステム構成及び処理を実現するためのコンピュータ・ハードウェアのブロック図が示されている。図１において、システム・バス１０２には、ＣＰＵ１０４と、主記憶（ＲＡＭ）１０６と、ハードディスク・ドライブ（ＨＤＤ）１０８と、キーボード１１０と、マウス１１２と、ディスプレイ１１４が接続されている。ＣＰＵ１０４は、好適には、３２ビットまたは６４ビットのアーキテクチャに基づくものであり、例えば、インテル社のＰｅｎｔｉｕｍ（商標）４、インテル社のＣｏｒｅ（商標） ２ ＤＵＯ、ＡＭＤ社のＡｔｈｌｏｎ（商標）などを使用することができる。主記憶１０６は、好適には、４ＧＢ以上の容量、より好ましくは、１６ＧＢ以上の容量をもつものである。

【0029】

ハードディスク・ドライブ１０８には、オペレーティング・システムが、格納されている。オペレーティング・システムは、Ｌｉｎｕｘ（商標）、マイクロソフト社のＷｉｎｄｏｗｓ ７、Ｗｉｎｄｏｗｓ ＸＰ（商標）、Ｗｉｎｄｏｗｓ（商標）２０００、アップルコンピュータのＭａｃ ＯＳ（商標）などの、ＣＰＵ１０４に適合する任意のものでよい。

【0030】

ハードディスク・ドライブ１０８にはさらに、図２及び図３に関連して後述する、道路のグラフ・データ２０２と相対旅行時間全サンプルのデータ２０４が格納されている。道路のグラフ・データ２０２は、基本的には対象とする道路網のノードとリンクの長さ及び方向を含むデータであり、地図情報のデータから変換して作成することができる。相対旅行時間全サンプルのデータ２０４は、好適には、複数の自動車から収集したプローブカーデータから作成されたものであり、自動車が実際にグラフ・データ２０２のリンク（エッジ）を通過した際の、当該リンクの通過に要した相対旅行時間の値を、当該リンクに関連付けて保存している。ここで、相対旅行時間とは、実際に自動車があるリンクの走行に要した絶対時間を分子とし、法定速度により同一リンクを走行した場合の参照時間を分母とした値である。グラフ・データ２０２及び相対旅行時間全サンプルのデータ２０４は、好適には、既知の行列表現をもちいて記述し、ハードディスク・ドライブ１０８に保存されている。

【0031】

ハードディスク・ドライブ１０８にはさらに、図２に関連して後述する、メイン・ルーチン２０６、基本推定部ルーチン２０８、パラメータ群のデータ２１０、及び交叉検定によるパラメータ選択部ルーチン２１２が保存されている。基本推定部ルーチン２０８及び交叉検定によるパラメータ選択部ルーチン２１２は、Ｃ、Ｃ＋＋、Ｃ＃、Ｊａｖａ(R)などの、既知のプログラミング言語でコンパイルして実行可能ファイルとして、ハードティスク・ドライブ１０８に保存され、必要に応じて、ユーザの操作に応じて、オペレーティング・システムの作用により主記憶１０６に呼び出され、実行される。

【0032】

キーボード１１０及びマウス１１２は、オペレーティング・システムが提供するグラフィック・ユーザ・インターフェースに従い、ディスプレイ１１４に表示されたアイコン、タスクバー、ウインドウなどのグラフィック・オブジェクトを操作するために使用される。キーボード１１０及びマウス１１２はまた、ユーザによる始点・終点の入力操作や本発明の実施例に係るプログラムを開始または終了する操作を行うためにも使用される。

【0033】

ディスプレイ１１４は、これには限定されないが、好適には、１０２４×７６８以上の解像度をもち、３２ビットtrue colorのＬＣＤモニタである。ディスプレイ１１４は、例えば、指定したリンクの旅行時間に関する統計値を表示するために使用される。

【0034】

次に、図２を参照して、本発明の機能論理ブロック図の構成を説明する。図２において、メイン・ルーチン２０６は、キーボード１１０及びマウス１１２などの操作に応答して操作画面（図示しない）をディスプレイ１１４を表示するインターフェース機能や、グラフ・データ２０２及び相対旅行時間全サンプルのデータ２０４を読み込む機能や、基本推定部ルーチン２０８及び交叉検定によるパラメータ選択部ルーチン２１２を動作させる制御機能などをもつ。

【0035】

パラメータ群２１０は、基本推定部ルーチン２０８及び交叉検定によるパラメータ選択部ルーチン２１２で使用されるパラメータ群を含む。パラメータ群には、メタパラメータL、メタパラメータr、メタパラメータβ、メタパラメータp、メタパラメータ種類数T、フォールド数Fなどが含まれる。パラメータ群には、オペレータによって予め妥当な値がセットされるが、好適には、交叉検定によるパラメータ選択部ルーチン２１２によって、より最適なパラメータが計算されて置換される。

【0036】

基本推定部ルーチン２０８は、グラフ・データ２０２、相対旅行時間全サンプルのデータ２０４、及びパラメータ群２１０を使用して、適合された確率密度関数２１４を設定する機能をもつ。適合された確率密度関数２１４は好適にはハードディスク・ドライブ１０８に保存され、所望のリンクにおける旅行時間の分布などの統計量を計算したりする目的に使用される。基本推定部ルーチン２０８の処理は、図４のフローチャートを参照して、後で詳細に説明する。

【0037】

交叉検定によるパラメータ選択部ルーチン２１２は、基本推定部ルーチン２０８の機能を呼び出して交叉検定を行うことで、より最適なパラメータを求める処理を実行する。交叉検定によるパラメータ選択部ルーチン２１２の処理は、図５のフローチャートを参照して、後で詳細に説明する。

【0038】

図３は、グラフ・データ２０２と相対旅行時間全サンプルのデータ２０４の例を示すための図である。グラフ・データ２０２は、ノードn₁、・・・、ノードn₁₄で示されるようなノードと、リンクe₁、・・・、リンクe₁₆で示されるようなリンク（エッジ）をもつ。なお、便宜上、図３に示すグラフ・データ２０２のノードとリンクには、参照記号を省略されているものがあるが、実際は、全てのノードと全てのリンクにＩＤが付与されて、識別される。

【0039】

相対旅行時間全サンプルのデータ２０４は、個別のリンクに関連して、例えば、プローブカーデータにより取得された相対旅行時間を記録している。このとき、グラフ・データ２０２のすべてのリンクに相対旅行時間が記録されているとは限らないことを理解されたい。図３では特に、相対旅行時間が記録されているリンクが相対的に太い線で示され、相対旅行時間が記録されていない、すなわちプローブカーデータが取得されていないリンクが相対的に細い線で示されている。

【0040】

ここで、次のような経路で相対旅行時間が記録されたとする。
パス1: ノードn₃から出発して、e₇ → e₈ → e₁₁ → e₁₄という経路でノードn₁₂に達する。
パス2: ノードn₃から出発して、e₇ → e₉ → e₁₀ → e₁₄という経路でノードn₁₂に達する。
パス3: ノードn₃から出発して、e₇ → e₈ → e₁₂ → e₁₅ → e₁₆という経路でノードn₁₄に達する。
パス4: ノードn₉から出発して、e₁₃ → e₁₅ → e₁₆という経路でノードn₁₄に達する。

【0041】

すると例えば、e₇を、パス1とパス2とパス3が通過したことが見て取れる。その際、パス1とパス2とパス3のe₇上の相対旅行時間はそれぞれ異なり得るので、単一のリンクには、複数の 異なる相対旅行時間が関連づけられることになる。相対旅行時間全サンプルのデータ２０４は、このような個々のリンクに関連づけられた複数の相対旅行時間のデータ集合を含む。複数の相対旅行時間のデータ集合が存在すると、個々のリンク毎に相対旅行時間の確率密度関数を計算することができる。このような相対旅行時間の確率密度関数については、図４のフローチャートを参照して、説明する。

【0042】

次に、図４のフローチャートを参照して、基本推定部ルーチン２０８の処理を説明する。基本推定部ルーチン２０８は、参照番号４０２で示すメタパラメータL、参照番号４０２で示すメタパラメータr、相対旅行時間全サンプル２０４、グラフ・データ２０２、及び参照番号４０６で示すメタパラメータβ,pを入力して、処理を行う。

【0043】

メタパラメータLは、確率密度関数ψ_iの個数を示し、メタパラメータrは、推定された分布のバンド幅を調整する正のスカラーであり、メタパラメータβ,pは、カーネル行列Kの計算に使用される。

【0044】

処理を開始すると、ステップ４０８では、基本推定部ルーチン２０８は、グラフ・データ２０２の全リンクの相対旅行時間のサンプルを集めたデータセットYを構成する。すなわち、基本推定部ルーチン２０８は、全てのリンクについて相対旅行時間のサンプルを、リンクの区別なく１つの集合にかき集めて、相対旅行時間でソートしたものがデータセットYである。そして、メタパラメータrに応じた重なりを許しつつ、YからL個のデータセットY₁,...,Y_Lを生成する。すなわち、データセットYがN個のサンプルからなるとすると、ソートされた順にデータセットY₁,...,Y_Lの各々にrN/Lにもっとも近い整数個のサンプルを分ける。rの値は例えば、1、1.5あるいは2などが適当である。交叉検定により決められる値であるが、L = 100の場合はr = 1.5が最適であることが実験的に確認されている。

【0045】

ステップ４１０からステップ４１４までは、l = 1からLまでのループである。そのループ内のステップ４１２では、基本推定部ルーチン２０８は、データセットY_lに対して確率密度関数ψ_lを、最尤推定により、適合する。ここでψ_lは、好適にはガンマ分布か対数正規分布のどちらかである。この推定は、指数分布族に対する推定であるから凸最適化可能なので、パラメータが一意的に決まる。

【0046】

こうして、ステップ４１０からステップ４１４まででl = 1からLまでのループが完了すると、確率密度関数ψ₁,...,ψ_Lが揃う。

【0047】

ステップ４１６では、基本推定部ルーチン２０８は、データセットY[0],Y[1],...,Y[m]を集計する。ここでmは、相対旅行時間サンプルが存在するリンク（図３での太線のリンク）の数である。特にY[0] = Yとし、i >= 1の場合、Y[i]を、サンプルが存在する中で、i番目のリンクに与えられたデータセットとする。

【0048】

ステップ４１８からステップ４２２までは、i = 0からmまでのループである。そのループ内のステップ４２０では、基本推定部ルーチン２０８は、Y[i]の分布に関する確率密度関数または確率質量関数を確率密度関数ψ₁,...,ψ_Lの線形結合であらわしたときの結合定数θ_i1,...,θ_iLを、クラスタリング技法を使用することにより、基底関数φ_iを適合する。式で書くと下記のとおりである。

【数2】

【0049】

ここでのクラスタリングに、混合ガンマ分布、あるいは混合対数正規分布を推定する任意の方法が適用可能である。例えば、ＥＭアルゴリズム、階層ディリクレ過程における変分ベイズ法が適用可能である。この実施例では、より好適なクラスタリング技法として、Rikiya Takahashi, "Sequential minimal optimization in convex clustering repetitions", Statistical Analysis and Data Mining, volume 5, issue 1 (Special Issue: Best Papers of SDM'11), pages 70-89, 2012に記述されている高速凸クラスタリング技法を使用する。この技法は、特願２０１０−２４１０６５号明細書にも記述されている。

【0050】

より具体的に述べると、次のとおりである。
先ず、データセットY[i]に含まれる旅行時間の相対値サンプルをy[i,j] ( j = 1,2,...,n[i])と記述する。ここでn[i]は、データセットY[i]に含まれる旅行時間の相対値サンプルの数である。

【0051】

データセットY[i]について、当該高速凸クラスタリング技法で必要となるカーネルベクトルをk[i,l](l = 1,2,...,L)とする。その成分を書き下すと次のとおりである。

【数3】

このようにカーネルベクトルk[i,l]を与えて高速凸クラスタリング技法を適用すると、特願２０１０−２４１０６５号明細書の記述内容におけるλ₁,...,λ_Lが得られるので、θ_i1 = λ₁,...,θ_iL = λ_Lとみなすことで、結合定数θ_i1,...,θ_iLが決定される。

【0052】

ステップ４１８からステップ４２２までのループがi = 0からmまで完了すると、i = 0からmまで、結合定数θ_i1,...,θ_iLが決定される。

【0053】

一方で、基本推定部ルーチン２０８は、ステップ４２４で、グラフ・データ２０２からグラフの隣接行列Aを次のようにして用意する。
(1) 大きさが(リンクの数)×(リンクの数)であるような行列A = {a_ij}の配列を用意する。
(2) リンク(交差点と交差点を結ぶ道路の最小単位で方向をもつ)e_iとリンクe_jが物理的につながっているときに、a_ij = 0.5 + 0.5*そのリンク同士のコサイン類似度とし、e_jが物理的につながっていないとき、a_ij = 0とする。コサイン類似度とは、内積をベクトルの大きさで割った値である。これを式で記述すると、次のとおりである。

【数4】

但し、この式で、リンクe_iの始点をu_i、終点をv_iとし、リンクe_jの始点をu_j、終点をv_jとする。また、Δ(e_i)は、リンクe_iの方向をもつ２次元ベクトルである。

【0054】

このようにしてグラフの隣接行列Aが求まると、基本推定部ルーチン２０８は、ステップ４２６で、下記の式で定義される行列Dを計算する。

【数5】

ここで、diag()は対角行列をあらわし、|E|はリンクの数をあらわす。

【0055】

そして次の式で、グラフの負のラプラシアン行列Hを計算する。

【数6】

【0056】

こうしてステップ４２６でグラフの負のラプラシアン行列Hが計算されると、基本推定部ルーチン２０８は、ステップ４２８で、行列Hと、パラメータβ,pを用いて、下記の式でカーネル行列Kを計算する。

【数7】

このようにして計算したカーネル行列Kの各成分が、エッジ間類似度スカラーK(e,e_π[i])として使用される。この式から見て取れるように、元のエッジ間隣接グラフの負のラプラシアン行列Hとパラメータβ,pを用いて疎行列のべき乗として与えることで、直接隣接していないエッジ間でも補間が働くようになされる。なお、pの好適な値は8であり、βの好適な値は1〜5などである。pは8〜20程度の値の範囲で選ぶことができるが、pが大きくなるにつれて、行列の非ゼロ成分数が増えて多くのメモリが必要になるため、コンピュータの可用なリソースの範囲でpを選ぶ。βは整数である必要はなく0以上の実数値が指定可能であるが、近似精度確保のためにはβはpより小さい値であることが望ましい。

【0057】

ステップ４３０では、基本推定部ルーチン２０８は、ステップ４１６で用意されたリンク別データセットY[0],Y[1],...,Y[m]と、ステップ４１８〜４２２で用意された基底関数φ₀,φ₁,...,φ_mと、ステップ４２８で用意されたカーネル行列Kを用いて、リンク重要度スカラーλ₀,λ₁,...,λ_mを推定する。これは好適には、Kullback-Leibler Importance Estimation Procedure (KLIEP)で定義される目的関数を最大化するように最適化する。KLIEPは、密度推定を経ることなく２つの密度関数の比を直接推定するためのアルゴリズムであって、Sugiyama, M., Suzuki, T., Nakajima, S., Kashima, H., von Bunau, P. & Kawanabe, M. Direct importance estimation for covariate shift adaptation. Annals of the Institute of Statistical Mathematics, vol.60, no.4, pp.699-746, 2008.に記述されている。上述のRikiya Takahashi, "Sequential minimal optimization in convex clustering repetitions", Statistical Analysis and Data Mining, volume 5, issue 1 (Special Issue: Best Papers of SDM'11), pages 70-89, 2012に記載された高速凸クラスタリング技法が、この計算を現実的な時間で行うために使用することができる。

【0058】

ステップ４３０でリンク重要度スカラーλ₀,λ₁,...,λ_mが推定されると、ステップ４３２で基本推定部ルーチン２０８は、下記の式により、確率密度関数（あるいは確率質量関数）P(y_e = y)を構成し、ステップ４３４で適合された確率密度関数２１４として別のプログラムまたはルーチンに利用可能とする。

【数8】

尚ここで、e_π[1],...,e_π[m]は観測値の存在するエッジである。

【0059】

次に、図５のフローチャートを参照して、交叉検定によるパラメータ選択部ルーチン２１２の処理を説明する。交叉検定によるパラメータ選択部ルーチン２１２は、参照番号５０２で示すメタパラメータ種類数T、参照番号５０４で示すフォールド数F、相対旅行時間全サンプル２０４、グラフ・データ２０２を入力して、基本推定部ルーチン２０８で使用されるメタパラメータL,r,β,pを最適化する処理を行う。

【0060】

交叉検定によるパラメータ選択部ルーチン２１２は、ステップ５０６で、種類数T個分だけ、調べたいパラメータの組(L₁,r₁,β₁,p₁),...,(L_T,r_T,β_T,p_T)を生成する。

【0061】

また、交叉検定によるパラメータ選択部ルーチン２１２は、ステップ５０８で、フォールド数F分だけ、(学習データ,テストデータ)の組(Y_{train 1},Y_{test 1}),...,(Y_{train F},Y_{test F})をランダム分割により生成する。そして、交叉検定によるパラメータ選択部ルーチン２１２の主要な処理は、tに関するループの中で、fに関するループを回すことにより、tとfの組み合わせで、パラメータの組と、(学習データ,テストデータ)の組を基本推定部ルーチン２０８に与えて対数尤度を計算し、その結果に応じて最適パラメータを置き換える処理となる。

【0062】

交叉検定によるパラメータ選択部ルーチン２１２は、外側のtに関するループを、ステップ５１０から開始し、そのループの最初に、ステップ５１２で、S_max = -∞とセットする。

【0063】

交叉検定によるパラメータ選択部ルーチン２１２は、次のステップ５１４で、S_t = 0とセットする。

【0064】

ステップ５１６からは、fに関する内側のループに入る。

【0065】

ステップ５１８では、交叉検定によるパラメータ選択部ルーチン２１２は、tとfの値に応じて、参照番号５２０で示すパラメータの組(L_t,r_t,β_t,p_t)と、参照番号５２２で示す学習データY_{train f}を、グラフ・データ２０２とともに基本推定部ルーチン２０８に提供して、モデルを推定させる。

【0066】

交叉検定によるパラメータ選択部ルーチン２１２は、そうして出力されたモデル(f,t)を用いて、ステップ５２４で、テストデータY_{test f}の対数尤度S_ftを計算し、S_t := S_t + S_ft/Fにより、S_tに値を積み上げる。

【0067】

ステップ５１６からステップ５２６までが1からFまで終わると、交叉検定によるパラメータ選択部ルーチン２１２はステップ５２８で、S_t > S_maxであるかどうか判断し、もしそうなら、ステップ５３０で、交叉検定によるパラメータ選択部ルーチン２１２は、S_max := S_tとセットし、さらに(L^*,r^*,β^*,p^*) = (L_t,r_t,β_t,p_t)して、tに関する外側のループの終端であるステップ５３２に進む。

【0068】

このようにステップ５１０からステップ５３２までをt=1からTまで繰り返すと、結果の(L^*,r^*,β^*,p^*)を最適パラメータ５３４として設定して、交叉検定によるパラメータ選択部ルーチン２１２は処理を終わる。

【0069】

次にメイン・ルーチン２０６は、最適パラメータ(L^*,r^*,β^*,p^*)で基本推定部２０８の推定処理を行い、適合された確率密度関数２１４を更新する。

【0070】

なお、交叉検定によるパラメータ選択部ルーチン２１２の処理は必須ではなく、予め妥当なパラメータが分かっているなら、そのパラメータで基本推定部ルーチン２０８を一度動作させるだけでよい。交叉検定は、妥当なパラメータが不明で、且つデータに基づいて最適化したいときに行う。

【0071】

適合された確率密度関数２１４が得られると、あるリンクに沿った平均の旅行時間、α％タイルの値など、確率密度関数を使用した任意の値を計算することができる。

【0072】

あるリンクに沿った平均の旅行時間は、例えば次のようにして計算される。すなわち、

【数9】

における基底関数φ_i(y) (i=0,1,...,m)について、

【数10】

を計算して、各μ_iを代入した下記の式

【数11】

は、あるリンクにおける旅行時間の相対値に関する平均をあらわす。

【0073】

また、基底関数φ_i(y)を、φ_i(y)をyで積分した累積分布関数で置き換えることにより、あるリンクに沿ったα％タイルの値を計算することができる。

【0074】

尚、上記実施例では、リンクの適合された確率密度関数を重要度スカラーとの混合で与えるための基底関数群φ_iを、別の確率密度関数群ψ_iの線形結合として用意するようにしたが、これには限定されず、場合により直接基底関数群φ_iを適合するようにしてもよい。

【0075】

ただ、基底関数群φ_iを、より少ない数の確率密度関数群ψ_iの線形結合としてあらわすことは、コンピュータの計算量的にも望ましい。

【0076】

m個の基底関数群φ_iを、より少ない数(L < m)の確率密度関数群ψ_iの線形結合としてあらわすことで、ψ_iとしてガンマ関数、対数正規分布などの比較的少数の計算しやすいパラメータの確率密度関数を使用しても、基底関数φ_iは 混合ガンマ関数、あるいは混合対数正規分布となって、基底関数φ_iの表現力が増すのである。

【0077】

しかし、確率密度関数ψ_iとして使える確率密度関数は、ガンマ関数や対数正規分布に限定されず、指数分布族全般が使える他、t分布、パレート分布、あるいは正規分布も使用が可能である。

【0078】

以上、主として交通路のリンクの旅行時間に関する確率密度関数を求める実施例について説明してきたが、消費電力、燃費、二酸化炭素排出量などリンクに関連付けられている任意の量の確率密度あるいは確率質量の推定に適用可能である。また、図４のフローチャートのステップ４２０などにも記述したように、確率密度関数でなく確率質量関数であっても、同様に処理することが可能である。

【0079】

また、発明の実装に使用されるコンピュータ・システムは、特定のハードウェア・アーキテクチャやオペレーティング・システムなどのプラットフォームに限定されず、任意のプラットフォームで実装することが可能である。

【符号の説明】

【0080】

１０２ システム・バス
１０４ ＣＰＵ
１０６ 主記憶
１０８ ハードディスク・ドライブ
１１０ キーボード
１１２ マウス
１１４ ディスプレイ
２０２ グラフ・データ
２０４ 相対旅行時間全サンプル
２０６ メイン・ルーチン
２０８ 基本推定部ルーチン
２１０ パラメータ群
２１２ 交叉検定によるパラメータ選択部
２１４ 適合された確率密度関数

【図1】

【図2】

【図3】

【図4】

【図5】