JP7298684B2 - 混雑度推定装置、混雑度推定方法及びプログラム - Google Patents

Description

本発明は、混雑度推定装置、混雑度推定方法及びプログラムに関する。

テーマパーク問題と呼ばれる問題が従来から知られている。テーマパーク問題とは、テーマパークの混雑状況をマルチエージェントシミュレーション（ＭＡＳ：Multi-Agent Simulation）で再現することにより、来園者のアトラクション選択行動の分析や混雑状況の推定、混雑を緩和するための制御策の検討・評価等を行うための問題である。このようなテーマパーク問題を対象として、これまでに様々な従来技術が提案されている。

例えば、特許文献１には、テーマパークにおけるアトラクションの待ち時間を予測する技術が記載されている。また、例えば、特許文献２には、テーマパーク等の混雑を緩和するための最適な制御策を探索する技術が記載されている。

特開２０１８－０７３３６１号公報 特開２０１８－１４７０８７号公報

ところで、テーマパークで来園者がアトラクションを選択する際には、来園者間で選好に個人差があることが一般的である。すなわち、アトラクションを選択する際には、各来園者は自身の好みに応じてアトラクションを選択するのが一般的である。しかしながら、従来技術では来園者間の選好の個人差が考慮されていなかった。このため、アトラクションの待ち時間等の混雑度を推定する場合、その精度が高くないことがあった。

本発明の実施の形態は、上記の点に鑑みてなされたもので、混雑度を高い精度で推定することを目的とする。

上記目的を達成するため、本発明の実施の形態におけるパラメータ推定装置は、総数Ｎの選択主体ｎ（ｎ＝１，・・・，Ｎ）の各々が総数Ｍの選択対象ｍ（ｍ＝１，・・・，Ｍ）のいずれかを選択するものとして、各選択主体ｎが許容する各選択主体ｍの混雑度の限界を表す許容限界α_ｎ，ｍを入力して、前記選択主体ｎの許容限界α_ｎ，ｍの前記選択対象ｍに関する最大値を表す忍耐力φ_ｎと、前記選択主体ｎが前記選択対象ｍを選択する際の選好を表す選好ベクトルψ_ｎとを算出し、算出した前記忍耐力φ_ｎと前記選好ベクトルψ_ｎとを用いて、総数Ｉ（＞Ｎ）の選択主体ｉ（ｉ＝１，・・・，Ｉ）の各々の許容限界α_ｉ，ｍを得るためのモデルのパラメータを推定するパラメータ推定手段、を有することを特徴とする。

また、本発明の実施の形態における混雑度推定装置は、総数Ｎの選択主体ｎ（ｎ＝１，・・・，Ｎ）の各々が総数Ｍの選択対象ｍ（ｍ＝１，・・・，Ｍ）のいずれかを選択するものとして、各選択主体ｎが許容する各選択主体ｍの混雑度の限界を表す許容限界α_ｎ，ｍを入力して、前記選択主体ｎの許容限界α_ｎ，ｍの前記選択対象ｍに関する最大値を表す忍耐力φ_ｎと、前記選択主体ｎが前記選択対象ｍを選択する際の選好を表す選好ベクトルψ_ｎとを算出し、算出した前記忍耐力φ_ｎと前記選好ベクトルψ_ｎとを用いて、総数Ｉ（＞Ｎ）の選択主体ｉ（ｉ＝１，・・・，Ｉ）の各々の許容限界α_ｉ，ｍを得るためのモデルのパラメータを推定するパラメータ推定手段と、推定した前記パラメータを用いた前記モデルにより、前記選択主体ｉの各々の許容限界α_ｉ，ｍを算出する算出手段と、前記許容限界α_ｉ，ｍとシミュレーション条件とを入力して、各時刻ｔ（ｔ＝１，・・・・，Ｔ）において前記選択主体ｉが前記選択対象ｍのいずれかを選択することをシミュレーションすることで、前記各時刻ｔにおける前記選択対象ｍの各々の混雑度を推定する混雑度推定手段と、を有することを特徴とする。

混雑度を高い精度で推定することができる。

本発明の実施の形態における混雑度推定装置の全体構成の一例を示す図である。 許容限界モデルのパラメータを推定する処理の一例を示すフローチャートである。 パラメータ推定用の許容限界データの一例を説明するための図である。 忍耐力φ_ｎ及びアトラクション選好ψ_ｎ，ｍの一例を説明するための図である。 シミュレーション用の来園者データを作成する処理の一例を示すフローチャートである。 許容限界ベクトルα_ｉの生成の一例を説明するための図である。 シミュレーション用の来園者データの一例を説明するための図である。 入園時刻Ｉ_ｉ及び退園時刻Ｏ_ｉの設定の一例を説明するための図である。 予定数ｋ_ｉの設定の一例を説明するための図である。 アトラクションの配置の一例を説明するための図である。 シミュレーション用の移動時間データの一例を説明するための図である。 シミュレーション用のアトラクションデータの一例を説明するための図である。 来園者の状態遷移の一例を説明するための図である。 シミュレーションにより混雑度を推定する処理の一例を示すフローチャートである。 シミュレーション結果の一例を説明するための図である。 本発明の実施の形態における混雑度推定装置のハードウェア構成の一例を示す図である。

以下、本発明の実施の形態について説明する。本発明の実施の形態では、複数のアトラクションが配置されたテーマパークにおいて、各来園者のアトラクションに対する選好を考慮して、当該テーマパークの各アトラクションの混雑度（例えば、各アトラクションの待ち時間）をシミュレーションにより推定する混雑度推定装置１０について説明する。ここで、テーマパークとは、何等かのテーマに基づいて施設の一部又は全体を演出した観光施設のことあり、具体例としては遊園地等が挙げられる。なお、テーマパークはレジャーランド等と称されることもある。

ただし、テーマパークにおける各アトラクションの混雑度をシミュレーションにより推定することは一例であって、本発明の実施の形態は、選択主体（例：来園者）が選択可能な複数の対象（例：アトラクション）が配置されている場合に、各対象の混雑度をシミュレーションにより推定することに同様に適用可能である。例えば、来場者が選択可能な複数のイベントブースが配置されたイベント会場において、各イベントブースの混雑度をシミュレーションにより推定することにも同様に適用可能である。

＜理論的構成＞
以降では、本発明の実施の形態の理論的な構成について説明する。

≪許容限界モデルのパラメータ推定≫
アトラクション数をＭ、各アトラクションのインデックスをｍ（ｍ＝１，・・・，Ｍ）として、インデックスｍのアトラクションを「アトラクションｍ」と表す。同様に、各来園者のインデックスをｎとして、インデックスｎの来園者を「来園者ｎ」と表す。また、アトラクションｍに対する来園者ｎの待ち時間の許容限界（つまり、アトラクションｍに対して来園者ｎが許容可能な待ち時間を表すスカラー値）をα_ｎ，ｍとする。

このとき、本発明の実施の形態では、例えばアンケート等によりＮ人の来園者から予め取得しておいた許容限界α_ｎ，ｍ（ｎ＝１，・・・，Ｎ， ｍ＝１，・・・，Ｍ）を用いて、以下の式（１）及び式（２）に示すモデル（このモデルを「許容限界モデル」とも表す。）のパラメータを推定する。

ここで、φ_ｎは来園者ｎの待ち時間の許容限界の最大値を表し、以降では「忍耐力」とも表す。α_ｎは来園者ｎの許容限界α_ｎ，ｍ（ｍ＝１，・・・，Ｍ）をｍ番目の要素とするＭ次元ベクトルであり、以降では「許容限界ベクトル」とも表す。ψ_ｎは、来園者ｎのアトラクションｍに対する相対的な選好を表すスカラー値ψ_ｎ，ｍをｍ番目の要素とするＭ次元ベクトルである。以降では、ψ_ｎ，ｍを「アトラクション選好」、ψ_ｎを「アトラクション選好ベクトル」とも表す。

本発明の実施の形態では、忍耐力φ_ｎが対数正規分布に従うものとして、つまり、ｌｏｇ（φ_ｎ）～Ｎ（μ，σ^２）であるものとして、パラメータμ及びσ^２を最尤法等により推定する。なお、μ及びσ^２はそれぞれ正規分布の平均及び分散である。

また、本発明の実施の形態では、アトラクション選好ベクトルψ_ｎがディリクレ分布Ｄｉｒ（β）に従うものとして、つまり、ψ_ｎ～Ｄｉｒ（β）であるものとして、パラメータβを最尤法により推定する。なお、βはディリクレ分布のパラメータであり、Ｍ次元ベクトル（β_１，・・・，β_Ｍ）で表される。

≪シミュレーション用の来園者データの作成≫
各アトラクションの混雑度をシミュレーションする際に用いる来園者数をＩ、各来園者のインデックスをｉとして、インデックスｉの来園者を「来園者ｉ」と表す。なお、シミュレーションに用いる来園者数Ｉは、許容限界モデルのパラメータ推定に用いた来園者数Ｎよりも非常に大きいものとする。

このとき、本発明の実施の形態では、上記で推定したパラメータμ、σ^２及びβを用いて忍耐力φ_ｉ（ｉ＝１，・・・，Ｉ）とアトラクション選好ベクトルψ_ｉ（ｉ＝１，・・・，Ｉ）とを生成した上で、以下の式（３）により許容限界ベクトルα_ｉ（ｉ＝１，・・・，Ｉ）を生成する。

そして、これらの許容限界ベクトルα_ｉを用いて、シミュレーション用の来園者データを作成する。なお、シミュレーション用の来園者データには、アトラクションｍに対する来園者ｉの許容限界α_ｉ，ｍと、来園者ｉの入園時刻Ｉ_ｉ及び退園時刻Ｏ_ｉと、来園者ｉがアトラクションを体験する予定数ｋ_ｉ（つまり、来園者ｉが何個のアトラクションの体験を予定しているかを示す個数）とが含まれる。

≪シミュレーション≫
シミュレーション用の来園者データと、シミュレーション用のアトラクションデータと、シミュレーション用の移動時間データとを用いて、各アトラクションｍのシミュレーション時刻ｔ毎の待ち時間（混雑度の一例）を推定する。なお、シミュレーション用のアトラクションデータとは、例えば、シミュレーションにおけるアトラクションの処理能力（つまり、単位時間あたりにアトラクションを体験可能な人数）を表すデータである。また、シミュレーション用の移動時間データとは、例えば、シミュレーションにおけるアトラクション間（及びテーマパークの出入口とアトラクションとの間）の移動に要する時間を表すデータである。

このとき、本発明の実施の形態では、来園者ｉのシミュレーション時刻ｔにおけるアトラクションｍの選択確率θ_{ｉ，ｍ，ｔ}を以下の式（４）に示すモデル（このモデルを「多項線形モデル」とも表す。）より算出する。

ここで、Ａ_{ｉ，ｍ，ｔ}＝ｍａｘ（０，α_ｉ，ｍ－Ｗ_ｍ，ｔ）であり、Ｗ_ｍ，ｔは時刻ｔにおけるアトラクションｍの待ち時間である。なお、α_ｉ，ｍは許容限界ベクトルα_ｉのｍ番目の要素（つまり、来園者ｉのアトラクションｍに対する許容限界）である。

これにより、シミュレーションの結果として、各来園者ｉのアトラクション選好ψ_ｉ，ｍを考慮した各アトラクションｍの待ち時間（混雑度の一例）が推定される。

なお、上記の式（４）に示す多項線形モデルは、従来から知られている線形を非負に拡張して、論理的整合性の条件を満たすようにしたモデルのことである。また、論理的整合性の条件とは、有限個の選択肢の中から或る選択肢を選択する場合に、何れかの選択肢を必ず選択するという条件の下では全ての選択肢の選択確率の和が１となり、かつ、全ての選択肢の選択確率が非負となる、という制約条件のことである。

＜混雑度推定装置１０の全体構成＞
次に、本発明の実施の形態における混雑度推定装置１０の全体構成について、図１を参照しながら説明する。図１は、本発明の実施の形態における混雑度推定装置１０の全体構成の一例を示す図である。

図１に示すように、本発明の実施の形態における混雑度推定装置１０は、パラメータ推定部１０１と、来園者データ作成部１０２と、シミュレーション部１０３と、記憶部１０４とを有する。なお、パラメータ推定部１０１と、来園者データ作成部１０２と、シミュレーション部１０３とは、例えば、混雑度推定装置１０にインストールされた１以上のプログラムがプロセッサ等に実行させる処理により実現される。また、記憶部１０４は、例えば、混雑度推定装置１０の補助記憶装置や記録媒体等の任意の記憶装置を用いて実現可能である。

記憶部１０４は、各種データを記憶する。記憶部１０４に記憶されるデータとしては、例えば、後述するパラメータ推定用の許容限界データ、許容限界モデルのパラメータμ、σ^２及びβ、シミュレーション用の来園者データ、シミュレーション用のアトラクションデータ、シミュレーション用の移動時間データ、シミュレーションに用いる来園者数Ｉ等が挙げられる。また、記憶部１０４には、各アトラクションｍのシミュレーション時刻ｔにおける待ち時間Ｗ_ｍ，ｔも記憶される。

パラメータ推定部１０１は、パラメータ推定用の許容限界データを入力として、上記の式（１）及び式（２）に示す許容限界モデルのパラメータμ、σ^２及びβを推定する。ここで、パラメータ推定用の許容限界データとは、例えばアンケート等によりＮ人の来園者から予め取得しておいた許容限界α_ｎ，ｍ（ｎ＝１，・・・，Ｎ， ｍ＝１，・・・，Ｍ）を示すデータである。

来園者データ作成部１０２は、パラメータ推定部１０１により推定されたパラメータμ、σ^２及びβと、シミュレーションに用いる来園者数Ｉとを入力として、許容限界ベクトルα_ｉ（ｉ＝１，・・・，Ｉ）を生成する。そして、来園者データ作成部１０２は、これらの許容限界ベクトルα_ｉ（ｉ＝１，・・・，Ｉ）を用いて、シミュレーション用の来園者データを作成する。

シミュレーション部１０３は、シミュレーション用の来園者データと、シミュレーション用のアトラクションデータと、シミュレーション用の移動時間データとを入力として、各アトラクションｍのシミュレーション時刻ｔ毎の待ち時間Ｗ_ｍ，ｔを推定する。ここで、本発明の実施の形態では、シミュレーション時刻ｔは非負整数値の「分」であるものとして、シミュレーション開始からの経過時間［分］を表すものとする。具体的には、シミュレーションの終了時刻をＴ［分］として、ｔ＝０，１，２，・・・，Ｔ［分］と表されるものとする。ただし、これに限られず、シミュレーション時刻ｔは、任意の時間毎（例えば、３０分毎や１時間毎等）のインデックスを表すものであってもよい。

なお、図１に示す混雑度推定装置１０の全体構成は一例であって、他の構成であってもよい。例えば、パラメータ推定部１０１、来園者データ作成部１０２及びシミュレーション部１０３のそれぞれを異なる装置が有していてもよい。すなわち、図１に示す混雑度推定装置１０は、例えば、パラメータ推定部１０１を有するパラメータ推定装置と、来園者データ作成部１０２を有する来園者データ作成装置と、シミュレーション部１０３を有するシミュレーション装置とに分けることができてもよい。

＜許容限界モデルのパラメータ推定処理＞
以降では、許容限界モデルのパラメータを推定する処理について、図２を参照しながら説明する。図２は、許容限界モデルのパラメータを推定する処理の一例を示すフローチャートである。

ステップＳ１０１：パラメータ推定部１０１は、パラメータ推定用の許容限界データを入力する。上述したように、パラメータ推定用の許容限界データとは、許容限界α_ｎ，ｍ（ｎ＝１，・・・，Ｎ， ｍ＝１，・・・，Ｍ）を示すデータ（例えば、Ｎ×Ｍの行列で表され、（ｎ，ｍ）成分をα_ｎ，ｍとするデータ）である。なお、パラメータ推定部１０１は、例えば、記憶部１０４に記憶されているパラメータ推定用の許容限界データを入力してもよいし、通信ネットワークを介して接続される他の装置から送信されたパラメータ推定用の許容限界データを入力してもよい。

ここで、一例として、Ｎ＝５、Ｍ＝３である場合におけるパラメータ推定用の許容限界データの具体例を図３に示す。図３では、許容限界α_１，１＝２４．０、許容限界α_１，２＝９７．８、許容限界α_１，３＝４９．７である場合を示している。同様に、許容限界α_２，１＝５８．７、許容限界α_２，２＝２８．７、許容限界α_２，３＝２７．６である場合を示している。ｎ＝３，４，５の場合も図３に示す通りである。

ステップＳ１０２：パラメータ推定部１０１は、パラメータ推定用の許容限界データが示す許容限界α_ｎ，ｍ（ｎ＝１，・・・，Ｎ， ｍ＝１，・・・，Ｍ）を用いて、上記の式（１）により忍耐力φ_ｎを計算すると共に、上記の式（２）によりアトラクション選好ベクトルψ_ｎを計算する。すなわち、パラメータ推定部１０１は、各来園者ｎに対して、各アトラクションｍの待ち時間の許容限界α_ｎ，ｍのうちの最大値を忍耐力φ_ｎとする。また、パラメータ推定部１０１は、各来園者ｎに対して、全アトラクションに対する各アトラクションｍの待ち時間の許容限界α_ｎ，ｍの割合をアトラクション選好ψ_ｎ，ｍとする。なお、アトラクション選好ψ_ｎ，ｍは０以上かつ１以下であり、その値が１に近いほど該当のアトラクションが選択されることを表す。

ここで、一例として、Ｎ＝５、Ｍ＝３である場合における忍耐力φ_ｎ及びアトラクション選好ベクトルψ_ｎの各要素であるアトラクション選好ψ_ｎ，ｍの具体例を図４に示す。図４では、アトラクション選好ψ_１，１＝０．１４、アトラクション選好ψ_１，２＝０．５７、アトラクション選好ψ_１，３＝０．２９である場合を示している。同様に、アトラクション選好ψ_２，１＝０．５１、アトラクション選好ψ_２，２＝０．２５、アトラクション選好ψ_２，３＝０．２４である場合を示している。ｎ＝３，４，５の場合も図４に示す通りである。

また、図４では、忍耐力φ_１＝９７．８、忍耐力φ_２＝５８．７、忍耐力φ_３＝１９．３、忍耐力φ_４＝７．５、忍耐力φ_５＝２５．９である場合を示している。

ステップＳ１０３：パラメータ推定部１０１は、上記のステップＳ１０２で計算された忍耐力φ_ｎ（ｎ＝１，・・・，Ｎ）が対数正規分布に従う（つまり、ｌｏｇ（φ_ｎ）～Ｎ（μ，σ^２））ものとして、パラメータμ及びσ^２を最尤法により推定する。この推定には、例えば、『C.M.ビショップ著「パターン認識と機械学習（上）ベイズ理論による統計的予測」p.24（1.2.4 ガウス分布）』に記載されている手法を用いればよい。

ステップＳ１０４：パラメータ推定部１０１は、上記のステップＳ１０２で計算されたアトラクション選好ベクトルψ_ｎがディリクレ分布Ｄｉｒ（β）に従う（つまり、ψ_ｎ～Ｄｉｒ（β））ものとして、パラメータβを最尤法により推定する。この推定には、例えば、『Thomas P. Minka, "Estimating a Dirichlet distribution", <URL: https://tminka.github.io/papers/dirichlet/minka-dirichlet.pdf>』に記載されている手法を用いればよい。なお、上述したように、パラメータβは、β＝（β_１，・・・，β_Ｍ）と表されるＭ次元ベクトルである。

以上により、許容限界モデルのパラメータμ、σ^２及びβが推定される。これらのパラメータμ、σ^２及びβは、例えば、パラメータ推定部１０１によって記憶部１０４に保存される。

＜シミュレーション用の来園者データ作成処理＞
以降では、シミュレーション用の来園者データを作成する処理について、図５を参照しながら説明する。図５は、シミュレーション用の来園者データを作成する処理の一例を示すフローチャートである。

ステップＳ２０１：来園者データ作成部１０２は、シミュレーションに用いる来園者数Ｉと、許容限界モデルのパラメータμ、σ^２及びβとを入力する。これらのパラメータμ、σ^２及びβは、パラメータ推定部１０１によって推定されたパラメータである。なお、来園者データ作成部１０２は、例えば、記憶部１０４に記憶されているパラメータμ、σ^２及びβを入力してもよいし、通信ネットワークを介して接続される他の装置から送信されたパラメータμ、σ^２及びβを入力してもよい。また、来園者データ作成部１０２は、例えば、記憶部１０４に記憶されている来園者数Ｉを入力してもよいし、通信ネットワークを介して接続される他の装置から送信された来園者数Ｉを入力してもよいし、キーボード等の入力装置で指定された来園者数Ｉを入力してもよい。

ステップＳ２０２：来園者データ作成部１０２は、正規分布Ｎ（μ，σ^２）に従う乱数ｒ_ｉ（ｉ＝１，・・・，Ｉ）を生成した上で、これらの乱数ｒ_ｉ（ｉ＝１，・・・，Ｉ）から忍耐力φ_ｉ（ｉ＝１，・・・，Ｉ）を生成する。すなわち、来園者データ作成部１０２は、各ｉ＝１，・・・，Ｉに対して、φ_ｉ＝ｅ^ｒｉとして、忍耐力φ_ｉを生成する。なお、ｅはネイピア数である。

ステップＳ２０３：来園者データ作成部１０２は、ディリクレ分布Ｄｉｒ（β）に従うＭ次元ベクトルをＩ個ランダムに生成し、これらのＭ次元ベクトルをアトラクション選好ベクトルψ_ｉ（ｉ＝１，・・・，Ｉ）とする。

ステップＳ２０４：来園者データ作成部１０２は、忍耐力φ_ｉとアトラクション選好ベクトルψ_ｉとを用いて、上記の式（３）により許容限界ベクトルα_ｉ（ｉ＝１，・・・，Ｉ）を生成する。上記の式（３）に示すように、ｉ＝１，・・・，Ｉに対して、来園者データ作成部１０２は、アトラクション選好ψ_ｉ，ｍの最大値が１となるように各アトラクション選好ψ_ｉ，ｍ（ｍ＝１，・・・，Ｍ）を正規化した上で、正規化後のψ_ｉ，ｍと忍耐力φ_ｉとの積が許容限界α_ｉ，ｍとなるように許容限界ベクトルα_ｉを生成する。

ここで、一例として、Ｍ＝３である場合における許容限界ベクトルα_ｉの生成の具体例を図６に示す。図６では、ｉ＝１の場合、ｍａｘψ_１＝０．５７であり、φ_１＝９７．８である。このため、許容限界ベクトルα_１の要素である許容限界α_１，１、α_１，２及びα_１，３は、α_１，１＝９７．８×（０．１４／０．５７）＝２４．０、α_１，２＝９７．８×（０．５７／０．５７）＝９７．８、α_１，３＝９７．８×（０．２９／０．５７）＝４９．７となる。

同様に、図６では、ｉ＝２の場合、ｍａｘψ_２＝０．５１であり、φ_２＝５８．７である。このため、許容限界ベクトルα_２の要素である許容限界α_２，１、α_２，２及びα_２，３は、α_２，１＝５８．７×（０．５１／０．５１）＝５８．７、α_２，２＝５８．７×（０．２５／０．５１）＝２８．７、α_２，３＝５８．７×（０．２４／０．５１）＝２７．６となる。ｉ＝３以降についても同様に、許容限界ベクトルα_ｉの各要素α_ｉ，１、α_ｉ，２及びα_ｉ，３が計算される。

ステップＳ２０５：来園者データ作成部１０２は、許容限界ベクトルα_ｉ（ｉ＝１，・・・，Ｉ）を用いて、シミュレーション用の来園者データを作成する。上述したように、シミュレーション用の来園者データには、アトラクションｍに対する来園者ｉの許容限界α_ｉ，ｍと、来園者ｉの入園時刻Ｉ_ｉ及び退園時刻Ｏ_ｉと、来園者ｉがアトラクションを体験する予定数ｋ_ｉとが含まれる。なお、入園時刻Ｉ_ｉ及び退園時刻Ｏ_ｉの組（ペア）は「滞在時間」と称されてもよい。

ここで、一例として、Ｍ＝３の場合におけるシミュレーション用の来園者データを図７に示す。図７では、ｉ＝１の場合、許容限界α_１，１＝２４．０、許容限界α_１，２＝９７．８、許容限界α_１，３＝４９．７、入園時刻Ｉ_１＝８：００、退園時刻Ｏ_１＝１４：００、予定数ｋ_１＝２である。同様に、ｉ＝２の場合、許容限界α_２，１＝５８．７、許容限界α_２，２＝２８．７、許容限界α_２，３＝２７．６、入園時刻Ｉ_２＝８：３０、退園時刻Ｏ_２＝１６：００、予定数ｋ_２＝２である。ｉ＝３以降についても図７に示す通りである。

ここで、来園者ｉの入園時刻Ｉ_ｉ及び退園時刻Ｏ_ｉには任意の時刻を設定することが可能であるが、一般に、実際のテーマパークは昼間に滞在人数がピークとなることが多いと考えられる。そこで、本発明の実施の形態では、テーマパーク内の来園者ｉの滞在人数が昼間にピークとなるように各来園者ｉの入園時刻Ｉ_ｉ及び退園時刻Ｏ_ｉを設定する。具体的には、テーマパークの開園時刻が「８：００」、閉園時刻が「２１：００」、１日あたり累計で３０００人が来園する場合を想定して、図８に示すように、シミュレーション時刻ｔをｔ＝０［分］～ｔ＝Ｔ＝７８０［分］まで、Ｉ＝３０００として、ｔ＝３００からｔ＝４００の間に滞在人数がピークとなるように各来園者ｉの入園時刻Ｉ_ｉ及び退園時刻Ｏ_ｉを設定した。なお、上述したように、Ｔはシミュレーション終了時刻である。

また、来園者ｉの予定数ｋ_ｉには０以上の任意の整数値を設定することが可能であるが、本発明の実施の形態では、平均３のポアソン分布（ただし、０は除く）に従うように各予定数ｋ_ｉを設定した。Ｉ＝３０００である場合の予定数ｋ_ｉ（ｉ＝１，・・・，３０００）のヒストグラムを図９に示す。

以上により、シミュレーション用の来園者データが作成される。この来園者データは、例えば、来園者データ作成部１０２によって記憶部１０４に保存される。

＜シミュレーションの条件＞
シミュレーションにより混雑度を推定する処理を説明する前に、その前提として、シミュレーション時における各種条件について説明する。

≪移動時間≫
シミュレーション時におけるアトラクション間及びテーマパークの出入口とアトラクションとの間を各来園者ｉが移動する際に要する時間は、シミュレーション用の移動時間データで与えられる。本発明の実施の形態では、Ｍ＝５であるものとして、各アトラクションｍの配置と、アトラクション間及び出入口－アトラクション間の移動経路とが図１０に示すものであるとする。なお、図１０中において、アトラクションを示す丸の中の数字はアトラクションのインデックスである。

このとき、図１１に示すシミュレーション用の移動時間データが与えられるものとする。図１１に示すシミュレーション用の移動時間データは、移動元から移動先への移動時間がマトリックス形式で表現されている。例えば、移動元が「出入口」、移動先が「ｍ＝１のアトラクション」である場合、その移動時間は３０［分］となる。同様に、移動元が「ｍ＝１のアトラクション」、移動先が「ｍ＝５のアトラクション」である場合、その移動時間は１０［分］となる。

≪アトラクションの処理能力≫
シミュレーション時におけるアトラクションの処理能力（つまり、単位時間あたりにアトラクションを体験可能な人数）は、シミュレーション用のアトラクションデータで与えられる。本発明の実施の形態では、Ｍ＝５であるものとして、図１２に示すシミュレーション用のアトラクションデータが与えられるものとする。図１２に示すシミュレーション用のアトラクションデータは、アトラクションの体験時間［分］と、アトラクションの定員［人］と、アトラクションの処理能力［人／分］と、アトラクションの稼働周期［分］とがマトリックス形式で表現されている。なお、処理能力は、「定員／体験時間」となる。

例えば、ｍ＝１のアトラクションは、体験時間が「５分］、定員が「１２人」、処理能力が「２．４」、稼働周期が「５」である。これは、ｍ＝１のアトラクションは、５分毎に稼働し、１回の稼働で１２人が同時に体験可能であり、１回の体験時間は５分間であることを表している。

なお、本発明の実施の形態では、シミュレーション用のアトラクションデータには、「処理能力」が含まれるものとしたが、これに限られず、少なくとも「処理能力を特定可能な情報」が含まれていればよい。処理能力を特定可能な情報とは、「定員」と「体験時間」との組であってもよいし、「処理能力」そのものであってもよい。

≪各来園者ｉの状態遷移等≫
シミュレーション時には、各来園者ｉには、図１３に示す各状態（「入園」、「退園」、「アトラクション選択」、「移動」、「待ち行列」、「体験」、「待機」）のうちのいずれか１つの状態やテーマパーク内の位置等が対応付けられている。つまり、例えば、来園者のインデックスｉと、この来園者の状態や位置等とが対応付けて記憶部１０４に記憶される。なお、初期状態においては、各来園者ｉには状態「入園」、位置「出入口」が対応付けられているものとする。

このとき、各シミュレーション時刻ｔにおいて、以下の（Ｃ１）～（Ｃ９）に示す条件に従って各来園者ｉの状態や位置等を更新するものとする。

（Ｃ１）状態「入園」の来園者ｉについて
入園時刻Ｉ_ｉがシミュレーション時刻ｔ以前である場合、来園者ｉの状態を「アトラクション選択」に更新する。このことは、シミュレーション時刻ｔが入園時刻Ｉ_ｉになった場合、来園者ｉがテーマパークに入園し、アトラクションの選択を開始することを意味する。

（Ｃ２）予定数ｋ_ｉ≠０かつ状態「アトラクション選択」の来園者ｉについて
上記の式（４）に示す多項線形モデルによりアトラクションｍの選択確率θ_{ｉ，ｍ，ｔ}を算出し、来園者ｉが次に体験するアトラクションｍを確率的に選択する。具体的には、待ち時間Ｗ_ｍ，ｔ＜α_ｉ，ｍを満たすアトラクションｍを候補として、選択確率θ_{ｉ，ｍ，ｔ}に基づいていずれか１つのアトラクションｍを選択する。

そして、いずれかのアトラクションｍが選択された場合は、来園者ｉの状態を「移動」に更新する。このとき、現在の位置から「選択されたアトラクションｍ」までの移動時間をシミュレーション時刻ｔに加算した「移動完了時刻」を当該来園者ｉに対応付けると共に、当該来園者ｉの位置を「選択されたアトラクションｍ」に更新する。移動時間は、上述したシミュレーション用の移動時間データから取得される。なお、各来園者ｉ間での移動時間の個人差を表現するために、シミュレーション用の移動時間データから取得された移動時間に対して乱数が加減乗除されてもよい。

一方で、いずれのアトラクションｍも選択されなかった場合（つまり、全てのアトラクションｍに対して、Ｗ_ｍ，ｔ≧α_ｉ，ｍであるため、θ_{ｉ，ｍ，ｔ}＝０となった場合）、来園者ｉの状態を「待機」に更新する。このとき、予め決められた待機時間（例えば、「３０分」等）をシミュレーション時刻ｔに加算した「待機終了時刻」を当該来園者ｉに対応付ける。なお、待機時間は予め決められていてもよいし、乱数を待機時間としてもよい。

（Ｃ３）予定数ｋ_ｉ＝０かつ状態「アトラクション選択」の来園者ｉについて
来園者ｉの状態を「退園」に更新する。また、当該来園者ｉの位置を「出入口」に更新する。このことは、当該来園者ｉが予定数のアトラクションを体験したため、テーマパークを退園することを意味する。

（Ｃ４）状態「移動」の来園者ｉについて
来園者ｉに対応付けられている移動完了時刻がシミュレーション時刻ｔ以前である場合、当該来園者ｉの状態を「待ち行列」に更新する。また、当該来園者ｉが体験するアトラクションｍの待ち時間Ｗ_ｍ，ｔをシミュレーション時刻ｔに加算した「体験開始時刻」と、当該アトラクションｍの体験時間を体験開始時刻に加算した「体験終了時刻」とを当該来園者ｉに対応付ける。なお、当該来園者ｉに対応付けられている移動完了時刻は削除される（又は、閉園時刻よりも後の時刻に更新されてもよい。）。

（Ｃ５）状態「待ち行列」の来園者ｉについて
来園者ｉに対応付けられている体験開始時刻がシミュレーション時刻ｔ以前である場合、当該来園者ｉの状態を「体験」に更新する。なお、当該来園者ｉに対応付けられている体験開始時刻は削除される（又は、閉園時刻よりも後の時刻に更新されてもよい。）。

（Ｃ６）状態「体験」の来園者ｉについて
来園者ｉに対応付けられている体験終了時刻がシミュレーション時刻ｔ以前である場合、当該来園者ｉの状態を「アトラクション選択」に更新する。また、当該来園者ｉの予定数ｋ_ｉを１減ずる。

（Ｃ８）状態「待機」の来園者ｉについて
来園者ｉに対応付けられている待機終了時刻がシミュレーション時刻ｔ以前である場合、当該来園者ｉの状態を「アトラクション選択」に更新する。

（Ｃ９）状態「入園」及び「退園」以外の全ての状態の来園者ｉについて
退園時刻Ｏ_ｉがシミュレーション時刻ｔ以前である場合、来園者ｉの状態を「退園」に更新する。また、当該来園者ｉの位置を「出入口」に更新する。このことは、シミュレーション時刻ｔが退園時刻Ｏ_ｉになった場合、来園者ｉがテーマパークを退園することを意味する。

＜シミュレーションによる混雑度推定処理＞
以降では、上記の「シミュレーションの条件」で述べた各種条件の下で、シミュレーションにより混雑度を推定する場合の処理について、図１４を参照しながら説明する。図１４は、シミュレーションにより混雑度を推定する処理の一例を示すフローチャートである。

ステップＳ３０１：シミュレーション部１０３は、シミュレーション用の来園者データと、シミュレーション用のアトラクションデータと、シミュレーション用の移動時間データとを入力する。なお、シミュレーション部１０３は、例えば、記憶部１０４に記憶されているこれらのデータを入力してもよいし、通信ネットワークを介して接続される他の装置から送信されたこれらのデータを入力してもよい。

ステップＳ３０２：シミュレーション部１０３は、シミュレーション時刻ｔをシミュレーション開始時刻に初期化すると共に、各アトラクションｍの待ち時間Ｗ_ｍ，ｔを０に初期化する。なお、シミュレーション開始時刻としてはｔ＝０とすればよいが、これに限られず、任意の時刻をシミュレーション開始時刻としてもよい。

ステップＳ３０３：シミュレーション部１０３は、上記の（Ｃ１）～（Ｃ９）に示す条件に従って各来園者ｉの状態や位置等を更新する。

ステップＳ３０４：次に、シミュレーション部１０３は、シミュレーション時刻ｔを更新する。すなわち、シミュレーション部１０３は、シミュレーション時刻ｔに１を加算する。

ステップＳ３０５：シミュレーション部１０３は、シミュレーション時刻ｔがシミュレーション終了時刻Ｔであるか否かを判定する。シミュレーション時刻ｔがシミュレーション終了時刻Ｔでないと判定された場合はステップＳ３０６に進む。一方で、シミュレーション時刻ｔがシミュレーション終了時刻Ｔであると判定された場合はステップＳ３０７に進む。

ステップＳ３０６：シミュレーション部１０３は、各アトラクションｍの待ち時間Ｗ_ｍ，ｔを算出し、記憶部１０４に保存する。ここで、待ち時間Ｗ_ｍ，ｔは、「シミュレーション時刻ｔにおいてアトラクションｍに並んでいる来園者ｉの人数／アトラクションｍの処理能力」で算出される。なお、シミュレーション時刻ｔにおいてアトラクションｍに並んでいる来園者ｉの人数とは、シミュレーション時刻ｔにおける状態が「待ち行列」で、かつ、位置が「アトラクションｍ」である来園者ｉの数のことである。

なお、シミュレーション部１０３は、待ち時間Ｗ_ｍ，ｔの算出及び保存を行った後、ステップＳ３０３に戻る。これにより、シミュレーション時刻ｔがシミュレーション終了時刻Ｔとなるまで、ステップＳ３０３～ステップＳ３０６が繰り返し実行される。

ステップＳ３０７：シミュレーション部１０３は、全ての来園者ｉの状態を「退園」、位置を「出入口」に更新する。この場合、シミュレーション時刻ｔ＝Ｔであり、テーマパークの閉園時刻に相当するためである。

以上により、シミュレーションの結果として、各シミュレーション時刻ｔにおける各アトラクションｍの待ち時間Ｗ_ｍ，ｔが得られる。これらの待ち時間Ｗ_ｍ，ｔが、各シミュレーション時刻ｔにおける各アトラクションｍの混雑度の推定結果である。

＜シミュレーション結果＞
次に、Ｍ＝５、Ｉ＝１０００，２０００，３０００，４０００，５０００，６０００，７０００，８０００，９０００とした場合におけるシミュレーション結果（つまり、各シミュレーション時刻ｔにおけるアトラクションｍの待ち時間Ｗ_ｍ，ｔの推移）を図１５に示す。なお、図１５に示すｍ＝１～ｍ＝５の待ち時間は、それぞれ、図５のステップＳ２０２及びステップＳ２０３で忍耐力φ_ｉ及びアトラクション選好ベクトルψ_ｉを生成する際の乱数のシードを変えて１０回ずつシミュレーションを行った結果得られた待ち時間Ｗ_ｍ，ｔの平均値である。

図１５に示すように、来園者数Ｉが増加するに従って全てのアトラクションｍで待ち時間が増加していることがわかる。このため、現実のテーマパークにおける各アトラクションの混雑度（待ち時間）を良く再現できており、本発明の実施の形態によって混雑度を高い精度で推定できていることがわかる。

したがって、本発明の実施の形態によれば、例えば、来園者の一部からアンケート等によって許容限界α_ｎ，ｍ（ｎ＝１，・・・，Ｎ， ｍ＝１，・・・，Ｍ）を事前に得ておくことで、各来園者のアトラクション選好を考慮した上で、Ｉ（＞Ｎ）人が来園した際の混雑度を容易に、かつ、高い精度で予測することが可能となる。

＜混雑度推定装置１０のハードウェア構成＞
最後に、本発明の実施の形態における混雑度推定装置１０のハードウェア構成について、図１６を参照しながら説明する。図１６は、本発明の実施の形態における混雑度推定装置１０のハードウェア構成の一例を示す図である。

図１６に示すように、本発明の実施の形態における混雑度推定装置１０は、入力装置２０１と、表示装置２０２と、外部Ｉ／Ｆ２０３と、ＲＡＭ（Random Access Memory）２０４と、ＲＯＭ（Read Only Memory）２０５と、プロセッサ２０６と、通信Ｉ／Ｆ２０７と、補助記憶装置２０８とを有する。これら各ハードウェアは、それぞれがバスＢを介して通信可能に接続されている。

入力装置２０１は、例えばキーボードやマウス、タッチパネル等であり、ユーザが各種操作を入力するのに用いられる。表示装置２０２は、例えばディスプレイ等であり、混雑度推定装置１０の処理結果等を表示する。なお、混雑度推定装置１０は、入力装置２０１及び表示装置２０２の少なくとも一方を有していなくてもよい。

外部Ｉ／Ｆ２０３は、外部装置とのインタフェースである。外部装置には、記録媒体２０３ａ等がある。混雑度推定装置１０は、外部Ｉ／Ｆ２０３を介して、記録媒体２０３ａ等の読み取りや書き込みを行うことができる。記録媒体２０３ａには、混雑度推定装置１０が有する各機能部（例えば、パラメータ推定部１０１、来園者データ作成部１０２及びシミュレーション部１０３）を実現する１以上のプログラムが記録されていてもよい。なお、記録媒体２０３ａには、例えば、ＣＤ（Compact Disc）やＤＶＤ（Digital Versatile Disk）、ＳＤメモリカード、ＵＳＢメモリカード等がある。

ＲＡＭ２０４は、プログラムやデータを一時保持する揮発性の半導体メモリである。ＲＯＭ２０５は、電源を切ってもプログラムやデータを保持することができる不揮発性の半導体メモリである。ＲＯＭ２０５には、例えば、ＯＳ（Operating System）に関する設定情報や通信ネットワークに関する設定情報等が格納されている。

プロセッサ２０６は、例えばＣＰＵ（Central Processing Unit）等であり、ＲＯＭ２０５や補助記憶装置２０８等からプログラムやデータをＲＡＭ２０４上に読み出して各種処理を実行する演算装置である。

通信Ｉ／Ｆ２０７は、混雑度推定装置１０を通信ネットワークに接続するためのインタフェースである。混雑度推定装置１０が有する各機能部を実現する１以上のプログラムは、通信Ｉ／Ｆ２０７を介して、所定のサーバ装置等から取得（ダウンロード）されてもよい。

補助記憶装置２０８は、例えばＨＤＤ（Hard Disk Drive）やＳＳＤ（Solid State Drive）等であり、プログラムやデータを格納している不揮発性の記憶装置である。補助記憶装置２０８に格納されているプログラムやデータには、例えば、ＯＳ、当該ＯＳ上で各種機能を実現するアプリケーションプログラム、混雑度推定装置１０が有する各機能部を実現する１以上のプログラム等がある。

本発明は、具体的に開示された上記の実施の形態に限定されるものではなく、請求の範囲の記載から逸脱することなく、種々の変形や変更が可能である。

１０ 混雑度推定装置
１０１ パラメータ推定部
１０２ 来園者データ作成部
１０３ シミュレーション部
１０４ 記憶部

Claims (4)

1. 総数Ｎの選択主体ｎ（ｎ＝１，・・・，Ｎ）の各々が総数Ｍの選択対象ｍ（ｍ＝１，・・・，Ｍ）のいずれかを選択するものとして、各選択主体ｎが許容する各選択対象ｍの混雑度の限界を表す許容限界α_ｎ，ｍを入力して、前記選択主体ｎの許容限界α_ｎ，ｍの前記選択対象ｍに関する最大値を表す忍耐力φ_ｎと、前記選択主体ｎが前記選択対象ｍを選択する際の選好を表す選好ベクトルψ_ｎとを算出し、算出した前記忍耐力φ_ｎと前記選好ベクトルψ_ｎとを用いて、総数Ｉ（＞Ｎ）の選択主体ｉ（ｉ＝１，・・・，Ｉ）の各々の許容限界α_ｉ，ｍを得るためのモデルのパラメータを推定するパラメータ推定手段と、
   推定した前記パラメータを用いた前記モデルにより、前記選択主体ｉの各々の許容限界α_ｉ，ｍを算出する算出手段と、
   前記許容限界α_ｉ，ｍと、前記選択対象ｍの各々の単位時間あたりの処理能力を特定可能な情報が含まれるシミュレーション条件とを入力して、各時刻ｔ（ｔ＝１，・・・・，Ｔ）において前記選択主体ｉが前記選択対象ｍのいずれかを選択することをシミュレーションすることで、前記各時刻ｔにおける前記選択対象ｍの各々の混雑度を推定する混雑度推定手段と、
   を有し、
   前記パラメータ推定手段は、
   前記忍耐力φ _ｎ が対数正規分布に従うものとして、前記対数正規分布に対応する正規分布Ｎ（μ，σ ^２ ）のパラメータμ及びσ ^２ を前記モデルのパラメータとして最尤法により推定し、
   前記選好ベクトルψ _ｎ がディリクレ分布Ｄｉｒ（β）に従うものとして、前記ディリクレ分布Ｄｉｒ（β）のパラメータβを前記モデルのパラメータとして最尤法により推定し、
   前記混雑度推定手段は、
   前記許容限界α _ｉ，ｍ から、各時刻ｔ（ｔ＝１，・・・・，Ｔ）において前記選択主体ｉが前記選択対象ｍのいずれかを選択する確率θ _{ｉ，ｍ，ｔ} を多項線形モデルにより算出し、
   算出した前記確率θ _{ｉ，ｍ，ｔ} と、前記選択対象ｍの各々の単位時間あたりの処理能力を特定可能な情報とに基づいて、前記各時刻ｔにおける前記選択対象ｍの各々の混雑度を推定する、ことを特徴とする混雑度推定装置。
2. 前記混雑度推定手段は、
   時刻ｔにおける前記選択対象ｍの各々の混雑度をＷ_ｍ，ｔとして、前記選択対象ｍのいずれかを選択する状態の各選択主体ｉについて、Ｗ_ｍ，ｔ＜α_ｉ，ｍを満たす選択対象ｍの中から前記確率θ_{ｉ，ｍ，ｔ}に従って１つの選択対象ｍを選択し、選択した選択対象ｍと、前記選択した選択対象ｍの単位時間当たりの処理能力を特定可能な情報とに基づいて、前記選択した選択対象ｍの混雑度を推定する、ことを特徴とする請求項１に記載の混雑度推定装置。
3. 総数Ｎの選択主体ｎ（ｎ＝１，・・・，Ｎ）の各々が総数Ｍの選択対象ｍ（ｍ＝１，・・・，Ｍ）のいずれかを選択するものとして、各選択主体ｎが許容する各選択対象ｍの混雑度の限界を表す許容限界α_ｎ，ｍを入力して、前記選択主体ｎの許容限界α_ｎ，ｍの前記選択対象ｍに関する最大値を表す忍耐力φ_ｎと、前記選択主体ｎが前記選択対象ｍを選択する際の選好を表す選好ベクトルψ_ｎとを算出し、算出した前記忍耐力φ_ｎと前記選好ベクトルψ_ｎとを用いて、総数Ｉ（＞Ｎ）の選択主体ｉ（ｉ＝１，・・・，Ｉ）の各々の許容限界α_ｉ，ｍを得るためのモデルのパラメータを推定するパラメータ推定手順と、
   推定した前記パラメータを用いた前記モデルにより、前記選択主体ｉの各々の許容限界α_ｉ，ｍを算出する算出手順と、
   前記許容限界α_ｉ，ｍと、前記選択対象ｍの各々の単位時間あたりの処理能力を特定可能な情報が含まれるシミュレーション条件とを入力して、各時刻ｔ（ｔ＝１，・・・・，Ｔ）において前記選択主体ｉが前記選択対象ｍのいずれかを選択することをシミュレーションすることで、前記各時刻ｔにおける前記選択対象ｍの各々の混雑度を推定する混雑度推定手順と、
   をコンピュータが実行し、
   前記パラメータ推定手順は、
   前記忍耐力φ _ｎ が対数正規分布に従うものとして、前記対数正規分布に対応する正規分布Ｎ（μ，σ ^２ ）のパラメータμ及びσ ^２ を前記モデルのパラメータとして最尤法により推定し、
   前記選好ベクトルψ _ｎ がディリクレ分布Ｄｉｒ（β）に従うものとして、前記ディリクレ分布Ｄｉｒ（β）のパラメータβを前記モデルのパラメータとして最尤法により推定し、
   前記混雑度推定手順は、
   前記許容限界α _ｉ，ｍ から、各時刻ｔ（ｔ＝１，・・・・，Ｔ）において前記選択主体ｉが前記選択対象ｍのいずれかを選択する確率θ _{ｉ，ｍ，ｔ} を多項線形モデルにより算出し、
   算出した前記確率θ _{ｉ，ｍ，ｔ} と、前記選択対象ｍの各々の単位時間あたりの処理能力を特定可能な情報とに基づいて、前記各時刻ｔにおける前記選択対象ｍの各々の混雑度を推定する、ことを特徴とする混雑度推定方法。
4. コンピュータを、請求項１又は２に記載の混雑度推定装置における各手段、として機能させるためのプログラム。

Images