JP2011243029A - シミュレーション・システム、方法及びプログラム - Google Patents

Abstract

【課題】実測値に対するシミュレーションの近似精度を高める。
【解決手段】実測値に基づき動作の異なる複数のシミュレーションが並列に起動される。時刻Ｔ_０で実測値が得られた時、シミュレーションの計算結果は、時刻Ｔ＞Ｔ_０まで進んでいる。シミュレーションの遡った値である、時刻Ｔ_０での計算結果が、時刻Ｔ_０での実測値と比較される。複数のシミュレーションから、重複を許し実測値から近いほど大きい値をとるような第１の重み付けの離散分布が得られ、離散分布に従い複数のシミュレーションが生成され、シミュレーションの計算値相互間が近いほど大きい値をとる第２の重み付けの離散分布、即ち、シミュレーション間の移行尤度が得られる。この移行尤度に従い、再スタートされる前の複数のシミュレーションから個別に次のシミュレーションが生成され実行される。
【選択図】図３

Description

この発明は、コンピュータ・ベースのシミュレーション技法に関し、より詳しくは、交通シミュレーションなどの短期予測に適合するシミュレーション技法に関する。

近年、都市計画の一環として、いわゆる、高度道路交通システム（Intelligent Transportation System）が、特に１９９０年代後半以降、活発に研究されている。

最近特に要望されてきているのは、自動車流入調整などの道路交通管制業務支援、混雑予測などのカーナビへの情報提供などの、リアルタイムの道路情報処理技術である。このような道路情報処理のベースとなるデータとして、例えば、実際に走行している自動車に取り付けられたセンサから得られたデータである、プローブカーデータがある。プローブカーデータが得られるようになると、それを用いてどのようにしてリアルタイムに、交通量、道路の混み具合などの交通情報を適切に処理するかが課題になりつつある。

シミュレーション一般の背景技術として、実測データに基づきシミュレーションを行う従来技術として、下記のものがある。

特開２０００−１６３４５１号公報は、シミュレーション計算を行う演算装置と、シミュレーション計算範囲の所定位置の観測データを求める観測機器を有する環境予測システムを開示する。このシステムにおいて、演算装置は、シミュレーション計算範囲を二次元座標マトリクスに区切り、それぞれの座標位置に対して、シミュレーション計算のための初期値を設定され、且つ境界条件が入力され、前記設定された初期値及び、入力された境界条件に基づき、所定期間経過後の各座標位置における観察項目の値をシミュレーション計算により求める。観測機器は、所定期間経過後の所定の座標位置に実際の観測データを求め、前記演算装置は、更に所定の座標位置のシミュレーション計算により求めた値と、観測機器による実際の観測データを比較し、差分があれば、前記初期値及び又は、境界条件を補正して再度シミュレーション計算を行うことを開示する。

特開２００２−２８６８６０号公報は、地下水の流れの解析を行う地下水流れ解析処理と、この地下水流れ解析処理にて得られる流れの条件下で移流分散遅延モデルにより流れ方向分散係数と遅延係数を物質輸送パラメータとして地下水中の物質移動解析とを行うシミュレーション装置において、前記移流分散遅延モデルの要素である地下水中の対象物質の地下水流れ方向に対する分散を表す流れ方向分散係数ａLおよび地下水の流れに対する地下水中の対象物質の移動の遅れを表す遅延係数Ｒdのうち、前記遅延係数Ｒdは所望の固定値ｒ0を用いると共に、流れ方向分散係数ａLは複数種を用い、これら複数種それぞれにおける前記物質の濃度分布または特定地点の濃度を経過時間対応に算出する物質輸送解析手段と、所要の経過時間ｍをおいて解析対象地での実測により得られた物質の濃度分布を実測値情報として用い、この実測値情報と前記物質輸送解析手段の算出結果とを比較して、前記算出した結果のうち、前記実測値情報に対する類似度が高い算出結果を選定し、当該算出結果の取得にあたって用いた流れ方向分散係数を以て流れ方向分散係数として決定し、かつ、この算出結果の該当する経過時間ｎと、前記実測での経過時間ｍから遅延係数をr0×ｍ／ｎとして決定する遅延係数決定手段と、を具備することを開示する。

また、交通状況監視に関する従来技術として、下記のものがある。

特開平９−１２８６７７号公報は、走行軌跡を計測する車載機を搭載した車両を走行させ、それによって得た走行軌跡データを用いて交通状況を監視し、また、信号タイミングデータの実績値および、車両走行モデルを用いて走行軌跡を実測の軌跡に近似させ、そのモデルを用いてシミュレーションを行い、信号タイミングを最適化することを開示する。

特開２００４−１１８７３５号公報は、現在から比較的長い将来にわたった交通状況を高い精度で予測するために、交通状況計算部で、道路ネットワークデータ、道路ネットワーク上の信号データ、及びＯＤ交通量を用いた交通流シミュレーションを行なって、交通量の計算値を計算し、将来予測処理部で、過去の交通量の観測値と現在の交通量の観測値とに基づいて、将来の交通量の予測値を予測し、ＯＤ交通量推定部で計算値が現在の交通量の観測値及び予測値に一致するように交通状況計算部で使用するＯＤ交通量を修正することを開示する。

特開２０００−１６３４５１号公報 特開２００２−２８６８６０号公報 特開平９−１２８６７７号公報 特開２００４−１１８７３５号公報

上記従来技術は、観測結果の実測値に基づき、実測値を近似するようにシミュレータを補正しつつ、シミュレーション動作を行うことは開示するが、シミュレータの精度を高めることについては、特に示唆すものではない。

従って、本発明の目的は、実測値に対する近似精度を高めることを可能ならしめるシミュレーション技法を提供することにある。

本発明の更なる目的は、コンピュータ・システムのマルチスレッド環境を利用して、高速且つ、実測値に対して高精度でシミュレーションを実行する技法を提供することにある。

本発明によれば、上記目的は、コンピュータ・システムにより、モンテカルロ・シミュレーションの原理を利用して、実測値に基づき、乱数的に複数のシミュレーションのスレッドを生成し、後の時点で別の実測値が得られたとき、その実測値から近いシミュレーションを複数残して、そこから更に複数のシミュレーションを再スタートさせるようにすることによって達成される。

本発明によれば、最初に、実測値に基づき、少しずつ動作の異なる複数のシミュレーションが並列に起動される。

シミュレーションの計算は、一般的に実世界の動きよりも高速なので、時刻T₀で実測値が得られたとき、シミュレーションの計算の計算結果は、時刻T > T₀まで進んでいる。ここで時刻Tは、コンピュータ上のシミュレーション時刻である。そこで、本発明によれば、既に時刻Tまで進んでいるシミュレーションの遡った値である、時刻T₀での計算結果が、時刻T₀で実測値と比較される。

これを可能とするため、生成された各シミュレーションにおいて、各時刻毎の計算結果が、後に参照可能にログとして保存される。

次に、再スタートされる前の複数のシミュレーションから、好適には重複を許して、実測値から近いほど大きい値をとるような第１の重み付けの離散分布が得られる。

次に、第１の重み付けの離散分布に従い複数のシミュレーションが生成され、そのシミュレーションの計算値相互間の近さを勘案して、近いほど大きい値をとるような第２の重み付けの離散分布、すなわち、シミュレーション間の移行尤度が得られる。

そして、この移行尤度に従い、再スタートされる前の複数のシミュレーションから個別に、次のシミュレーションが生成されて実行される。

本発明の１つの側面によれば、生成される各シミュレーションは、好適には、ある１つのシミュレーション動作関数で、それぞれ異なる引数またはパラメータを付与ことでバリエーションを与えられる。こうしておくと、シミュレーションを再スタートさせるとき、パラメータをコピーすることで、シミュレーションの動作を、再スタートするシミュレーションに継承することができる。

本発明の別の側面によれば、再スタートされるシミュレーションの数は、再スタートされる前のシミュレーションの数と同数で、発生された乱数と、上記第１及び第２の重み付け係数の組み合せにより、重複を許して、可能なシミュレーションから選ばれる。この際、重み付け係数の存在により、実測値に近く且つ移行尤度の高い傾向を示すシミュレーションがより選ばれやすくなる。

本発明の更に別の側面によれば、コンピュータ・システムは、マルチスレッド環境をもち、各々のシミュレーションは、個別のスレッドとして生成される。その際、好適には、コンピュータ・システムは、マルチコアまたはマルチプロセッサのシステムとして実装され、個別のコアまたはプロセッサに、シミュレーションのスレッドを割り当てることにより、シミュレーションがより高速に実行される。

この発明によれば、乱数的に複数生成されたシミュレーションのスレッドから、実測値に近いものがより高い尤度になるように重みをつけてシミュレーションを選び実行されることで、より高い精度のシミュレーションを実行することが可能となる。

また、マルチコアまたはマルチプロセッサのコンピュータ・システム上で、個別のコアまたはプロセッサに、シミュレーションのスレッドを割り当てることにより、シミュレーションがより高速に実行される。このことは、非常に複雑なシミュレーション動作を、実時間の計測よりも高速に実行させることを可能ならしめる。

シミュレーションを実施するための一例のハードウェアのブロック図である。 シミュレーションを実施するための一例の機能ブロック図である。 シミュレーション全体の処理のフローチャートを示す図である。 初期化の処理のフローチャートを示す図である。 シミュレーションの再配置処理のフローチャートを示す図である。 離散分布に従い、シミュレーションの複数のサンプルを生成する処理のフローチャートを示す図である。 シミュレーション処理を模式的に説明する図である。 シミュレーション処理を模式的に説明する図である。

以下、図面を参照して、本発明の一実施例の構成及び処理を説明する。以下の記述では、特に断わらない限り、図面に亘って、同一の要素は同一の符号で参照されるものとする。なお、ここで説明する構成と処理は、一実施例として説明するものであり、本発明の技術的範囲をこの実施例に限定して解釈する意図はないことを理解されたい。

先ず、図１を参照して、本発明を実施するために使用されるコンピュータのハードウェアについて説明する。図１において、ホスト・バス１０２には、複数のＣＰＵ１ １０４ａ、ＣＰＵ２ １０４ｂ、ＣＰＵ３ １０４ｃ、・・・ＣＰＵｎ １０４ｎが接続されている。ホスト・バス１０２にはさらに、ＣＰＵ１ １０４ａ、ＣＰＵ２ １０４ｂ、ＣＰＵ３ １０４ｃ、・・・ＣＰＵｎ １０４ｎの演算処理のためのメイン・メモリ１０６が接続されている。

一方、Ｉ／Ｏバス１０８には、キーボード１１０、マウス１１２、ディスプレイ１１４及びハードティスク・ドライブ１１６が接続されている。Ｉ／Ｏバス１０８は、Ｉ／Ｏブリッジ１１８を介して、ホスト・バス１０２に接続されている。キーボード１１０及びマウス１１２は、オペレータが、コマンドを打ち込んだり、メニューをクリックするなどして、操作するために使用される。ディスプレイ１１４は、必要に応じて、後述する本発明に係るプログラムをＧＵＩで操作するためのメニューを表示するために使用される。

この目的のために使用される好適なコンピュータ・システムのハードウェアとして、ＩＢＭ（Ｒ）Ｓｙｓｔｅｍ Ｘがある。その際、ＣＰＵ１ １０４ａ、ＣＰＵ２ １０４ｂ、ＣＰＵ３ １０４ｃ、・・・ＣＰＵｎ １０４ｎは、例えば、インテル（Ｒ）Ｘｅｏｎ（Ｒ）であり、オペレーティング・システムは、Ｗｉｎｄｏｗｓ（商標）Ｓｅｒｖｅｒ ２００３である。オペレーティング・システムは、ハードティスク・ドライブ１１６に格納され、コンピュータ・システムの起動時に、ハードティスク・ドライブ１１６からメイン・メモリ１０６に読み込まれる。

本発明を実施するためには、必須ではないが、マルチプロセッサ・システムを用いることが望ましい。十分な数のコアまたはプロセッサがあれば、１つのコアまたはプロセッサに１つのシミュレーションのスレッドを割り当てて、高速に実行できるからである。ここでマルチプロセッサ・システムとは、一般に、独立に演算処理し得るプロセッサ機能のコアを複数もつプロセッサを用いるシステムを意図しており、従って、マルチコア・シングルプロセッサ・システム、シングルコア・マルチプロセッサ・システム、及びマルチコア・マルチプロセッサ・システムのどれかでよいことを理解されたい。

なお、本発明を実施するために使用可能なコンピュータ・システムのハードウェアは、ＩＢＭ（Ｒ）Ｓｙｓｔｅｍ Ｘに限定されず、本発明のシミュレーション・プログラムを走らせることができるものであれば、通常のパーソナル・コンピュータなどの任意のコンピュータ・システムを使用することができる。オペレーティング・システムも、Ｗｉｎｄｏｗｓ（Ｒ）に限定されず、Ｌｉｎｕｘ（Ｒ）、Ｍａｃ ＯＳ（Ｒ）など、任意のオペレーティング・システムを使用することができる。さらに、シミュレーション・プログラムを高速で動作させるために、ＰＯＷＥＲ（商標）６ベースで、オペレーティング・システムがＡＩＸ（商標）のＩＢＭ（Ｒ）Ｓｙｓｔｅｍ Ｐなどのコンピュータ・システムを使用してもよい。

ハードティスク・ドライブ１１６にはさらに、図２で後述する本発明に係るシミュレータ２０４、重み計算モジュール２０６、シミュレータ選択モジュール２０８、シミュレーションを実行するためのパラメータ２１０などが格納されており、オペレータのキーボードやマウス操作に応答して、メイン・メモリ１０６にロードされて実行される。

Ｉ／Ｏバス１０８にはさらに、通信インターフェース１２０により、インターネット１２２を介して、観測装置１２４に接続されている。観測装置１２４は、この実施例では、プローブカー・データの観測装置であり、図１のシステムは、通信インターフェース１２０を介して、一定周期、あるいは不定周期にプローブカー・データからデータを取得することができる。尚ここで、プローブカー・データでなく、定点観測データを用いるようにしてもよい。

通信インターフェース１２０は、例えば、イーサネット・カードであり、この場合、図示しないが、特定のウェブ・サーバを介して、観測装置１２４に有線接続される。

あるいは通信インターフェース１２０は、データ通信カードを用いて、サーバにダイアルアップ接続するタイプのものでもよく、さらには公衆無線ＬＡＮによるタイプのものでもよい。本発明において、観測装置１２４に対する接続形態に制限はなく、任意の接続構成を用いることができることを理解されたい。

次に、図２の機能ブロックを参照して、この実施例で使用されるプログラム・モジュールとその機能について説明する。

図２において、制御モジュール２０２は、ハードティスク・ドライブ１１６に保存されており、キーボード１１０やマウス１１２などの操作に応答して、オペレーティング・システムの働きにより、メイン・メモリ１０６にロードされて、実行される。

制御モジュール２０２は、シミュレーション全体の動作を管理する役割を果たすものであり、これには限定されないが、複数のＣＰＵ１ １０４ａ、ＣＰＵ２ １０４ｂ、ＣＰＵ３ １０４ｃ、・・・ＣＰＵｎ １０４ｎのうちの１つに割り当てられる管理スレッドとして実装されてもよい。

本発明によれば、一般的に、多くのシミュレーションを並列に走らせるほど、精度が高まる。並列に走らせることができるシミュレーションの最大数は、システムの可用なリソース、オペレーティング・システムの機能、システムのハードウェア・アーキテクチャ、必要なシミュレーションの動作速度などで決定される。

シミュレータ２０４は、この実施例では、道路上の自動車の動きをシミュレートするものであり、関数Fで以って、下記のように記述される。
x_n+1 = x_n + Δt
y_n+1 = F(x_n,y_n,θ_i)

この式で、n → n+1 を繰り返して、x_nが所定の時刻に到達するまで、シミュレーションのスレッドが走ることになる。
ここで、x_nはシミュレーション時間であり、Δtはクロック1ステップ分の時間刻みである。
また、y_nは、一般的に、複数次元をもつベクトル量であり、シミュレーションの状態量をあらわす。θ_iは、シミュレータ２０４のパラメータであり、このθ_iとして異なる値を与えることにより、シミュレータ２０４は、異なる挙動を示す。θ_iは一般的に、複数の要素からなるベクトル量である。パラメータθ_iの要素には、シミュレータ２０４で使用される乱数の種も含まれてもよい。

シミュレータ２０４は、動作中に、シミュレーション時刻T = x_nと、
計算結果y_nの特徴量S(T)を、例えば、ハードディスク・ドライブ１１６の所定領域に記録する。すると、制御モジュール２０２は、任意の時点で、所与の時刻のシミュレータ２０４の計算結果を取得できる。なお、実測値が得られる時間が予め分かっていれば、計算結果y_nの特徴量S(T)は、その時間付近だけ記録するようにしてもよい。

制御モジュール２０２は、起動した任意のシミュレーションのスレッドを停止、あるいは消滅させることができる。このような機能は、Ｗｉｎ３２ＡＰＩ、Ｌｉｎｕｘ ＡＰＩなど、通常のオペレーティング・システムのＡＰＩとして標準で備わっている機能である。

制御モジュール２０２は、メイン・メモリ１０６などのリソースが許す限り、任意の数のシミュレータ２０４を起動して、個別のスレッドとして並列に走らせることができる。

シミュレーションのためのパラメータ２１０がハードディスク・ドライブ１１６に保存されて、そこには、シミュレータ２０４に引数として与えるためのパラメータθ_i(i = 1,...,k)が含まれる。

制御モジュール２０２が参照しやすくするため、パラメータθ_i(i = 1,...,k)は、システムの起動時に、メイン・メモリ１０６の所定の領域にロードしておいてもよい。

制御モジュール２０２は、シミュレータ２０４を複数個、個別のスレッドとして起動する際、パラメータの引数として、異なるθ_i(i = 1,...,k)を順次与えて、少しずつ挙動の異なるシミュレーションを走らせることができる。

制御モジュール２０２がシミュレータ２０４を複数起動する際、利用可能なプロセッサまたはコアが、起動するシミュレータ２０４の数より多い場合、制御モジュール２０２は、個々のプロセッサまたはコアに、起動された個々のシミュレータ２０４のスレッドが１つずつアサインされる。

重み計算モジュール２０６は、制御モジュール２０２からコールされて、観測値と起動されたシミュレータ２０４の計算結果に基づき、後述するシミュレーションの再配置のために、重みの係数を計算する機能をもつ。

シミュレータ選択モジュール２０８は、制御モジュール２０２からコールされて、重み計算モジュール２０６から得られた重み係数と、起動されたシミュレータ２０４の計算結果に基づき、重複を許して、ｋ個のシミュレータ２０４を生成する機能をもつ。

重み計算モジュール２０６とシミュレータ選択モジュール２０８のより詳しい機能は、後で、図３以下のフローチャートを参照して説明する。

表示モジュール２１２は、制御モジュール２０２により制御されて、シミュレーションの途中経過を、ディスプレイ１１４に表示する機能を有する。

データ観測モジュール２１４は、制御モジュール２０２により制御されて、観測装置１２４から、観測データと、その観測時刻を取得し、メイン・メモリ１０６の所定領域に配置する機能をもつ。

尚、図２に示すモジュールは、使用されているオペレーティング・システム及びシステム環境に適合する限り、Ｃ、Ｃ＋＋、Ｃ＃、Ｊａｖａ(R)などの既知の任意のプログラミング言語により作成することができる。

次に、図３〜図６のフローチャートを参照して、この実施例のシミュレーション実行処理を説明する。図３のステップ３０２で、制御モジュール２０２は、シミュレーション開始に当たっての初期化を行う。

図４のフローチャートは、初期化の処理をより詳細に示す。図４において、制御モジュール２０２は、ステップ４０２で、データ観測モジュール２１４の機能を利用して、実状態s(0)を観測する。これは、時間T = 0での観測値に対応する。s(0)は、シミュレータ２０４の関数Fの計算結果に対応する次元と状態量をもつ。s(0)は、観測値自体であってもよいし、観測値に基づき計算された特徴量であってもよい。この実施例では、s(0)は、観測値に基づき計算された特徴量であるとする。

ステップ４０４では、制御モジュール２０２は、パラメータθ_i(i = 1,...,k)を、メイン・メモリ１０６または、ハードディスク・ドライブ１１６の所定の領域から読み込む。

図３に戻って、制御モジュール２０２は、観測した初期値と、パラメータθ_i(i = 1,...,k)により、下記のようにk個のシミュレーションを並列に走らせる。
ここで、y⁽ⁱ⁾ ₀ = s(0) (i = 1,...,k)とする。
x⁽¹⁾ _n+1 = x⁽¹⁾ _n + Δt
y⁽¹⁾ _n+1 = F(x⁽¹⁾ _n,y⁽¹⁾ _n,θ₁)
---------------------------------------------------------
x⁽²⁾ _n+1 = x⁽²⁾ _n + Δt
y⁽²⁾ _n+1 = F(x⁽²⁾ _n,y⁽²⁾ _n,θ₂)
---------------------------------------------------------
x⁽³⁾ _n+1 = x⁽³⁾ _n + Δt
y⁽³⁾ _n+1 = F(x⁽³⁾ _n,y⁽³⁾ _n,θ₃)
---------------------------------------------------------
・・・・・・・・・・・・・・・・・・・・・・・・・・
・・・・・・・・・・・・・・・・・・・・・・・・・・
---------------------------------------------------------
x^(k) _n+1 = x^(k) _n + Δt
y^(k) _n+1 = F(x^(k) _n,y^(k) _n,θ_k)

ステップ３０６では、制御モジュール２０２は、所定時間（例えば、T₀)が経過したかどうか、判断する。もしそうでなければ、ステップ３０８で、k個のシミュレーションは、そのまま並列実行される。

ステップ３０６で、所定時間が経過したと判断されると、ステップ３１０で、シミュレーションの再配置が行われ、ステップ３１２では、再配置された結果のk個のシミュレーションが並列実行される。尚、ステップ３１０の再配置処理は、図５のフローチャートを参照して、後でより詳細に説明する。

ステップ３０６の判断は、例えば、カウンタとしての所定の変数に値をセットして、デクリメントし、その値がゼロに達したことに応答して行う。シミュレーションの再配置が行われた後は、カウンタとしての所定の変数に値をセットして、デクリメントする処理が再開される。こうして、シミュレーションの再配置は、所定時間毎に行われることになる。

ステップ３０６の判断の段階で、所定時間が経過しても、何かの事情で、データ観測モジュールで実測値が得られていない、という場合がありえる。その場合は、所定時間が経過してもステップ３０６の判断を直ちに真にせず、後で実測値が更新されていたことに応答して、ステップ３０６の判断を真にしてもよい。

ステップ３１４では、制御モジュール２０２は、終了条件に達したかどうかを判断する。ここで終了条件とは例えば、シミュレーション開始後、所定の時間が経過したかどうか、あるいは、観測結果が所定の値になったかどうか、などである。

さらに、観測結果がある程度定期的に得られるなら、ステップ３０６の判断を所定時間の経過とせず、観測結果が入来したタイミングとすることもできる。

次に、図５のフローチャートを参照して、ステップ３１０の再配置処理を、より詳細に説明する。

ステップ５０２では、制御モジュール２０２は、実状態sを観測した時間T₀と、現在のシミュレーション時間Tを読み取る。実状態sを観測した時間T₀は、データ観測モジュール２１４が、観測値とその観測時間を、メモリ上の所定の記憶領域に配置するので、読み取ることができる。また、シミュレーション時間Tは、システムのクロックから読み取ることができる。このとき、シミュレーションは、実世界よりも高速なので、T > T₀である。

ステップ５０４では、制御モジュール２０２が、重み計算モジュール２０６を呼び出して、重み値w⁽¹⁾(i)を、i = 1,...,kに亘って、下記の式により、計算する。

すなわち、k個のシミュレーションの時間T₀における値をそれぞれ、
S1(T₀),S2(T₀),...,Sk(T₀)とする。ここで、Si(T₀)というのは、上記の式で、y⁽ⁱ⁾ _nと同じものである。T₀がステップ数nに対応する。
また、時間T₀における実測値の特徴量をs(T₀)とする。すなわち、シミュレーションの計算結果は、特徴量に対応している。
そこで、ある定義された距離関数d1を以って、
d1(Si(T₀),s(T₀))が大きいほど、w⁽¹⁾(i)が小さくなるようにw⁽¹⁾(i)が選ばれる。ここで、d1は、これには限定されないが、例えば、Si(T₀),s(T₀)をベクトルとしたとき、ユークリッド距離でよい。あるいは、用途に応じて、マンハッタン距離なども使用することができる。
また、これには限定されないが、例えば、 w⁽¹⁾(i) 〜 exp(-d1(Si(T₀),s(T₀)))のように選んでよい。w⁽¹⁾(i) は、実測値とシミュレーションの計算値の単調減少関数になっていればよい。
こうしておいて、下記のように、w⁽¹⁾(i)をスケーリングする。

重み値w⁽¹⁾(i)は、実測値と近い状態値をもつシミュレーションが選ばれやすくなるように分布される。

ステップ５０６では、制御モジュール２０２が、P(i) = w⁽¹⁾(i) i = 1, ...,kの離散分布に従って、シミュレータ選択モジュール２０８を呼び出し、k個のシミュレーション・サンプル{j}を得る。尚、シミュレータ選択モジュール２０８の詳細な処理は、図６のフローチャートを参照して、後で説明する。

ステップ５０８では、制御モジュール２０２が、ある定義された距離関数d2を以って、ステップ５０６で生成されたk個のシミュレーション・サンプル{j}の各要素jにおいて、
d2((Si(T),θ_i),(Sj(T),θ_j))が大きいほど、w⁽²⁾ _j(i)が小さくなるようにw⁽²⁾ _j(i)を選ぶ。ここで、d2は、これには限定されないが、例えば、(Si(T),θ_i),(Sj(T),θ_j)をベクトルとしたとき、ユークリッド距離でよい。d2の式は、d1と同様でもよく、異なっていてもよい。さらに、(Si(T),θ_i)という記法は、ベクトルSi(T)と、ベクトルθ_iを並べたベクトルを意味する。前述のように、Si(T)は特徴量なので、パラメータθ_iも一種の特徴量と解釈してもよい。
また、これには限定されないが、例えば、 w⁽²⁾ _j(i) 〜 exp(-d2((Si(T),θ_i),(Sj(T),θ_j))))のように選んでよい。w⁽²⁾ _j(i)も、実測値とシミュレーションの計算値の単調減少関数になっていればよい。
こうしておいて、下記のように、w⁽²⁾ _j(i)をスケーリングする。

重み値w⁽²⁾ _j(i)は、シミュレーションjからシミュレーションiへの移行尤度を示す。

ステップ５１０では、制御モジュール２０２が、各jについて、P_j(i) = w⁽²⁾ _j(i) i = 1, ...,kの離散分布に従って、シミュレータ選択モジュール２０８を呼び出し、k個のシミュレーション・サンプル{i*}を得る。すなわち、このステップでは、各jについて、
P_j(i) = w⁽²⁾ _j(i) i = 1, ...,kの離散分布に従って、その都度シミュレータ選択モジュール２０８を呼び出し、１つずつシミュレーション・サンプルが選ばれ、j = 1,..,kまで進んだとき、k個のシミュレーション・サンプル{i*}が選ばれている、という次第である。

ステップ５１２では、制御モジュール２０２が、ステップ５１０で生成されたk個のシミュレーション・サンプル{i*}を新たなシミュレータ集合とするように、{i*}のシミュレーション時刻Tでのシミュレーション状態S^*i(T)と、パラメータθ^* _iを、元のシミュレータ集合に下記のように順次コピーする。
for (i = 1;i <= k; i++ ) {
Si(T) = S^*i(T);
θ_i = θ^* _i;
}
ここでパラメータθ^* _iとあるのは、実際は、θ_i ( i = 1,...,k)から選ばれたものであることに留意されたい。

次に図６のフローチャートを参照して、シミュレータ選択モジュール２０８の処理を説明する。この処理は、離散分布P(i) = w(i) i = 1,...,kと、選択するシミュレータの数Nを引数とする。

ステップ６０２では、シミュレータ選択モジュール２０８は、ローカル変数m = 1とおく。

ステップ６０４では、シミュレータ選択モジュール２０８は、既知のアルゴリズムを用いて、範囲(0,1)から一様乱数を生成し、それをxとおく。

ステップ６０４では、シミュレータ選択モジュール２０８は、下記の式をみたすjを選択する。

但し、

とする。

これによってjが決まるので、ステップ６０８で、m番目のシミュレータとしてシミュレータjが選択される。これは言い換えると、m番目のシミュレータは、パラメータとして、θ_jを使用することを意味する。

ステップ６１０では、mが1増分され、ステップ６１２でm > Nかどうかが判断される。そうでなければ、処理はステップ６０４に戻り、そうであるなら、所定の数のシミュレータが選ばれたので、処理を終って、シミュレータ選択モジュール２０８は、もとの処理に戻る。

尚、xの乱数性により、ステップ６０８で、同一のシミュレータが重複して選ばれることがあることを理解されたい。

図７は、この実施例のシミュレータの動作を模式的に示す図である。ここで、縦軸は、シミュレーションの状態値、横軸は時間である。また、説明の便宜上、並列に走るシミュレータは３個だけだとする。

この図は、実測値がt = T₀で得られたとき、シミュレーション時刻は、T > T₀まで進んでいることを示す。ここで、図３のステップ３０６により、シミュレーションの再配置のステップ３１０が呼び出される。

すると、本発明によれば、t = T₀での実測値とシミュレーションの計算値の距離に従い、ステップ５０４で、重みw⁽¹⁾(i)が計算される。

ステップ５０６で、シミュレータ選択モジュール２０８が呼び出され、重みw⁽¹⁾(i)に従い、k個のシミュレーション・サンプルが選ばれる。

そのようにして選ばれたシミュレーション・サンプルに基づき、ステップ５０８において、重みw⁽²⁾ _j(i)が計算され、それに基づき、各j毎にシミュレータ選択モジュール２０８が呼ばれて、改めて結果的にk個のシミュレーション・サンプルが選ばれる。

そのk個のシミュレーション・サンプルの状態値とパラメータが、ステップ５１２で、元のシミュレータ集合の引数としてコピーされ、ステップ３１２で、その結果再配置されたシミュレーションが実行される。

図８は、その結果選ばれたシミュレータS'1, S'2及びS'3がシミュレーション時刻Tから始まって、先に進んでいる様子を示す。

このようなシミュレーションの挙動は、表示モジュール２１２の機能により、ディスプレイ１１４に表示される。

このようにシミュレーション・スレッドが複数あると、それらの状態値を複数のシミュレーションで平均した結果を、最終的にシミュレーションされた結果とみなすことができる。

上記実施例では、第１の重み係数としてのw⁽¹⁾(i)と、第２の重み係数としてのw⁽²⁾ _j(i)の両方を用いて、シミュレーションの離散分布を生成するようにした。この実施例は、精度の点で好ましいが、第１の重み係数のみを用いてシミュレーションの離散分布を生成しても、シミュレーションを実現することができることを理解されたい。このとき、生成したj番目のシミュレーションがもとのシミュレーションのi番目に対応するなら、生成したj番目のシミュレーションに、もとのi番目のシミュレーションの計算結果とパラメータがコピーされる。

以上、本発明を、プローブカーデータまたは定点観測データの測定結果に基づく、交通シミュレーションを例にとって説明したが、本発明は、これには限定されず、金融工学のリスク・シミュレーション、商品価格や株価の変動シミュレーション、天気予報など、実測値に基づき離散分布に従う複数のシミュレーションを確率的に発生されることが有効であるさまざまな分野に適用可能であることを理解されたい。

また、例示されているハードウェア構成は一例であって、本発明は、マルチタスクあるいははマルチスレッドをサポートする環境であるなら、特定のハードウェア、オペレーティングシステムに限定されるものではない。

１０４ａ・・・１０４ｎ ＣＰＵ
１０６ メイン・メモリ
１０８ バス
１１０ キーボード
１１２ マウス
１１４ ディスプレイ
１１６ ハードティスク・ドライブ
１２４ 観測装置
２０２ 制御モジュール
２０４ シミュレータ
２０６ 計算モジュール
２０８ シミュレータ選択モジュール
２１０ パラメータ
２１２ 表示モジュール
２１４ データ観測モジュール

Claims (18)

1. コンピュータの処理によって、シミュレーションを実行する方法であって、
   前記コンピュータの処理によって、所定の観測時点から出発して、実際の事象の実測値を観測して観測値を得るステップと、
   前記コンピュータの処理によって、前記所定の観測時点から出発して、動作が、異なるパラメータによって特徴づけられる、複数のシミュレーション処理を並列に走らせるステップと、
   前記コンピュータの処理によって、前記所定の観測時点後の所定のタイミングで、そのタイミングまでに得られた前記実測値と、前記各シミュレーション処理の計算結果の間の距離に従う離散分布を生成するステップと、
   前記コンピュータの処理によって、前記離散分布を用いて、複数のシミュレーション処理を生成するステップを有する、
   シミュレーション方法。
2. 前記コンピュータの処理によって、前記離散分布に従い複数のシミュレーション処理を生成した際の、任意の前記シミュレーション処理と該生成されたシミュレーション処理の間の移行尤度を生成するステップをさらに有し、前記複数のシミュレーション処理を生成するステップは、該移行尤度を利用してシミュレーション処理を生成する、請求項１に記載の方法。
3. 前記所定のタイミングが、前記所定の観測時点から所定時間経過後である、請求項１に記載の方法。
4. 前記所定のタイミングが、前記実測値が得られたタイミングである、請求項１に記載の方法。
5. 前記移行尤度に従い、複数のシミュレーション処理を生成するステップが、それまでに実行されていたシミュレーション処理から、計算結果と前記パラメータを受け取るステップを有する、請求項１に記載の方法。
6. 前記シミュレーションが交通シミュレーションであり、前記実測値が、プローブカーデータまたは定点観測データである、請求項１に記載の方法。
7. コンピュータの処理によって、シミュレーションを実行するプログラムであって、
   前記コンピュータに、
   所定の観測時点から出発して、実際の事象の実測値を観測して観測値を得るステップと、
   前記所定の観測時点から出発して、動作が、異なるパラメータによって特徴づけられる、複数のシミュレーション処理を並列に走らせるステップと、
   前記所定の観測時点後の所定のタイミングで、そのタイミングまでに得られた前記実測値と、前記各シミュレーション処理の計算結果の間の距離に従う離散分布を生成するステップと、
   前記離散分布を用いて、複数のシミュレーション処理を生成するステップを実行させる、
   シミュレーション・プログラム。
8. 前記離散分布に従い複数のシミュレーション処理を生成した際の、任意の前記シミュレーション処理と該生成されたシミュレーション処理の間の移行尤度を生成するステップをさらに有し、前記複数のシミュレーション処理を生成するステップは、該移行尤度を利用してシミュレーション処理を生成する、請求項７に記載のプログラム。
9. 前記所定のタイミングが、前記所定の観測時点から所定時間経過後である、請求項７に記載のプログラム。
10. 前記所定のタイミングが、前記実測値が得られたタイミングである、請求項７に記載のプログラム。
11. 前記移行尤度に従い、複数のシミュレーション処理を生成するステップが、それまでに実行されていたシミュレーション処理から、計算結果と前記パラメータを受け取るステップを有する、請求項７に記載のプログラム。
12. 前記シミュレーションが交通シミュレーションであり、前記実測値が、プローブカーデータまたは定点観測データである、請求項７に記載のプログラム。
13. コンピュータの処理によって、シミュレーションを実行するシステムであって、
    所定の観測時点から出発して、実際の事象の実測値を観測して観測値を得る手段と、
    前記所定の観測時点から出発して、動作が、異なるパラメータによって特徴づけられる、複数のシミュレーション処理を並列に走らせる手段と、
    前記所定の観測時点後の所定のタイミングで、そのタイミングまでに得られた前記実測値と、前記各シミュレーション処理の計算結果の間の距離に従う離散分布を生成する手段と、
    前記離散分布を用いて、複数のシミュレーション処理を生成する手段とを有する、
    シミュレーション・システム。
14. 前記離散分布に従い複数のシミュレーション処理を生成した際の、任意の前記シミュレーション処理と該生成されたシミュレーション処理の間の移行尤度を生成する手段をさらに有し、前記複数のシミュレーション処理を生成する手段が、該移行尤度に基づき複数のシミュレーション処理を生成する、請求項１３に記載のシステム。
15. 前記所定のタイミングが、前記所定の観測時点から所定時間経過後である、請求項１３に記載のシステム。
16. 前記所定のタイミングが、前記実測値が得られたタイミングである、請求項１３に記載のシステム。
17. 前記移行尤度に従い、複数のシミュレーション処理を生成する手段が、それまでに実行されていたシミュレーション処理から、計算結果と前記パラメータを受け取る処理を行う、請求項１３に記載のシステム。
18. 前記シミュレーションが交通シミュレーションであり、前記実測値が、プローブカーデータまたは定点観測データである、請求項１３に記載のシステム。

Images