WO2019138443A1

WO2019138443A1 - 移動経路推定装置、移動経路推定システムおよび移動経路推定方法

Info

Publication number: WO2019138443A1
Application number: PCT/JP2018/000183
Authority: WO
Inventors: 響介小西; 洋志亀田
Original assignee: 三菱電機株式会社
Priority date: 2018-01-09
Filing date: 2018-01-09
Publication date: 2019-07-18
Also published as: WO2019138638A1; JPWO2019138638A1; JP6632779B2

Abstract

移動経路推定装置（４）は、観測対象候補の消失および出現が反映された観測対象候補の存在確率の予測値および運動諸元確率分布の予測値に基づいて、観測時刻における観測対象候補の存在確率および運動諸元確率分布を更新する。移動経路推定装置（４）は、更新した観測対象候補の存在確率および運動諸元確率分布に基づいて、観測対象の移動経路を示す運動諸元推定値を算出する。

Description

移動経路推定装置、移動経路推定システムおよび移動経路推定方法

　この発明は、観測対象の移動経路を推定する移動経路推定装置、移動経路推定システムおよび移動経路推定方法に関する。

　従来から、観測対象の移動経路を推定する技術が知られている。例えば、非特許文献１には、観測領域に設置されたセンサによって得られた観測対象の検出データに基づいて、観測対象の移動経路を逐次推定する技術が記載されている。

Ｂ．Ｒｉｓｔｉｃ，　Ｂａ－Ｔｕｏｎｇ　Ｖｏ，　Ｂａ－Ｎｇｕ　Ｖｏ，　Ａ．Ｆａｒｉｎａ，　"Ａ　Ｔｕｔｏｒｉａｌ　ｏｎ　Ｂｅｒｎｏｕｌｌｉ　Ｆｉｌｔｅｒｓ：　Ｔｈｅｏｒｙ、　Ｉｍｐｌｅｍｅｎｔａｔｉｏｎ　ａｎｄ　Ａｐｐｌｉｃａｔｉｏｎｓ，　"　ＩＥＥＥ　Ｔｒａｎｓ．　ｏｎ　Ｓｉｇｎａｌ　Ｐｒｏｃ．　ｖｏ．６１，　ｎｏ．１３，　ｐｐ．　３４０６－３４３０，　Ｊｕｌ．　２０１３．

　非特許文献１に記載された従来の技術では、観測領域に同時に存在する観測対象の個数が既知であることを前提としているため、個数が未知の観測対象が観測領域に同時に存在する場合、観測対象のそれぞれの移動経路を推定できないという課題があった。

　この発明は上記課題を解決するものであり、個数が未知の観測対象の移動経路を推定することができる移動経路推定装置、移動経路推定システムおよび移動経路推定方法を得ることを目的とする。

　この発明に係る移動経路推定装置は、消失予測部、出現予測部、移動予測部、データ更新部および推定部を備える。消失予測部は、観測領域内の検出範囲を通過した観測対象候補が通過方向ごとに検出された入退検出データ、入退検出データが得られた入退検出センサ位置、観測時刻よりも一定時間だけ過去の時刻において観測対象候補が観測領域に存在する可能性を示す存在確率、および過去の時刻における観測対象候補の運動諸元確率分布に基づいて、過去の時刻から観測時刻までの間の観測対象候補の消失が反映された観測対象候補の存在確率の予測値および運動諸元確率分布の予測値を算出する。出現予測部は、入退検出データおよび入退検出センサ位置に基づいて、過去の時刻から観測時刻までの間に出現した観測対象候補の存在確率の予測値および運動諸元確率分布の予測値を算出する。移動予測部は、消失予測部と出現予測部によって算出された観測対象候補の存在確率の予測値および運動諸元確率分布の予測値に基づいて、過去の時刻から観測時刻までの間の観測対象候補の移動が反映された観測対象候補の存在確率の予測値および運動諸元確率分布の予測値を算出する。データ更新部は、検出範囲に存在する観測対象候補が検出された存在検出データ、存在検出データが得られた存在検出センサ位置および移動予測部によって算出された予測値に基づいて、観測時刻における観測対象候補の存在確率および運動諸元確率分布を算出する。推定部は、データ更新部によって算出された観測対象候補の存在確率および運動諸元確率分布に基づいて、観測時刻に観測対象が観測領域に存在するか否かを推定し、観測時刻に存在すると推定した観測対象の運動諸元推定値を算出する。

　この発明によれば、移動経路推定装置が、観測対象候補の消失および出現が反映された観測対象候補の存在確率の予測値および運動諸元確率分布の予測値に基づいて、観測時刻における観測対象候補の存在確率および運動諸元確率分布を更新する。移動経路推定装置は、更新した観測対象候補の存在確率および運動諸元確率分布に基づいて、観測対象の移動経路を示す運動諸元推定値を算出する。このように、移動経路推定装置は、観測対象候補の消失および出現を予測することにより、個数が未知の観測対象の移動経路を推定することができる。

この発明の実施の形態１に係る移動経路推定システムの構成を示すブロック図である。実施の形態１における観測領域の例を示す図である。観測領域におけるセンサ動作の例を示す図である。図４Ａは、実施の形態１に係る移動経路推定システムの機能を実現するハードウェア構成を示すブロック図である。図４Ｂは、実施の形態１に係る移動経路推定システムの機能を実現するソフトウェアを実行するハードウェア構成を示すブロック図である。実施の形態１に係る移動経路推定方法を示すフローチャートである。実施の形態１に係る移動経路推定装置の動作を示すフローチャートである。図６のステップＳＴ７の詳細な処理を示すフローチャートである。この発明の実施の形態２に係る移動経路推定システムの構成を示すブロック図である。実施の形態２に係る移動経路推定装置の動作を示すフローチャートである。この発明の実施の形態３に係る移動経路推定システムの構成を示すブロック図である。実施の形態３における観測領域の例を示す図である。実施の形態３に係る移動経路推定装置の動作を示すフローチャートである。この発明の実施の形態４に係る移動経路推定システムの構成を示すブロック図である。実施の形態４に係る移動経路推定装置の動作を示すフローチャートである。図１４のステップＳＴ７ｄの詳細な処理を示すフローチャートである。この発明の実施の形態５に係る移動経路推定システムの構成を示すブロック図である。実施の形態５における観測領域の例を示す図である。実施の形態５に係る移動経路推定方法を示すフローチャートである。実施の形態５に係る移動経路推定装置の動作を示すフローチャートである。図１９のステップＳＴ７ｆの詳細な処理を示すフローチャートである。

　以下、この発明をより詳細に説明するため、この発明を実施するための形態について、添付の図面に従って説明する。
実施の形態１．
　図１は、この発明の実施の形態１に係る移動経路推定システム１の構成を示すブロック図である。図１に示すように、移動経路推定システム１は、Ｎ（Ｎは２以上の整数）個の入退検出センサ２－１～２－Ｎ、Ｍ（Ｍは２以上の整数）個の存在検出センサ３－１～３－Ｍ、移動経路推定装置４、および移動経路保存部５を備えて構成される。移動経路推定装置４は、消失予測部４０、出現予測部４１、移動予測部４２、データ更新部４３および推定部４４を備える。

　以降の説明では、入退検出センサ２－１～２－Ｎおよび存在検出センサ３－１～３－Ｍによって検出された物体の消失または出現を考慮して、移動経路推定装置４によって移動経路が推定されるまでの物体を“観測対象候補”と呼ぶ。すなわち、センサによって検出された物体のうち、移動経路が推定されるまでの物体は、移動経路推定装置４において、仮の観測対象であるため、これを観測対象候補と定義する。

　入退検出センサ２－１～２－Ｎは、観測領域に設けられて、入退検出データを検出するセンサであり、光学カメラ、赤外線センサ、および測距センサなどのセンサによって実現される。また、入退検出データは、観測領域内の検出範囲を通過した観測対象候補が通過方向ごとに検出されたデータである。入退検出センサ２－１～２－Ｎは、入退検出データおよびセンサ位置を含む検出信号を観測時刻ごとに消失予測部４０および出現予測部４１に送信する。センサ位置は、入退検出センサが設置された位置であり、入退検出データが得られた入退検出センサ位置に相当する。

　存在検出センサ３－１～３－Ｍは、観測領域に設けられて、存在検出データを検出するセンサであって、光学カメラ、赤外線センサ、測距センサ、振動センサ、マイクロフォンおよび熱検知センサなどのセンサによって実現される。また、存在検出データは、センサの検出範囲に存在する観測対象候補が検出されたデータである。存在検出センサ３－１～３－Ｍは、存在検出データおよびセンサ位置を含む検出信号を、データ更新部４３に送信する。センサ位置は、存在検出センサが設置された位置であり、存在検出データが得られた存在検出センサ位置に相当する。
　入退検出センサ２－１～２－Ｎおよび存在検出センサ３－１～３－Ｍによって検出された複数の観測対象候補には、観測対象候補のそれぞれに固有な識別番号が付与される。

　図２は、実施の形態１における観測領域１０の例を示す図である。図２に示す観測領域１０は、出入り口を有した壁に囲まれた領域であり、出入り口を通ってのみ出入りが可能である。入退検出センサ２－ｎ（ｎは、１からＮのうちのいずれかの整数）は、例えば、観測領域１０の出入り口に設けられて、検出範囲１１－ｎが出入り口の領域に重なっている。出入り口を通って観測領域１０に出入りする観測対象候補は、検出範囲１１－ｎを通過したときに入退検出センサ２－ｎによって検出される。

　存在検出センサ３－ｍ（ｍは、１からＭのうちのいずれかの整数）は、図２に示すように、観測領域１０の任意の箇所に設けられる。存在検出センサ３－ｍの検出範囲１２－ｍは、例えば、破線で示すように、センサ位置を中心とした円形の領域を形成する。観測領域１０内の観測対象は、検出範囲１２－ｍに進入したときに存在検出センサ３－ｍに検出される。なお、検出範囲１２－ｍが円形である場合を示したが、楕円形、矩形などの円形以外の形状であってもよい。

　図３は、観測領域１０におけるセンサ動作の例を示す図であり、図２に示した観測領域１０におけるセンサの動作を示している。図３において、観測対象候補を検出したセンサは灰色で記載されている。人Ａ１が観測領域１０に進入する方向（矢印ａ１で示す方向）に移動する。このとき、入退検出センサ２－ｎは、出入り口における検出範囲１１－ｎを通過して観測領域１０内に進入した人Ａ１を検出する。人Ａ１が検出された入退検出データは、入退検出センサ２－ｎから移動経路推定装置４へ送信される。

　観測領域１０内を人Ａ１が移動して検出範囲１２－ｍを通過する度に、存在検出センサ３－ｍは、人Ａ１を検出する。人Ａ１が検出された存在検出データは、存在検出センサ３－ｍから移動経路推定装置４へ送信される。
　人Ａ１が、観測領域１０から退出する方向（矢印ａ２で示す方向）に移動すると、入退検出センサ２－ｎは、出入り口の検出範囲１１－ｎを通過して観測領域１０外に退出した人Ａ１を検出する。人Ａ１が検出された入退検出データは、入退検出センサ２－ｎから移動経路推定装置４へ送信される。

　同様に、人Ａ２が観測領域１０に進入する方向（矢印ｂ１で示す方向）に移動すると、入退検出センサ２－ｎは、出入り口における検出範囲１１－ｎを通過して観測領域１０内に進入した観測対象の人Ａ２を検出する。人Ａ２が検出された入退検出データは、入退検出センサ２－ｎから移動経路推定装置４へ送信される。

　観測領域１０内を人Ａ２が移動して検出範囲１２－ｍを通過する度に、存在検出センサ３－ｍは、人Ａ２を検出する。人Ａ２が検出された存在検出データは、存在検出センサ３－ｍから移動経路推定装置４へ送信される。
　人Ａ２が、観測領域１０から退出する方向（矢印ｂ２で示す方向）に移動すると、入退検出センサ２－ｎは、出入り口の検出範囲１１－ｎを通過して観測領域１０外に退出した人Ａ２を検出する。人Ａ２が検出された入退検出データは、入退検出センサ２－ｎから移動経路推定装置４へ送信される。

　人Ａ１および人Ａ２を検出した存在検出センサ３－ｍのそれぞれを辿っていくことで、観測領域１０内での人Ａ１の移動経路１３ａと人Ａ２の移動経路１３ｂとを推定することができる。なお、入退検出センサ２－ｎおよび存在検出センサ３－ｍは、図３に示すように、観測対象候補の誤検出または検出失敗が起こり得る。“観測対象候補の誤検出”とは、検出範囲に観測対象候補が存在しないにもかかわらず、センサによって観測対象候補が検出されることを意味する。また、“観測対象候補の検出失敗”とは、検出範囲に存在する観測対象候補が検出されないことを意味する。

　図１の説明に戻る。
　移動経路推定装置４は、入退検出センサ２－１～２－Ｎから入力した入退検出データ、および、存在検出センサ３－１～３－Ｍから入力した存在検出データを用いて、観測対象候補の運動諸元推定値を算出する。運動諸元推定値は、観測対象の移動経路を示すデータとして、時刻フレームごとに、かつ観測対象ごとに移動経路保存部５に保存される。移動経路保存部５に保存された移動経路を示すデータは、移動経路推定システム１の利用者からの要求に応じて出力される。

　入退検出センサ２－１～２－Ｎおよび存在検出センサ３－１～３－Ｍが、検出データを移動経路推定装置４に出力する時間の区切りを“時刻フレーム”と呼び、センサそれぞれの時刻フレームは同期しているものとする。センサによって観測対象候補が検出された観測時刻に相当する時刻フレームを“現時刻”と呼び、現時刻の１時刻フレーム前の過去の時刻に相当する時刻フレームを“前時刻”と呼ぶ。
　前時刻を現時刻よりも１時刻フレーム前の時刻と定義したが、これは例示であり、実施の形態１では、現時刻よりも一定時間だけ過去の時刻であれば、２以上の時刻フレームを遡った時刻であってもよい。

　消失予測部４０は、入退検出センサ２－１～２－Ｎから受信した入退検出データおよびセンサ位置を含む検出信号、前時刻における観測対象候補の存在確率および運動諸元確率分布に基づいて、前時刻から現時刻までの間の観測対象候補の消失が反映された観測対象候補の存在確率の予測値および運動諸元確率分布の予測値を算出する。なお、上記センサ位置は、入退検出データを検出した入退検出センサが設置された位置である。
　消失予測部４０によって算出された存在確率の予測値および運動諸元確率分布の予測値は、移動予測部４２に出力される。

　存在確率は、観測対象候補が観測領域に存在する可能性を示す確率であって、例えば、“存在確率が高い”とは、観測対象候補が観測領域に真に存在する観測対象である可能性が高いことを示している。存在確率は、０以上１以下のスカラー値で表される。

　運動諸元確率分布は、観測対象候補が観測領域に存在する場合において観測対象候補がどのような運動諸元であるかを示す確率分布である。運動諸元は、観測対象候補の位置、速度および加速度などである。例えば、２次元平面内の位置および２次元速度空間内の速度成分を運動諸元とした場合、運動諸元確率分布は、観測対象候補の２次元の位置および２次元の速度成分を表す４次元空間内における確率分布として定義される。

　出現予測部４１は、入退検出センサ２－１～２－Ｎから受信した入退検出データおよびセンサ位置を含む検出信号に基づいて、前時刻から現時刻までの間に出現した観測対象候補の存在確率の予測値および運動諸元確率分布の予測値を算出する。上記センサ位置は、入退検出データを検出した入退検出センサが設置された位置である。出現予測部４１によって算出された存在確率の予測値および運動諸元確率分布の予測値は、移動予測部４２に出力される。

　移動予測部４２は、消失予測部４０と出現予測部４１とから入力した観測対象候補の存在確率の予測値および運動諸元確率分布の予測値に基づいて、前時刻から現時刻までの間の観測対象候補の移動が反映された、観測対象候補の存在確率の予測値および運動諸元確率分布の予測値を算出する。移動予測部４２によって算出された存在確率の予測値および運動諸元確率分布の予測値は、データ更新部４３に出力される。

　データ更新部４３は、存在検出センサ３－１～３－Ｍから受信した存在検出データおよびセンサ位置を含む検出信号、移動予測部４２によって算出された観測対象候補の存在確率の予測値および運動諸元確率分布の予測値に基づいて、現時刻における観測対象候補の存在確率および運動諸元確率分布を算出する。上記センサ位置は、存在検出データを検出した存在検出センサが設置された位置である。データ更新部４３は、現時刻における１つまたは複数の観測対象候補の存在確率および運動諸元確率分布を算出し、算出した存在確率および運動諸元確率分布を、推定部４４および次の時刻フレームの予測を行う消失予測部４０に出力する。

　推定部４４は、データ更新部４３から入力した現時刻の観測対象候補の存在確率および運動諸元確率分布に基づいて、現時刻に観測対象が観測領域に存在するか否かを推定し、現時刻に存在すると推定された観測対象の運動諸元推定値を算出する。
　推定部４４によって算出された運動諸元推定値は、移動経路保存部５に出力されて保存される。移動経路保存部５に保存された観測対象の運動諸元推定値が時系列に並べられたデータが、観測対象の移動経路を示すデータとなる。
　なお、移動経路保存部５に保存された観測対象の移動経路を示すデータは、移動経路推定装置４の利用者の要求に応じて表示装置などに出力される。

　図４Ａは、移動経路推定システム１の機能を実現するハードウェア構成を示すブロック図である。図４Ｂは、移動経路推定システム１の機能を実現するソフトウェアを実行するハードウェア構成を示すブロック図である。図４Ａおよび図４Ｂにおいて、センサ群１００は、入退検出センサ２－１～２－Ｎに相当するセンサ群である。センサ群１０１は、存在検出センサ３－１～３－Ｍに相当するセンサ群である。記憶装置１０２には、移動経路推定装置４による推定処理の実行中に得られる中間データが記憶され、さらに、移動経路保存部５が構築されている。

　移動経路推定装置４における、消失予測部４０、出現予測部４１、移動予測部４２、データ更新部４３および推定部４４のそれぞれの機能は、処理回路によって実現される。
　すなわち、移動経路推定装置４は、図５を用いて後述するステップＳＴ１からステップＳＴ５までの処理を実行するための処理回路を備える。
　処理回路は、専用のハードウェアであってもよいが、メモリに記憶されたプログラムを実行するＣＰＵ（Ｃｅｎｔｒａｌ　Ｐｒｏｃｅｓｓｉｎｇ　Ｕｎｉｔ）であってもよい。

　処理回路が図４Ａに示す専用のハードウェアの処理回路１０３である場合、処理回路１０３は、例えば、単一回路、複合回路、プログラム化したプロセッサ、並列プログラム化したプロセッサ、ＡＳＩＣ（Ａｐｐｌｉｃａｔｉｏｎ　Ｓｐｅｃｉｆｉｃ　Ｉｎｔｅｇｒａｔｅｄ　Ｃｉｒｃｕｉｔ）、ＦＰＧＡ（Ｆｉｅｌｄ－Ｐｒｏｇｒａｍｍａｂｌｅ　Ｇａｔｅ　Ａｒｒａｙ）またはこれらを組み合わせたものが該当する。
　消失予測部４０、出現予測部４１、移動予測部４２、データ更新部４３および推定部４４のそれぞれの機能を別々の処理回路で実現してもよいし、これらの機能をまとめて１つの処理回路で実現してもよい。

　処理回路が図４Ｂに示すプロセッサ１０４である場合、消失予測部４０、出現予測部４１、移動予測部４２、データ更新部４３および推定部４４のそれぞれの機能は、ソフトウェア、ファームウェアまたはソフトウェアとファームウェアとの組み合わせによって実現される。ソフトウェアまたはファームウェアは、プログラムとして記述されて、メモリ１０５に記憶される。

　プロセッサ１０４は、メモリ１０５に記憶されたプログラムを読み出して実行することによって、消失予測部４０、出現予測部４１、移動予測部４２、データ更新部４３および推定部４４のそれぞれの機能を実現する。
　すなわち、移動経路推定装置４は、プロセッサ１０４によって実行されるときに、図５に示すステップＳＴ１からステップＳＴ５までの処理が結果的に実行されるプログラムを記憶するためのメモリ１０５を備える。これらのプログラムは、消失予測部４０、出現予測部４１、移動予測部４２、データ更新部４３および推定部４４の手順または方法をコンピュータに実行させるものである。メモリ１０５は、コンピュータを、消失予測部４０、出現予測部４１、移動予測部４２、データ更新部４３および推定部４４として機能させるためのプログラムが記憶されたコンピュータ可読記憶媒体であってもよい。

　メモリ１０５には、例えば、ＲＡＭ（Ｒａｎｄｏｍ　Ａｃｃｅｓｓ　Ｍｅｍｏｒｙ）、ＲＯＭ（Ｒｅａｄ　Ｏｎｌｙ　Ｍｅｍｏｒｙ）、フラッシュメモリ、ＥＰＲＯＭ（Ｅｒａｓａｂｌｅ　Ｐｒｏｇｒａｍｍａｂｌｅ　Ｒｅａｄ　Ｏｎｌｙ　Ｍｅｍｏｒｙ）、ＥＥＰＲＯＭ（Ｅｌｅｃｔｒｉｃａｌｌｙ－ＥＰＲＯＭ）などの不揮発性または揮発性の半導体メモリ、磁気ディスク、フレキシブルディスク、光ディスク、コンパクトディスク、ミニディスク、ＤＶＤなどが該当する。

　消失予測部４０、出現予測部４１、移動予測部４２、データ更新部４３および推定部４４のそれぞれの機能について一部を専用のハードウェアで実現し、一部をソフトウェアまたはファームウェアで実現してもよい。
　例えば、消失予測部４０、出現予測部４１および移動予測部４２のそれぞれは、専用のハードウェアとしての処理回路で機能を実現する。データ更新部４３および推定部４４については、プロセッサ１０４が、メモリ１０５に記憶されたプログラムを読み出して実行することにより機能を実現してもよい。
　このように、処理回路は、ハードウェア、ソフトウェア、ファームウェアまたはこれらの組み合わせにより上記機能のそれぞれを実現することができる。

　次に動作について説明する。
　図５は、実施の形態１に係る移動経路推定方法を示すフローチャートである。
　消失予測部４０が、入退検出センサ２－１～２－Ｎから受信した入退検出データおよびセンサ位置、および、データ更新部４３から入力した前時刻における観測対象候補の存在確率および運動諸元確率分布に基づいて、前時刻から現時刻までの間の観測対象候補の消失が反映された観測対象候補の存在確率の予測値および運動諸元確率分布の予測値を算出する（ステップＳＴ１）。

　次に、出現予測部４１が、入退検出センサ２－１～２－Ｎから受信した入退検出データおよびセンサ位置に基づいて、前時刻から現時刻までの間に出現した観測対象候補の存在確率の予測値および運動諸元確率分布の予測値を算出する（ステップＳＴ２）。

　移動予測部４２が、消失予測部４０および出現予測部４１から入力した観測対象候補の存在確率の予測値および運動諸元確率分布の予測値に基づいて、前時刻から現時刻までの間の観測対象候補の移動が反映された観測対象候補の存在確率の予測値および運動諸元確率分布の予測値を算出する（ステップＳＴ３）。

　データ更新部４３が、存在検出センサ３－１～３－Ｍから受信した存在検出データおよびセンサ位置、および、移動予測部４２から入力した観測対象候補の存在確率の予測値および運動諸元確率分布の予測値に基づいて、現時刻における観測対象候補の存在確率および運動諸元確率分布を算出する（ステップＳＴ４）。

　推定部４４が、データ更新部４３から入力した観測対象候補の存在確率および運動諸元確率分布に基づいて、現時刻に観測対象が観測領域に存在するか否かを推定し、現時刻に存在すると推定した観測対象の移動経路を推定する（ステップＳＴ５）。例えば、推定部４４は、観測対象の移動経路を示す運動諸元推定値を算出し、算出した運動諸元推定値を移動経路保存部５に保存する。

　次に、観測対象の移動経路の推定処理について詳細に説明する。
　以降では、観測対象の運動諸元として２次元位置空間内の位置および２次元速度空間内の速度成分を推定する場合について説明する。
　２次元位置空間の座標系は、直交するＸＹ座標軸で表され、２次元速度成分もＸ軸方向の値とＹ軸方向の値があるものとする。ただし、これらの前提は例示であって、これらの運動諸元に限定されるものではない。
　例えば、３次元空間を移動する観測対象について３次元位置空間内の位置および３次元速度空間内の速度を運動諸元としてもよい。観測対象が加速する可能性を考慮して、加速度成分、加速度の微分成分を運動諸元に含めてもよい。さらに、直交座標系の代わりに、極座標形を利用してもよい。

　また、移動経路推定装置４では、観測対象が等速直進運動を行うことを前提として運動諸元が推定される。なお、移動経路推定装置４では、等速直進運動の代わりに、例えば、等加速度直進運動、等角速度旋回運動などを前提として運動諸元を推定してもよい。
　入退検出センサ２－１～２－Ｎは、入退検出データの誤検出があるものと仮定し、存在検出センサ３－１～３－Ｍについても、存在検出データの誤検出または存在検出データの検出失敗があるものと仮定する。

　さらに、運動諸元確率分布を、重み付けされた有限個数の“粒子”の集合で近似する。運動諸元確率分布の値が大きい領域には、重み係数が大きい粒子が密に分布し、運動諸元確率分布の値が小さい領域には、重み係数が小さな粒子が疎に分布する。なお、近似する方法は、前述した方法に限定されるものではなく、例えば、ガウス分布の確率密度関数の線形結合で運動諸元確率分布を近似してよい。

　現時刻の時刻フレーム番号を“ｋ”として、現時刻を“時刻ｋ”と記載する。前時刻の時刻フレーム番号を“ｋ－１”として前時刻を“時刻ｋ－１”と記載する。ｋの添え字が付く記号は、その記号が示すパラメータ値が時刻ｋにおける値であることを示している。ｋ－１の添え字についても同様である。

　観測対象の運動諸元を表す４行１列のベクトルｘは、下記式（１）で定義される。
　ただし、下記式（１）において、ｘ^（１）は位置のＸ軸成分を表し、ｘ^（２）は位置のＹ軸成分を表し、ｘ^（３）は速度のＸ軸成分を表し、ｘ^（４）は速度のＹ軸成分を表している。
　行列右上の添え字Ｔは行列の転置を表している。このような運動諸元を表すベクトルを“状態ベクトル”と呼ぶ。
　ｘ＝［ｘ^（１）　ｘ^（２）　ｘ^（３）　ｘ^（４）］^Ｔ　　　　　　　・・・（１）

　また、ｎ番目（ｎ＝１，２，・・・，Ｎ）の入退検出センサ２－ｎのセンサ位置を示す２行１列のベクトルｚ_ｎ ^{（ＩｎＯｕｔ）}を、下記式（２）で定義する。
　下記式（２）において、ｚ_ｎ ^{（ＩｎＯｕｔ、１）}は、入退検出センサ２－ｎのセンサ位置のＸ軸成分であり、ｚ_ｎ ^{（ＩｎＯｕｔ、２）}は、入退検出センサ２－ｎのセンサ位置のＹ軸成分である。
　ｚ_ｋ，ｎ ^{（ＩｎＯｕｔ）}＝［ｚ_ｋ，ｎ ^{（ＩｎＯｕｔ，１）}　ｚ_ｋ，ｎ ^{（ＩｎＯｕｔ，２）}］^Ｔ　　　・・・（２）

　さらに、ｍ番目（ｍ＝１，２，・・・，Ｍ）の存在検出センサ３－ｍのセンサ位置を示す２行１列のベクトルｚ_ｋ，ｍ ^{（Ｅｘｉｓｔ）}を、下記式（３）で定義する。
　下記式（３）において、ｚ_ｋ，ｍ ^{（Ｅｘｉｓｔ，１）}は、存在検出センサ３－ｍのセンサ位置のＸ軸成分であり、ｚ_ｋ，ｍ ^{（Ｅｘｉｓｔ，２）}は、存在検出センサ３－ｍのセンサ位置のＹ軸成分である。
　ｚ_ｋ，ｍ ^{（Ｅｘｉｓｔ）}＝［ｚ_ｋ，ｍ ^{（Ｅｘｉｓｔ，１）}　ｚ_ｋ，ｍ ^{（Ｅｘｉｓｔ，２）}］^Ｔ　　　・・・（３）

　時刻ｋにおける観測対象候補の識別番号の集合をＬ_ｋとし、時刻ｋにおける識別番号ｉの観測対象候補の存在確率をｒ_ｋ｜ｋ ^（ｉ）とする。時刻ｋにおける識別番号ｉの観測対象候補の存在確率の予測値ｒ_{ｋ｜ｋ－１} ^（ｉ）は、時刻ｋ－１における識別番号ｉの観測対象候補の存在確率をｒ_{ｋ－１｜ｋ－１} ^（ｉ）を用いて算出される。

　時刻ｋにおける、識別番号ｉの観測対象候補の運動諸元確率分布ｆ_ｋ｜ｋ ^（ｉ）（ｘ）は、時刻ｋにおいて識別番号ｉの観測対象候補の運動諸元が状態ベクトルｘである確率を示している。重み付けされた粒子の集合で運動諸元確率分布ｆ_ｋ｜ｋ ^（ｉ）（ｘ）を近似することで、識別番号ｉの観測対象候補の運動諸元確率分布ｆ_ｋ｜ｋ ^（ｉ）（ｘ）は、下記式（４）で表すことができる。下記式（４）において、δはディラックのデルタ関数であり、粒子の個数Ｊは事前設定パラメータとする。ｗ_ｋ｜ｋ ^{（ｉ，ｊ）}はｊ番目（ｊ＝１，２，・・・，Ｊ）の粒子の重み係数を示す０以上１以下のスカラー値であり、ｘ_ｋ｜ｋ ^{（ｉ，ｊ）}はｊ番目の粒子の状態を示す状態ベクトルである。

　時刻ｋおける、識別番号ｉの観測対象候補の運動諸元確率分布の予測値ｆ_{ｋ｜ｋ－１} ^（ｉ）（ｘ）は、時刻ｋ－１における識別番号ｉの観測対象候補の運動諸元確率分布を用いて算出される。重み付けされた粒子の集合でｆ_{ｋ｜ｋ－１} ^（ｉ）（ｘ）を近似することで、運動諸元確率分布の予測値ｆ_{ｋ｜ｋ－１} ^（ｉ）（ｘ）は、下記式（５）で表すことができる。
　ただし、ｗ_{ｋ｜ｋ－１} ^{（ｉ，ｊ）}は、ｊ番目（ｊ＝１，２，・・・，Ｊ）の粒子の重み係数を示す０以上１以下のスカラー値であり、ｘ_{ｋ｜ｋ－１} ^{（ｉ，ｊ）}は、ｊ番目の粒子の状態を示す状態ベクトルである。

　図６は、移動経路推定装置４の動作を示すフローチャートである。図６では、運動諸元確率分布ｆ_ｋ｜ｋ ^（ｉ）（ｘ）を、重み付けされた粒子の集合で近似することを前提として、時刻ｋにおける識別番号ｉの観測対象候補の運動諸元確率分布ｆ_ｋ｜ｋ ^（ｉ）（ｘ）は、Ｊ個の重み係数ｗ_ｋ｜ｋ ^{（ｉ，ｊ）}および粒子状態ｘ_ｋ｜ｋ ^{（ｉ，ｊ）}を用いて算出される。なお、運動諸元確率分布の予測値ｆ_{ｋ｜ｋ－１} ^（ｉ）（ｘ）についても同様に算出される。
　以降では、ある確率密度関数ｆ（ｘ）に従って生成された乱数がｘであることを示す式を、ｘ～ｆ（ｘ）と記載する。

　さらに、確率密度関数として、Ｎ（ｘ；μ，Ｐ）およびＵ（ｘ，ａ，ｂ）を使用する。
　Ｎ（ｘ；μ，Ｐ）は、確率変数ベクトルｘに関する平均値ベクトルμ、誤差共分散行列Ｐの多変量ガウス分布の確率密度関数である。
　Ｕ（ｘ，ａ，ｂ）は、確率変数スカラー値ｘに関する、最小値ａおよび最大値ｂの一様分布の確率密度関数である。

　消失予測部４０および出現予測部４１は、入退検出センサ２－１～２－Ｎから受信した時刻ｋにおけるＮ個の入退検出データのうち、観測領域に対する観測対象候補の進入または退出が検出されたデータを、任意の順番で選択して入力する（ステップＳＴ１ａ）。これにより入退検出データを処理するループが開始される。
　以降では、ステップＳＴ２ａからステップＳＴ４ａまででは、ｎ番目（ｎ＝１，２，・・・，Ｎのいずれか）の入退検出センサ２－ｎによって検出された入退検出データが選択されたものと仮定する。

　消失予測部４０および出現予測部４１は、ｎ番目の入退検出センサ２－ｎによって検出された入退検出データに基づいて、観測対象候補の移動方向が観測領域への進入方向であるか否かを判定する（ステップＳＴ２ａ）。
　消失予測部４０は、観測対象候補の移動方向が観測領域からの退出方向であると判定した場合（ステップＳＴ２ａ；ＮＯ）、ステップＳＴ３ａへ移行する。
　出現予測部４１は、観測対象候補の移動方向が観測領域への進入方向であると判定した場合（ステップＳＴ２ａ；ＹＥＳ）、ステップＳＴ４ａへ移行する。

　なお、入退検出センサが、観測領域への進入方向と観測領域からの退出方向とを同一の時刻フレームで検出可能である場合、入退検出センサによって検出された観測対象候補の移動方向に基づいて、消失予測部４０がステップＳＴ３ａの処理を実行し、出現予測部４１が、ステップＳＴ４ａの処理を実行してもよい。

　ステップＳＴ３ａにおいて、消失予測部４０は、入退検出センサ２－ｎによって検出された入退検出データ、上記式（２）で示した入退検出センサ２－ｎのセンサ位置、時刻ｋ－１の観測対象候補の存在確率および運動諸元確率分布に基づいて、時刻ｋ－１から時刻ｋまでの間の観測対象候補の消失が反映された、観測対象候補の存在確率の予測値および運動諸元確率分布の予測値を算出する。

　なお、消失予測部４０が、ステップＳＴ３ａの処理を１回目に行う場合、時刻ｋ－１における観測対象候補の存在確率および運動諸元確率分布を用いる。２回目以降にステップＳＴ３ａの処理を行う場合、消失予測部４０は、直前のステップＳＴ３ａで算出した存在確率の予測値および運動諸元確率分布の予測値を用いる。

　消失予測部４０は、下記式（６）から下記式（１１）までを用いて時刻ｋ－１から時刻ｋまでの間の観測対象候補の消失が反映された、識別番号ｉの観測対象候補の存在確率の予測値ｒ_{ｋ｜ｋ－１，Ｓ} ^（ｉ）を算出し、時刻ｋ－１から時刻ｋまでの間の観測対象候補の消失が反映された、識別番号ｉの観測対象候補の運動諸元確率分布の予測値を表す粒子である｛ｗ_{ｋ｜ｋ－１，Ｓ} ^{（ｉ，ｊ）}，ｘ_{ｋ｜ｋ－１，Ｓ} ^{（ｉ，ｊ）}｝_{ｊ＝１，・・・，Ｊ}を算出する。
　識別番号ｉは、時刻ｋ－１（前時刻）の観測対象候補と同一であるものとする。
　なお、｜｜ｘ｜｜_２は、ベクトルｘのユークリッドノルムである。β_{ｋ，ＦＬＳ} ^{（ＩｎＯｕｔ）}およびｄ_Ｔｈは事前に設定されたパラメータであり、特に、β_{ｋ，ＦＬＳ} ^{（ＩｎＯｕｔ）}は、入退検出センサが誤検出する確率である。

　上記式（６）から上記式（１１）までを用いて観測対象候補の存在確率の予測値および運動諸元確率分布の予測値を算出する処理は、識別番号ｉの観測対象候補の運動諸元確率分布の平均位置Ｈ_Ｐｏｓμ_{ｋ－１｜ｋ－１} ^（ｉ）が、観測対象候補の退出を検出したｎ番目の入退検出センサ２－ｎが設置されたセンサ位置ｚ_ｎ ^{（ＩｎＯｕｔ）}に近くなるにつれて、識別番号ｉの観測対象候補の存在確率をβ_{ｋ，ＦＬＳ} ^{（ＩｎＯｕｔ）}倍に減少させる処理に相当する。この処理により、入退検出センサ２－ｎによって検出された観測対象候補の退出を、観測対象候補の予測値に反映させることができる。

　識別番号ｉの観測対象候補の存在確率をβ_{ｋ，ＦＬＳ} ^{（ＩｎＯｕｔ）}倍することで、識別番号ｉの観測対象候補の存在確率に対して観測対象候補の退出（消失）を反映させたが、これに限定されるものではない。例えば、β_{ｋ，ＦＬＳ} ^{（ＩｎＯｕｔ）}＝０として、上記式（９）に示すｄ^{（ｉ，ｎ）}が閾値よりも小さい場合、消失予測部４０が、識別番号ｉの観測値候補の存在確率および運動諸元確率分布（粒子）を削除してもよい。

　なお、存在確率を減らす観測対象候補を、運動諸元確率分布の平均位置のみに基づいて選択する場合を示したが、これに限定されるものではない。
　例えば、消失予測部４０は、運動諸元確率分布の平均速度ベクトルの方位角度と、入退検出センサが、観測領域から退出する観測対象候補を検出する方向の方位角度と、の差分を用いて、存在確率を減らす観測対象候補を選択してもよい。

　ステップＳＴ４ａにおいて、出現予測部４１は、入退検出センサ２－ｎによって検出された入退検出データ、および上記式（２）で示した入退検出センサ２－ｎのセンサ位置に基づいて、時刻ｋ－１から時刻ｋまでの間に出現した観測対象候補の存在確率の予測値および運動諸元確率分布の予測値を算出する。

　出現予測部４１は、下記式（１２）から下記式（１６）までを用いて、時刻ｋ－１から時刻ｋまでの間に出現した観測対象候補の存在確率の予測値ｒ_{ｋ｜ｋ－１，Ｂ} ^（ｉｎ）を算出し、時刻ｋ－１から時刻ｋまでの間に出現した観測対象候補の運動諸元確率分布の予測値を表す粒子である｛ｗ_{ｋ｜ｋ－１，Ｂ} ^{（ｉｎ，ｊ）}，ｘ_{ｋ｜ｋ－１，Ｂ} ^{（ｉｎ，ｊ）}｝_{ｊ＝１，・・・，Ｊ}を算出する。
　識別番号ｉ_ｎは、ｎ番目の入退検出センサ２－ｎの周辺に出現したと予測された観測対象候補の識別番号であり、既存の観測対象候補の識別番号と異なる値である。
　β_{ｋ，ＦＬＳ} ^{（ＩｎＯｕｔ）}、Ｖ_Ｂ ^（３）およびＶ_Ｂ ^（４）は、事前設定パラメータであって、特に、β_{ｋ，ＦＬＳ} ^{（ＩｎＯｕｔ）}は、入退検出センサが観測対象候補を誤検出する確率を示している。

　上記式（１２）から上記式（１６）までを用いて存在確率の予測値および運動諸元確率分布の予測値を算出する処理は、ｎ番目の入退検出センサ２－ｎの設置位置を初期位置として新規の観測対象候補の粒子を生成し、生成した粒子のそれぞれの速度を一様分布の乱数として、存在確率と粒子の重み係数とを固有値とする識別番号ｉ_ｎの観測対象候補を追加する処理に相当する。この処理により、時刻ｋに入退検出センサ２－ｎによって観測領域への進入が検出された観測対象候補の予測値が得られる。

　なお、新たに出現した観測対象候補に関する粒子の初期速度を、上記式（１５）および上記式（１６）に示した０を中心とする一様分布の乱数としたが、ｎ番目の入退検出センサ２－ｎが観測領域に進入する観測対象候補を検出する方向に偏った分布の乱数を用いてもよい。

　ステップＳＴ５ａにおいて、消失予測部４０および出現予測部４１は、入退検出センサ２－１～２－Ｎから受信した時刻ｋにおけるＮ個の入退検出データのうち、観測領域に対する観測対象候補の進入または退出が検出された全てのデータを選択したか否かを判定する。未選択のデータがあると判定した場合（ステップＳＴ５ａ；ＮＯ）、消失予測部４０および出現予測部４１は、未選択の入退検出データを選択して、ステップＳＴ２ａの処理に戻る。全てのデータを選択していれば（ステップＳＴ５ａ；ＹＥＳ）、ステップＳＴ６ａの処理に移行する。

　ステップＳＴ６ａにおいて、移動予測部４２は、時刻ｋ－１から時刻ｋまでの間の観測対象候補の消失が反映された存在確率の予測値および運動諸元確率分布の予測値と、時刻ｋ－１から時刻ｋまでの間に出現した観測対象候補の存在確率の予測値および運動諸元確率分布の予測値に基づいて、時刻ｋ－１から時刻ｋまでの間の観測対象候補の移動が反映された観測対象候補の存在確率の予測値および運動諸元確率分布の予測値を算出する。

　時刻ｋ－１（前時刻）から時刻ｋ（現時刻）までの間の観測対象候補の消失が反映された、観測対象候補の存在確率の予測値をｒ_{ｋ｜ｋ－１，Ｓ} ^（ｉｓ）で表す。時刻ｋ－１から時刻ｋまでの間の観測対象候補の消失が反映された観測対象候補の運動諸元確率分布の予測値を｛ｗ_{ｋ｜ｋ－１，Ｓ} ^{（ｉｓ，ｊ）}，ｘ_{ｋ｜ｋ－１，Ｓ} ^{（ｉｓ，ｊ）}｝_{ｊ＝１，・・・，Ｊ}で表す。ｉ_Ｓは、消失予測部４０によって算出された観測対象候補の識別番号とする。

　さらに、時刻ｋ－１から時刻ｋまでの間に出現した観測対象候補の存在確率の予測値をｒ_{ｋ｜ｋ－１，Ｂ} ^（ｉＢ）で表す。時刻ｋ－１から時刻ｋまでの間に出現した観測対象候補の運動諸元確率分布の予測値を｛ｗ_{ｋ｜ｋ－１，Ｂ} ^{（ｉＢ，ｊ）}，ｘ_{ｋ｜ｋ－１，Ｂ} ^{（ｉＢ，ｊ）}｝_{ｊ＝１，・・・，Ｊ}で表す。ｉ_Ｂは、出現予測部４１によって算出された観測対象候補の識別番号とする。

　さらに、時刻ｋ－１から時刻ｋまでの間の観測対象候補の移動が反映された存在確率の予測値をｒ_{ｋ｜ｋ－１} ^（ｉ）で表す。時刻ｋ－１から時刻ｋまでの間の観測対象候補の移動が反映された運動諸元確率分布の予測値を｛ｗ_{ｋ｜ｋ－１} ^{（ｉ，ｊ）}，ｘ_{ｋ｜ｋ－１} ^{（ｉ，ｊ）}｝_{ｊ＝１，・・・，Ｊ}で表す。識別番号ｉは、識別番号ｉ_Ｓと識別番号ｉ_Ｂとの両方を含んだ識別番号とする。すなわち、識別番号ｉ_Ｓの集合をＬ_{ｋ｜ｋ－１，Ｓ}で表し、識別番号ｉ_Ｂの集合をＬ_{ｋ｜ｋ－１、Ｂ}で表した場合、識別番号ｉは、下記式（１７）で定義される。
　ｉ∈Ｌ_{ｋ｜ｋ－１}＝｛Ｌ_{ｋ｜ｋ－１，Ｓ}　∪　Ｌ_{ｋ｜ｋ－１，Ｂ}｝　　　・・・（１７）

　移動予測部４２は、下記式（１８）から下記式（２１）までを用いて、観測対象候補の移動が反映された観測対象候補の存在確率の予測値ｒ_{ｋ｜ｋ－１} ^（ｉ）および観測対象候補の運動諸元確率分布の予測値｛ｗ_{ｋ｜ｋ－１} ^{（ｉ，ｊ）}，ｘ_{ｋ｜ｋ－１} ^{（ｉ，ｊ）}｝_{ｊ＝１，・・・，Ｊ}を算出する。ただし、Φ_ｋは、等速直進運動の予測モデルに基づいて状態ベクトルを遷移させる行列である。τ_ｋは、時刻ｋ－１から時刻ｋまでの経過時間である。Ｑ_ｋは、等速直進運動の予測モデルの誤差を表すシステム雑音を示す行列であって、事前設定パラメータである。

　これまで、ガウス分布に基づいた乱数が生成可能であることを前提として、存在確率の予測値と運動諸元確率分布の予測値とを算出する場合を示したが、これに限定されるものではない。例えば、既存の粒子フィルタにおける運動予測処理の概念を利用して、任意の確率分布関数（提案分布）から生成された乱数を状態ベクトルとする。提案分布において生成された乱数の確率と、等速直進運動における遷移確率とを比較した結果に基づいて、予測後の重み係数を算出してもよい。

　ステップＳＴ７ａにおいて、データ更新部４３は、存在検出センサ３－１～３－Ｍから受信した存在検出データおよびセンサ位置、および、時刻ｋ－１から時刻ｋまでの間の観測対象候補の移動が反映された観測対象候補の存在確率の予測値および運動諸元確率分布の予測値に基づいて、時刻ｋにおける観測対象候補の存在確率および運動諸元確率分布を算出（更新）する。このとき、データ更新部４３は、複数の観測対象候補と複数の存在検出データとの複数通りの対応付けが反映された運動諸元確率分布を算出する。

　図７は、図６のステップＳＴ７の詳細な処理を示すフローチャートである。
　ステップＳＴ７ａ－１において、データ更新部４３は、識別番号ｉの観測対象候補の存在確率の予測値ｒ_{ｋ｜ｋ－１} ^（ｉ）のそれぞれの値に基づいて、識別番号の部分集合ごとの信頼度ｂ_{ｋ｜ｋ－１}（Ｉ^（λ））を算出する。λ＝１，２，・・・，Λである。
　識別番号の部分集合Ｉ^（λ）とは、０個以上の識別番号ｉから構成される集合である。
　例えば、ｉ＝１，２である場合、識別番号の部分集合Ｉ^（λ）は、λ＝１，２，３，４として、Ｉ^（１）＝φ、Ｉ^（２）＝｛１｝、Ｉ^（３）＝｛２｝、Ｉ^（４）＝｛１，２｝の４通りとなる。ここで、φは空集合を表す。

　データ更新部４３は、下記式（２２）を用いて、部分集合ごとの信頼度ｂ_{ｋ｜ｋ－１}（Ｉ^（λ））を算出する。下記式（２２）において、Ｌ_{ｋ｜ｋ－１}は、予測された観測対象候補の識別番号の集合である。下記式（２２）の右辺における、Ｌ_{ｋ｜ｋ－１}＼Ｉ^（λ）は、集合Ｌ_{ｋ｜ｋ－１}から集合Ｉ^（λ）の要素を除いた差集合を示している。

　ステップＳＴ７ａ－２において、データ更新部４３は、存在検出センサ３－１～３－Ｍのうち、時刻ｋにおいて観測対象候補を検出した存在検出センサが設置されたセンサ位置｛ｚ_ｋ，ｍ ^{（Ｅｘｉｓｔ）}｝_ｍ∈Ｍ＊、部分集合信頼度ｂ_{ｋ｜ｋ－１}（Ｉ^（λ））および識別番号ｉの運動諸元確率分布の予測値｛ｗ_{ｋ｜ｋ－１} ^{（ｉ，ｊ）}，ｘ_{ｋ｜ｋ－１} ^{（ｉ，ｊ）}｝_{ｊ＝１，・・・，Ｊ}に基づいて、写像ごとの信頼度ｂ’_ｋ｜ｋ（θ_ｋ ^{（λ，ｈ）}）および写像ごとの識別番号ｉの運動諸元確率分布｛ｗ’_ｋ｜ｋ ^{（ｉ，ｊ）}（θ_ｋ ^{（λ，ｈ）}），ｘ_ｋ｜ｋ ^{（ｉ，ｊ）}（θ_ｋ ^{（λ，ｈ）}）｝_{ｊ＝１，・・・，Ｊ}を算出する。なお、Ｍ^＊は、時刻ｋにおいて観測対象候補を検出した存在検出センサに付与された識別番号（以下、存在検出センサ番号と記載する）の集合である。
　また、信頼度および重み係数の記号に付けられたプライム（’）は、規格化される前の値であることを示している。規格化は、後述するステップＳＴ７ａ－３で説明する。

　写像θ_ｋ ^{（λ，ｈ）}は、部分集合Ｉ^（λ）の要素を、時刻ｋにおいて観測対象候補を検出した存在検出センサの存在検出センサ番号の少なくとも１つに対応付ける写像と定義する。これを数式で表現すると、下記式（２３）のようになる。下記式（２３）において、集合｛０｝への写像は、時刻ｋにおいて、観測対象候補がいずれの存在検出センサにも検出されなかったことを示す写像である。
　存在検出センサは、観測対象候補の個数に関わらず、検出範囲に観測対象候補が存在するか否かを示す存在検出データのみを出力するので、写像θ_ｋ ^{（λ，ｈ）}は必ずしも１対１の写像でなくともよい。
　θ_ｋ ^{（λ，ｈ）}：Ｉ^（λ）→［｛０｝　∪　Ｍ^＊］　　　・・・（２３）

　データ更新部４３は、下記式（２４）、下記式（２４ａ）および下記式（２５）を用いて、写像ごとの信頼度ｂ’_ｋ｜ｋ（θ_ｋ ^{（λ，ｈ）}）および写像ごとの識別番号ｉの運動諸元確率分布｛ｗ’_ｋ｜ｋ ^{（ｉ，ｊ）}（θ_ｋ ^{（λ，ｈ）}），ｘ_ｋ｜ｋ ^{（ｉ，ｊ）}（θ_ｋ ^{（λ，ｈ）}）｝_{ｊ＝１，・・・，Ｊ}を算出する。

　上記式（２４）、上記式（２４ａ）および上記式（２５）における変数は、下記式（２６）、下記式（２７）、下記式（２８）および下記式（２９）に従って算出される。
　下記式（２７）におけるｐ_ｋ，ＤＴは、時刻ｋにおいて存在検出センサが検出失敗せずに観測対象候補を検出する確率である。β_{ｋ，ＦＬＳ} ^{（Ｅｘｉｓｔ）}は、時刻ｋにおいて存在検出センサが観測対象候補を誤検出する確率である。ｄ_ＤＴは、存在検出センサの検出範囲の半径である。πは円周率である。

　ステップＳＴ７ａ－３において、データ更新部４３は、ステップＳＴ７ａ－２で算出した、写像ごとの信頼度ｂ’_ｋ｜ｋ（θ_ｋ ^{（λ，ｈ）}）を規格化し、写像ごとの識別番号ｉの運動諸元確率分布における重み係数｛ｗ’_ｋ｜ｋ ^{（ｉ，ｊ）}（θ_ｋ ^{（λ，ｈ）}）｝_{ｊ＝１，・・・，Ｊ}を規格化する。

　データ更新部４３は、下記式（３０）および下記式（３１）を用いて、信頼度ｂ’_ｋ｜ｋ（θ_ｋ ^{（λ，ｈ）}）を規格化し、写像ごとの識別番号ｉの運動諸元確率分布における重み係数｛ｗ’_ｋ｜ｋ ^{（ｉ，ｊ）}（θ_ｋ ^{（λ，ｈ）}）｝_{ｊ＝１，・・・，Ｊ}を規格化して、信頼度ｂ_ｋ｜ｋ（θ_ｋ ^{（λ，ｈ）}）および重み係数｛ｗ_ｋ｜ｋ ^{（ｉ，ｊ）}（θ_ｋ ^{（λ，ｈ）}）｝_{ｊ＝１，・・・，Ｊ}を得る。
　なお、Ｈ_λは、λ番目の部分集合Ｉ^（λ）における写像の総数である。

　ステップＳＴ７ａ－４において、データ更新部４３は、規格化された写像ごとの信頼度ｂ_ｋ｜ｋ（θ_ｋ ^{（λ，ｈ）}）および重み係数が規格化された写像ごとの識別番号ｉの運動諸元確率分布｛ｗ_ｋ｜ｋ ^{（ｉ，ｊ）}（θ_ｋ ^{（λ，ｈ）}），ｘ_ｋ｜ｋ ^{（ｉ，ｊ）}（θ_ｋ ^{（λ，ｈ）}）｝_{ｊ＝１，・・・，Ｊ}に基づいて、時刻ｋにおける識別番号ｉの観測対象候補の存在確率ｒ_ｋ｜ｋ ^（ｉ）および識別番号ｉの観測対象候補の運動諸元確率分布｛ｗ_ｋ｜ｋ ^{（ｉ，ｊ）}，ｘ_ｋ｜ｋ ^{（ｉ，ｊ）}｝_{ｊ＝１，・・・，Ｊ}を算出する。

　データ更新部４３は、下記式（３２）から下記式（３４）までを用いて、時刻ｋにおける識別番号ｉの観測対象候補の存在確率ｒ_ｋ｜ｋ ^（ｉ）および識別番号ｉの観測対象候補の運動諸元確率分布｛ｗ_ｋ｜ｋ ^{（ｉ，ｊ）}，ｘ_ｋ｜ｋ ^{（ｉ，ｊ）}｝_{ｊ＝１，・・・，Ｊ}を算出する。
　下記式（３３）の“ｆｏｒ　λ，ｈ　ｓ．ｔ．ｉ∈Ｉ^（λ），ｈ＝１，・・・，Ｈ_λ”という条件は、“λは要素ｉを持つＩ^（λ）を与え、ｈは１以上Ｈ_λ以下の整数であること”を示している。

　図６の説明に戻る。
　ステップＳＴ８ａにおいて、推定部４４は、時刻ｋ（現時刻）における観測対象候補のうちのいずれか１つを選択して入力する。これにより、観測対象候補を選択して処理するループが開始される。以降では、識別番号ｉの観測対象候補が選択された場合を例に挙げて説明する。

　次に、推定部４４は、下記式（３５）に従い、ステップＳＴ８ａで選択した識別番号ｉの観測対象候補の存在確率ｒ_ｋ｜ｋ ^（ｉ）が閾値ｒ_Ｔｈよりも大きいか否かを判定する（ステップＳＴ９ａ）。観測対象候補の存在確率の閾値ｒ_Ｔｈは、事前設定パラメータである。
　ｒ_ｋ｜ｋ ^（ｉ）＞ｒ_Ｔｈ　　　　・・・（３５）

　識別番号ｉの観測対象候補の存在確率ｒ_ｋ｜ｋ ^（ｉ）が閾値ｒ_Ｔｈよりも大きいと判定された場合（ステップＳＴ９ａ；ＹＥＳ）、ステップＳＴ１０ａの処理に移行する。
　一方、識別番号ｉの観測対象候補の存在確率ｒ_ｋ｜ｋ ^（ｉ）が閾値ｒ_Ｔｈ以下であれば（ステップＳＴ９ａ；ＮＯ）、ステップＳＴ１１ａの処理に移行する。

　ステップＳＴ１０ａにおいて、推定部４４は、時刻ｋにおける、識別番号ｉの観測対象候補の運動諸元確率分布に基づいて、当該観測対象候補に対応する観測対象の運動諸元を推定する。例えば、推定部４４は、識別番号ｉの観測対象候補に対応する観測対象の運動諸元推定値μ_ｋ ^（ｉ）を、下記式（３６）に従って算出する。

　移動経路保存部５には、異なる時刻フレームに保存された識別番号が同一の運動諸元推定値同士が時系列順に保存される。これにより、運動諸元推定値は、観測対象の推定された移動経路として移動経路保存部５に保存される。
　なお、移動経路保存部５には、観測対象の位置以外の運動諸元も保存され、移動経路保存部５に保存されたデータは、利用者からの要求に応じて抽出される。

　ステップＳＴ１１ａにおいて、推定部４４は、時刻ｋにおける、全ての観測対象候補を選択したか否かを確認する。このとき、未選択の観測対象候補があれば（ステップＳＴ１１ａ；ＮＯ）、ステップＳＴ９ａの処理に戻る。一方、全ての観測対象候補が選択されていた場合（ステップＳＴ１１ａ；ＹＥＳ）、図６の一連の処理が終了される。

　なお、運動諸元確率分布を有限個の粒子によって表す場合、粒子の個数が少ないほど、運動諸元確率分布の微小な勾配を表現することができず、粒子の分布が過剰に偏りやすくなることが知られている。過剰な粒子の偏りを軽減するため、図６に示した現時刻（時刻ｋ）における全てのステップの処理が完了した後に、観測対象候補の運動諸元確率分布の粒子を再構成する処理を追加してもよい。粒子を再構成する方法としては、リサンプリング、または、マルコフ連鎖モンテカルロ（ＭＣＭＣ）ムーブといった手法がある。

　以上のように、実施の形態１に係る移動経路推定装置４は、観測対象候補の消失および出現が反映された観測対象候補の存在確率の予測値および運動諸元確率分布の予測値に基づいて、観測時刻における観測対象候補の存在確率および運動諸元確率分布を更新する。移動経路推定装置４は、更新した観測対象候補の存在確率および運動諸元確率分布に基づいて、観測対象の移動経路を示す運動諸元推定値を算出する。このように、移動経路推定装置４は、観測対象候補の消失および出現を予測して、個数が未知の観測対象の移動経路を推定することができる。
　移動経路推定装置４が観測対象候補の消失および出現を予測することで、観測対象候補の個数の範囲が絞り込まれる。これにより、観測領域に真に存在する観測対象候補を判定するときの演算規模は、絞り込まれた観測対象候補の個数の範囲に依存する。
　例えば、非特許文献１に記載される技術を適用した移動経路推定装置を並列に接続した場合、移動経路の推定処理に要する演算規模は、観測対象候補の最大個数に依存することになる。これに対し、移動経路推定装置４では、絞り込まれた観測対象候補の個数の範囲に依存した演算規模になるため、個数が未知の観測対象の移動経路をより少ない演算規模で推定することが可能である。

　また、実施の形態１に係る移動経路推定装置４において、消失予測部４０が、入退検出データ、センサ位置、前時刻における観測対象候補の存在確率および運動諸元確率分布に基づいて、前時刻から現時刻までの間の観測対象候補の消失が反映された観測対象候補の存在確率の予測値および運動諸元確率分布の予測値を算出する。出現予測部４１が、入退検出データおよびセンサ位置に基づいて、前時刻から現時刻までの間に出現した観測対象候補の存在確率の予測値および運動諸元確率分布の予測値を算出する。
　このように、移動経路推定装置４は、前時刻よりも過去のデータを用いることなく、前時刻における処理結果、現時刻における入退検出データ、現時刻における存在検出データおよびセンサ位置に基づいて、現時刻における移動経路の推定が可能である。
　このように構成することで、消失予測部４０が、前時刻から現時刻までの間の観測対象候補の消失が反映された予測値を算出したときに、予測値の算出に用いた前時刻における観測対象候補の存在確率および運動諸元確率分布を削除することができるので、記憶装置１０２の使用容量を節約することができる。
　同様に、出現予測部４１が、前時刻から現時刻までの間に出現した観測対象候補の予測値を算出したときに、予測値の算出に用いた前時刻における観測対象候補の存在確率および運動諸元確率分布を削除することができるので、記憶装置１０２の使用容量を節約することができる。

　さらに、実施の形態１に係る移動経路推定装置４において、消失予測部４０は、前時刻から現時刻までの間の観測対象候補の消失が反映された観測対象候補の存在確率の予測値を、入退検出センサが観測対象候補を誤検出する確率β_{ｋ，ＦＬＳ} ^{（ＩｎＯｕｔ）}を用いて算出する。出現予測部４１は、現時刻に出現した観測対象候補の存在確率の予測値を、入退検出センサが観測対象候補を誤検出する確率β_{ｋ，ＦＬＳ} ^{（ＩｎＯｕｔ）}を用いて算出する。
　このように構成することで、消失予測部４０および出現予測部４１は、入退検出センサによって観測対象候補の観測領域からの退出または観測領域への進入が誤検出された場合であっても、観測対象候補が誤検出である可能性を含む存在確率の予測値を算出することができる。これにより、実際には存在しないが誤検出された観測対象候補について、移動経路が推定されることを防ぐことができる。

　さらに、実施の形態１に係る移動経路推定装置４において、データ更新部４３は、現時刻における観測対象候補の存在確率および運動諸元確率分布を、存在検出センサが観測対象候補を誤検出する確率β_{ｋ，ＦＬＳ} ^{（Ｅｘｉｓｔ）}、および存在検出センサが検出に失敗せずに観測対象候補を検出する確率ｐ_ｋ，ＤＴを用いて算出する。
　データ更新部４３は、予測値と整合する存在検出データが検出された場合であっても、存在検出データが誤検出である可能性を含む存在確率の予測値および運動諸元確率分布の予測値を算出する。これにより、推定した移動経路が誤検出された観測対象候補の方向にずれることが防止される。
　また、データ更新部４３が確率ｐ_ｋ，ＤＴを用いることによって、予測値に整合する存在検出データが検出されなかった場合であっても、観測対象候補の検出失敗の可能性を含む存在確率および運動諸元確率分布が算出される。これにより、推定した移動経路が、部分的に途切れる、または別の観測対象候補が検出された方向にずれることが防止される。

　さらに、実施の形態１に係る移動経路推定装置４において、移動予測部４２が、消失予測部４０および出現予測部４１によって算出された観測対象候補の存在確率の予測値および運動諸元確率分布の予測値に基づいて、前時刻から現時刻までの間の観測対象候補の移動が反映された観測対象候補の存在確率の予測値および運動諸元確率分布の予測値を算出する。このように、移動予測部４２が、観測対象候補の移動に関する予測モデルの誤差を含む運動諸元確率分布の予測値を算出するので、実際の観測対象候補の移動が予測モデルから外れた場合であっても、運動諸元確率分布の予測誤差を考慮することで、移動経路の推定誤差の増加を抑えることができる。

　例えば、等速直進運動の予測モデルに対して、実際の観測対象候補が停止または旋回を行った場合であっても、移動予測部４２は、予測モデルの誤差に応じた広がりをもつ確率密度関数から粒子を算出する。これにより、観測対象候補の運動諸元確率分布の予測値に対して予測モデルの誤差を反映させることができる。データ更新部４３は、等速直進運動の予測値から外れた位置で観測対象候補が検出されても、予測値から外れた粒子の重み係数を増加させて観測対象候補の存在確率および運動諸元確率分布を更新する。このため、推定部４４は、予測モデルの誤差が修正された移動経路の推定値（運動諸元推定値）を算出することができる。

　さらに、実施の形態１に係る移動経路推定装置４において、データ更新部４３が、複数の観測対象候補と複数の存在検出データとが対応付けられる可能性が反映された運動諸元確率分布を算出する。これにより、観測領域内の観測対象候補同士が接近または移動経路が交差する場合であっても、観測対象候補の運動諸元確率分布において、存在検出データと観測対象候補との整合性が低い対応付けには低い確率となる。このため、観測対象候補同士の移動経路が誤って入れ替わる頻度が低減する。

　実施の形態１に係る移動経路推定システム１は、入退検出センサ２－１～２－Ｎおよび存在検出センサ３－１～３－Ｍを備える。入退検出センサ２－１～２－Ｎは、観測領域内の検出範囲を通過した観測対象候補が通過方向ごとに検出された入退検出データを移動経路推定装置４に送信し、存在検出センサ３－１～３－Ｍは、検出範囲に存在する観測対象候補が検出された存在検出データを移動経路推定装置４に送信する。このように構成することで、センサと移動経路推定装置との通信容量の増加を抑えることができる。

　例えば、入退検出センサ２－１～２－Ｎがカメラであり、存在検出センサ３－１～３－Ｍが測距センサである場合、入退検出センサ２－１～２－Ｎは、カメラで撮影された映像を送信する代わりに、検出範囲を通過した観測対象候補が通過方向ごとに検出された入退検出データを送信すればよい。存在検出センサ３－１～３－Ｍは、測距データを送信する代わりに、検出範囲に存在する観測対象候補が検出された存在検出データを送信すればよい。これにより、センサの検出データをより少ない通信容量で移動経路推定装置４に送信することができる。また、個人情報保護の観点からセンサとして撮影装置を設置することができない場合であっても、移動経路推定システム１は、顔、服装といった個人情報を検出せずに、観測対象の移動経路を推定することが可能である。

　また、入退検出センサ２－１～２－Ｎおよび存在検出センサ３－１～３－Ｍは、入退検出データおよび存在検出データを検出可能なセンサであればよいため、高機能なセンサでなくてもよい。また、観測対象にセンサを設ける必要がない。
　例えば、データ更新部４３は、存在検出データの個数が多いほど、観測対象候補の運動諸元を絞り込むことができる。このため、移動経路推定システム１は、低精度の存在検出センサを多数備えることで、高精度のセンサを多数備えたシステムと同程度の推定精度を実現することができる。

実施の形態２．
　実施の形態１では、入退検出センサによって観測対象候補の個数が検出されないため、入退検出センサによる１回の観測（１時刻フレームでの観測）で誤った個数の観測対象候補が検出される場合がある。例えば、入退検出センサは、観測時間よりも短い時間で複数の観測対象候補が検出範囲に移動してくると、複数の観測対象候補のそれぞれを正確に検出できずに、１つの観測対象候補と誤検出されることがある。この場合、観測対象候補の個数が過小であるため、観測対象の移動経路の推定精度が劣化する。
　そこで、実施の形態２に係る移動経路推定装置は、観測対象候補の個数を考慮した移動経路の推定を行う。

　図８は、この発明の実施の形態２に係る移動経路推定システム１Ａの構成を示すブロック図である。図８において、図１と同一構成要素には同一符号を付して説明を省略する。移動経路推定システム１Ａは、Ｎ（Ｎは２以上の整数）個の入退検出センサ２Ａ－１～２Ａ－Ｎ、Ｍ（Ｍは２以上の整数）個の存在検出センサ３－１～３－Ｍ、移動経路推定装置４Ａ、および移動経路保存部５を備える。移動経路推定装置４Ａは、消失予測部４０Ａ、出現予測部４１Ａ、移動予測部４２、データ更新部４３および推定部４４を備える。

　入退検出センサ２Ａ－１～２Ａ－Ｎは、観測領域に設けられて、入退検出データを検出するセンサであり、光学カメラ、赤外線センサ、および測距センサなどのセンサによって実現される。入退検出データには、観測領域内の検出範囲を通過した観測対象候補が通過方向ごとに検出されたデータの他に、個数データも含まれる。個数データは、入退検出センサによる１回の観測で検出範囲を通過した観測対象候補の個数が検出されたデータである。例えば、入退検出センサは、測距データと閾値との比較結果に基づいて、検出範囲に存在した観測対象候補の個数を検出してもよく、検出範囲の撮影画像を画像解析して観測対象候補の個数を検出してもよい。

　消失予測部４０Ａは、個数データを含む入退検出データ、入退検出センサ位置、および前時刻における観測対象候補の存在確率および運動諸元確率分布に基づいて、前時刻から現時刻までの間の観測対象候補の消失が反映された観測対象候補の存在確率の予測値および運動諸元確率分布の予測値を算出する。
　出現予測部４１Ａは、個数データを含む入退検出データおよび入退検出センサ位置に基づいて、前時刻から現時刻までの間に出現した観測対象候補の存在確率の予測値および運動諸元確率分布の予測値を算出する。

　入退検出センサ２Ａ－１～２Ａ－Ｎは、図４Ａおよび図４Ｂに示したセンサ群１００によって実現される。また、移動経路推定装置４Ａにおける、消失予測部４０Ａ、出現予測部４１Ａ、移動予測部４２、データ更新部４３、および推定部４４のそれぞれの機能は、処理回路によって実現される。すなわち、移動経路推定装置４Ａは、図９を用いて後述するステップＳＴ１ｂからステップＳＴ１１ｂまでの処理を実行するための処理回路を備える。処理回路は、専用のハードウェアの処理回路１０３であってもよいが、メモリ１０５に記憶されたプログラムを実行するプロセッサ１０４であってもよい。

　次に動作について説明する。
　図９は、移動経路推定装置４Ａの動作を示すフローチャートである。図９のステップＳＴ１ｂ、ステップＳＴ２ｂ、ステップＳＴ５ｂからステップＳＴ１１ｂの処理は、図６のステップＳＴ１ａ、ステップＳＴ２ａ、ステップＳＴ５ａからステップＳＴ１１ａの処理と同様であるので説明を省略する。以降の説明においても、図６と同様に、運動諸元確率分布ｆ_ｋ｜ｋ ^（ｉ）（ｘ）を、重み付けされた粒子の集合で近似することを前提とする。
　記号および添え字のうち、説明がないものは、実施の形態１と同一の意味である。
　以降では、現時刻にｎ番目の入退検出センサ２Ａ－ｎによって検出された観測対象候補のうち、検出範囲に進入する方向に移動してきた観測対象候補の個数をν_ｎ ^＋とし、検出範囲から退出する方向に移動した観測対象候補の個数をν_ｎ ^－とする。

　ステップＳＴ３ｂにおいて、消失予測部４０Ａは、ｎ番目の入退検出センサ２Ａ－ｎによって検出された観測対象候補の数（個数データ）を含む入退検出データ、入退検出センサ２Ａ－ｎが設置されたセンサ位置、および時刻ｋ－１における観測対象候補の存在確率および運動諸元確率分布に基づいて、時刻ｋ－１から時刻ｋまでの間の観測対象候補の消失が反映された観測対象候補の存在確率の予測値および運動諸元確率分布の予測値を算出する。

　なお、消失予測部４０Ａが、ステップＳＴ３ｂの処理を１回目に行う場合、時刻ｋ－１における観測対象候補の存在確率および運動諸元確率分布を用いる。２回目以降にステップＳＴ３ｂの処理を行う場合、消失予測部４０Ａは、直前のステップＳＴ３ｂで算出した存在確率の予測値および運動諸元確率分布の予測値を用いる。

　消失予測部４０Ａは、下記式（３７）および上記式（７）から下記式（１１）までを用いて、時刻ｋ－１から時刻ｋまでの観測対象候補の消失が反映された、識別番号ｉの観測対象候補の存在確率の予測値ｒ_{ｋ｜ｋ－１，Ｓ} ^（ｉ）を算出し、時刻ｋ－１から時刻ｋまでの観測対象候補の消失が反映された識別番号ｉの観測対象候補の運動諸元確率分布の予測値を表す粒子である｛ｗ_{ｋ｜ｋ－１，Ｓ} ^{（ｉ，ｊ）}，ｘ_{ｋ｜ｋ－１，Ｓ} ^{（ｉ，ｊ）}｝_{ｊ＝１，・・・，Ｊ}を算出する。
　なお、識別番号ｉは、時刻ｋ－１（前時刻）の観測対象候補と同一であるものとする。｜｜ｘ｜｜_２は、ベクトルｘのユークリッドノルムであり、β_{ｋ，ＦＬＳ} ^{（ＩｎＯｕｔ）}は事前に設定されたパラメータであり、入退検出センサが誤検出する確率を示している。
　ｄ_Ｔｈ（ν_ｎ ^－）は、上記式（９）に示したｄ^{（ｉ，ｎ）}を識別番号ｉに関して昇順に並べたときのν_ｎ ^－番目の値である。観測対象候補の数がν_ｎ ^－よりも小さい場合、ｄ_Ｔｈ（ν_ｎ ^－）はｄ^{（ｉ，ｎ）}の最大値となる。

　上記式（７）から上記式（１１）と上記式（３７）を用いて観測対象候補の存在確率の予測値および運動諸元確率分布の予測値を算出する処理は、識別番号ｉの観測対象候補の運動諸元確率分布の平均位置Ｈ_Ｐｏｓμ_{ｋ－１｜ｋ－１} ^（ｉ）と、観測対象候補の退出を検出したｎ番目の入退検出センサ２Ａ－ｎが設置されたセンサ位置ｚ_ｎ ^{（ＩｎＯｕｔ）}との近さに関する、上位ν_ｎ ^－個の識別番号ｉの観測対象候補の存在確率をβ_{ｋ，ＦＬＳ} ^{（ＩｎＯｕｔ）}倍に減少させる処理に相当する。この処理により、入退検出センサ２Ａ－ｎによって検出された観測対象候補の退出を、観測対象候補の予測値に反映させることができる。

　これまで、識別番号ｉの観測対象候補の存在確率をβ_{ｋ，ＦＬＳ} ^{（ＩｎＯｕｔ）}倍することで、識別番号ｉの観測対象候補の存在確率に対して観測対象候補の退出（消失）が反映されたが、これに限定されるものではない。
　例えば、β_{ｋ，ＦＬＳ} ^{（ＩｎＯｕｔ）}＝０とした場合、上記式（９）で示すｄ^{（ｉ，ｎ）}が閾値よりも小さい場合、消失予測部４０Ａは、識別番号ｉの観測値候補の存在確率および運動諸元確率分布（粒子）を削除してもよい。

　なお、運動諸元確率分布の平均位置のみに基づいて、存在確率を減らす観測対象候補が選択される場合を示したが、これに限定されるものではない。
　例えば、消失予測部４０Ａは、運動諸元確率分布の平均速度ベクトルの方位角度と、入退検出センサが、観測領域から退出する観測対象候補を検出する方向の方位角度と、の差分を用いて、存在確率を減らす観測対象候補を選択してもよい。

　ステップＳＴ４ｂにおいて、出現予測部４１Ａは、入退検出センサ２Ａ－ｎによって検出された観測対象候補の数（個数データ）を含む入退検出データおよび入退検出センサ２Ａ－ｎのセンサ位置に基づいて、時刻ｋ－１から時刻ｋまでの間に出現した観測対象候補の存在確率の予測値および運動諸元確率分布の予測値を算出する。

　出現予測部４１Ａは、上記式（１２）から上記式（１６）までを用いて、時刻ｋ－１から時刻ｋまでに出現した観測対象候補の存在確率の予測値｛ｒ_{ｋ｜ｋ－１，Ｂ} ^（ｉｎ）｝_{ｉｎ∈Ｌｎ}を算出し、時刻ｋ－１から時刻ｋまでに出現した観測対象候補の運動諸元確率分布の予測値を表す粒子である｛ｗ_{ｋ｜ｋ－１，Ｂ} ^{（ｉｎ，ｊ）}，ｘ_{ｋ｜ｋ－１，Ｂ} ^{（ｉｎ，ｊ）}｝_{ｉｎ∈Ｌｎ，ｊ＝１，・・・，Ｊ}を算出する。集合Ｌ_ｎは、ｎ番目の入退検出センサ２Ａ－ｎの周辺に出現したと予測された観測対象候補の識別番号ｉ_ｎを要素とする集合であり、集合Ｌ_ｎの要素のそれぞれは既存の識別番号のいずれとも異なる値である。集合Ｌ_ｎの要素の個数はν_ｎ ^＋個である。

　上記式（１２）から上記式（１６）を用いて、識別番号ｉ_ｎの観測対象候補の存在確率の予測値｛ｒ_{ｋ｜ｋ－１，Ｂ} ^（ｉｎ）｝_{ｉｎ∈Ｌｎ}および識別番号ｉ_ｎの観測対象候補の運動諸元確率分布の予測値｛ｗ_{ｋ｜ｋ－１，Ｂ} ^{（ｉｎ，ｊ）}，ｘ_{ｋ｜ｋ－１，Ｂ} ^{（ｉｎ，ｊ）}｝_{ｉｎ∈Ｌｎ，ｊ＝１，・・・，Ｊ}を算出する処理は、ｎ番目の入退検出センサ２Ａ－ｎの設置位置を初期位置として新規の観測対象候補の粒子を生成し、生成した粒子のそれぞれの速度が一様分布の乱数であるものとして、存在確率と粒子の重み係数を固有値とした観測対象候補を、ν_ｎ ^＋個追加する処理に相当する。この処理により、時刻ｋにおいて入退検出センサ２－ｎによって観測領域への進入が検出された観測対象候補の予測値が得られる。

　なお、新たに出現した観測対象候補に関する粒子の初期速度を、上記式（１５）および上記式（１６）に示した０を中心とする一様分布の乱数としたが、入退検出センサ２Ａ－ｎが、観測領域に進入する観測対象候補を検出する方向に偏った分布の乱数を用いてもよい。

　以上のように、実施の形態２に係る移動経路推定装置４Ａにおいて、消失予測部４０Ａが、観測対象候補の数を含む入退検出データ、入退検出センサ位置、および前時刻における観測対象候補の存在確率および運動諸元確率分布に基づいて、前時刻から現時刻までの間の観測対象候補の消失が反映された観測対象候補の存在確率の予測値および運動諸元確率分布の予測値を算出する。出現予測部は４１Ａが、観測対象候補の数を含む入退検出データおよび入退検出センサ位置に基づいて、前時刻から現時刻までの間に出現した観測対象候補の存在確率の予測値および運動諸元確率分布の予測値を算出する。
　このように構成することで、移動経路推定装置４Ａは、入退検出センサ２Ａ－１～２Ａ－Ｎのそれぞれの１回の観測中に複数の観測対象候補が検出範囲に進入してきても、観測対象候補の個数が過大または過小に推定されることなく移動経路を推定することが可能となる。これにより、移動経路推定装置４Ａは、複数の観測対象候補が互いに異なる経路で移動しても、複数の観測対象候補のそれぞれの移動経路を正確に推定できる。
　また、消失予測部４０Ａが、観測領域から退出する観測対象候補の数分の存在確率を、β_{ｋ，ＦＬＳ} ^{（ＩｎＯｕｔ）}倍に減少させる。これにより、観測領域から退出する観測対象候補を検出した入退検出センサの検出範囲に実際には退出していない別の観測対象候補が存在する場合であっても、観測領域から退出していない観測対象候補が誤って削除される頻度が減少する。また、これにより、観測領域から実際に退出した観測対象候補が複数存在する場合であっても、一部の観測対象候補が削除されずに残存し続ける頻度が減少する。

実施の形態３．
　実施の形態１では、観測領域内の任意の位置に観測対象候補が存在できることを前提として、観測対象候補の移動経路を推定していた。
　しかしながら、実際には、障害物、境界の仕切り壁などによって観測領域内に観測対象候補が進入できない領域が存在する場合がある。この場合、実施の形態１に係る移動経路推定装置では、観測対象候補が進入できない領域への観測対象候補の移動も予測される。
　このため、観測対象候補が進入できない領域が観測領域内に存在しない場合に比べて、観測対象候補の移動が反映された予測値に基づいて推定された運動諸元と真の移動経路を示す運動諸元との間の誤差が増大する可能性がある。
　そこで、実施の形態３に係る移動経路推定装置は、観測対象候補が進入できない領域を考慮した移動経路の推定を行う。

　図１０は、この発明の実施の形態３に係る移動経路推定システム１Ｂの構成を示すブロック図である。図１０において、図１と同一構成要素には同一符号を付して説明を省略する。移動経路推定システム１Ｂは、Ｎ（Ｎは２以上の整数）個の入退検出センサ２－１～２－Ｎ、Ｍ（Ｍは２以上の整数）個の存在検出センサ３－１～３－Ｍ、移動経路推定装置４Ｂ、および移動経路保存部５を備える。移動経路推定装置４Ｂは、消失予測部４０、出現予測部４１、移動予測部４２Ａ、データ更新部４３、推定部４４、およびデータ保存部４５を備える。

　移動予測部４２Ａは、データ保存部４５に保存された領域データ、および、消失予測部４０および出現予測部４１によってそれぞれ算出された観測対象候補の存在確率の予測値および運動諸元確率分布の予測値に基づいて、前時刻から現時刻までの間の観測対象候補の移動が反映された観測対象候補の存在確率の予測値および運動諸元確率分布の予測値を算出する。データ保存部４５には、進入不可領域データが保存される。進入不可領域データは、観測領域のうち、観測対象候補が存在し得ない領域を示すデータである。

　図１１は、実施の形態３における観測領域１０Ａの例を示す図である。図１１に示す観測領域１０Ａには、仕切り壁で囲まれた進入不可領域１４が設けられている。進入不可領域１４には、観測対象候補が進入することができない。この進入不可領域１４を示す進入不可領域データは、Ｇ個のベクトル値γ_ｇから構成されており、ｇ＝１，２，・・・，Ｇである。

　ベクトルγ_ｇのそれぞれは、Ｘ軸に平行な辺とＹ軸に平行な辺とを有する四角形の進入不可領域１４の位置および大きさを示すベクトル値であり、下記式（３８）に示す４つの要素から構成されている。下記式（３８）において、γ_{ｇ，ｍｉｎＸ}は、進入不可領域１４の四隅の位置のうち、最小のＸ座標値を示しており、γ_{ｇ，ｍｉｎＹ}は、進入不可領域１４の四隅の位置のうち、最小のＹ座標値を示している。γ_ｇ，ΔＸは、進入不可領域１４のＸ軸に平行な辺の長さを示しており、γ_ｇ，ΔＹは、進入不可領域１４のＹ軸に平行な辺の長さを示している。すなわち、下記式（３８）では、進入不可領域１４が、１つの四角形の位置と大きさを示すγ_ｇがＧ個組み合わされて表現されている。
　γ_ｇ＝［γ_{ｇ，ｍｉｎＸ}　γ_{ｇ，ｍｉｎＹ}　γ_ｇ，ΔＸ　γ_ｇ，ΔＹ］　　　・・・・（３８）

　また、任意の状態ベクトルｘに対して、下記式（３９）で示す判定関数を定義する。
　ｘ^（１）およびｘ^（２）の定義は上記式（１）と同じである。判定関数Γ_ｇ（ｘ）は、状態ベクトルｘが進入不可領域１４の内部であれば“１”となり、進入不可領域１４の外部であれば“０”となる関数である。

　図１１では、進入不可領域１４の形状を四角形の組み合わせで表現したが、円形、三角形などの四角形以外の形状で表現してもよい。
　また、進入不可領域１４を位置と大きさで表現する代わりに、進入不可領域１４の四隅の位置ベクトルまたは進入不可領域１４の境界線で表現してもよい。

　データ保存部４５は、図４Ａおよび図４Ｂに示した記憶装置１０２により実現される。また、移動経路推定装置４Ｂにおける、消失予測部４０、出現予測部４１、移動予測部４２Ａ、データ更新部４３および推定部４４のそれぞれの機能は、処理回路により実現される。すなわち、移動経路推定装置４Ｂは、図１２を用いて後述するステップＳＴ１ｃからステップＳＴ１１ｃまでの処理を実行するための処理回路を備える。処理回路は、専用のハードウェアの処理回路１０３であってもよいが、メモリ１０５に記憶されたプログラムを実行するプロセッサ１０４であってもよい。
　なお、データ保存部４５は、移動経路推定装置４Ｂとは別に設けられた記憶装置に構築してもよい。この場合、移動予測部４２Ａは、移動経路推定装置４Ｂとは別に設けられた記憶装置と通信して、当該記憶装置から進入不可領域データを取得する。

　次に動作について説明する。
　図１２は、移動経路推定装置４Ｂの動作を示すフローチャートである。図１２のステップＳＴ１ｃからステップＳＴ５ｃおよびステップＳＴ７ｃからステップＳＴ１１ｃの処理は、図６のステップＳＴ１ａからステップＳＴ５ａおよびステップＳＴ７ａからステップＳＴ１１ａの処理と同様であるので説明を省略する。
　以降の説明においても、図６と同様に、運動諸元確率分布ｆ_ｋ｜ｋ ^（ｉ）（ｘ）を、重み付けされた粒子の集合で近似することを前提とする。また、記号および添え字のうち、説明がないものは、実施の形態１と同一の意味である。

　ステップＳＴ６ｃにおいて、移動予測部４２Ａは、進入不可領域データ、および、観測対象候補の消失が反映された存在確率の予測値および運動諸元確率分布の予測値に基づいて、時刻ｋ－１から時刻ｋまで間の観測対象候補の移動が反映された観測対象候補の存在確率の予測値および運動諸元確率分布の予測値を算出する。

　移動予測部４２Ａは、時刻ｋ－１から時刻ｋまでの観測対象候補の消失が反映された、観測対象候補の存在確率の予測値ｒ_{ｋ｜ｋ－１，Ｓ} ^（ｉｓ）と運動諸元確率分布の予測値｛ｗ_{ｋ｜ｋ－１，Ｓ} ^{（ｉｓ，ｊ）}，ｘ_{ｋ｜ｋ－１，Ｓ} ^{（ｉｓ，ｊ）}｝_{ｊ＝１，・・・，Ｊ}、時刻ｋ－１から時刻ｋまでに出現した観測対象候補の存在確率の予測値ｒ_{ｋ｜ｋ－１，Ｂ} ^（ｉＢ）と運動諸元確率分布の予測値｛ｗ_{ｋ｜ｋ－１，Ｂ} ^{（ｉＢ，ｊ）}，ｘ_{ｋ｜ｋ－１，Ｂ} ^{（ｉＢ，ｊ）}｝_{ｊ＝１，・・・，Ｊ}、および、進入不可領域データに基づいて、時刻ｋ－１から時刻ｋまでの間の観測対象候補の移動が反映された、観測対象候補の存在確率の予測値ｒ_{ｋ｜ｋ－１} ^（ｉ）および運動諸元確率分布の予測値｛ｗ_{ｋ｜ｋ－１} ^{（ｉ，ｊ）}，ｘ_{ｋ｜ｋ－１} ^{（ｉ，ｊ）}｝_{ｊ＝１，・・・，Ｊ}を算出する。ｉ_Ｓは、ステップＳＴ３ｃで予測値が算出された観測対象候補の識別番号である。ｉ_Ｂは、ステップＳＴ４ｃで予測値が算出された観測対象候補の識別番号である。識別番号ｉは、ｉ_Ｓとｉ_Ｂとの両方を含む識別番号である。ｉは、上記式（１７）で定義される。

　さらに、移動予測部４２Ａは、上記式（１８）、上記式（２０）、上記式（２１）および下記式（４０）を用いて、観測対象候補の存在確率の予測値ｒ_{ｋ｜ｋ－１} ^（ｉ）および運動諸元確率分布の予測値｛ｗ_{ｋ｜ｋ－１} ^{（ｉ，ｊ）}，ｘ_{ｋ｜ｋ－１} ^{（ｉ，ｊ）}｝_{ｊ＝１，・・・，Ｊ}を算出する。
　まず、移動予測部４２Ａは、上記式（１８）、上記式（２０）および上記式（２１）を用いて、観測対象候補の存在確率の予測値ｒ_{ｋ｜ｋ－１} ^（ｉ）および粒子状態ｘ_{ｋ｜ｋ－１} ^{（ｉ，ｊ）}を算出する。
　続いて、移動予測部４２Ａは、下記式（４０）を用いて、重み係数ｗ_{ｋ｜ｋ－１} ^{（ｉ，ｊ）}を算出する。Ｌ_{ｋ｜ｋ－１，Ｓ}は識別番号ｉ_Ｓの集合であり、Ｌ_{ｋ｜ｋ－１，Ｂ}は識別番号ｉ_Ｂの集合である。

　これまで、ガウス分布に基づいた乱数が生成可能であることを前提として、存在確率の予測値と運動諸元確率分布の予測値とを算出する場合を示したが、これに限定されるものではない。例えば、既存の粒子フィルタにおける運動予測処理の概念を利用して、任意の確率分布関数（提案分布）から生成された乱数を状態ベクトルとする。提案分布において生成された乱数の確率と、等速直進運動における遷移確率とを比較した結果に基づいて、予測後の重み係数を算出してもよい。
　また、乱数によって生成された粒子の状態ベクトルの位置成分が、進入不可領域１４の内部であった場合にその粒子の重み係数を０とした。この他に、進入不可領域１４の内部の粒子が生成されたときにその粒子を削除して、改めて、上記式（２０）および上記式（４０）によって粒子を生成してもよい。

　以上のように、実施の形態３に係る移動経路推定装置４Ｂにおいて、移動予測部４２Ａは、観測領域のうち、観測対象候補が存在し得ない進入不可領域１４を示す進入不可領域データ、および、消失予測部４０および出現予測部４１によってそれぞれ算出された観測対象候補の存在確率の予測値および運動諸元確率分布の予測値に基づいて、前時刻から現時刻までの間の観測対象候補の移動が反映された観測対象候補の存在確率の予測値および運動諸元確率分布の予測値を算出する。
　このように構成することで、移動予測部４２Ａが、観測対象候補の運動諸元確率分布が進入不可領域１４の外側に限定されるように観測対象候補の移動を予測する。
　これにより、観測対象候補の運動諸元確率分布の予測値がとり得る範囲の広がりが抑えられ、この予測値に基づいた移動経路の推定精度が向上する。

　例えば、観測領域が室内であり、室内の歩行者が観測対象である場合に、室内に歩行者が歩く通路が存在すれば、歩行者が時刻フレーム間に移動する範囲は、通路がない室内に比べて絞り込まれる。実施の形態３に係る移動経路推定装置４Ｂは、通路以外の領域を示す進入不可領域データを用いることで、通路内部に絞り込まれた観測対象候補の運動諸元確率分布の予測値を算出することができる。

　また、観測対象候補の運動諸元確率分布の予測値が絞り込まれると、存在検出データが移動経路の推定精度に与える影響が小さくなる。このため、同一の精度で移動経路を推定するときに必要な存在検出センサの個数が削減される。
　例えば、観測領域が室内であり、室内の歩行者が観測対象であり、室内にＴ字の通路が存在する場合、歩行者がＴ字の通路を曲がる動きは、速度および加速度の変動が大きく、等速直進運動または等角速度旋回運動を示す予測モデルを用いた予測が困難である。

　これに対し、実施の形態３に係る移動経路推定装置４Ｂは、通路以外の領域を示す進入不可領域データを用いることで、Ｔ字の通路の内部に絞り込まれた観測対象候補の運動諸元確率分布の予測値を算出する。これにより、運動諸元確率分布の予測値がＴ字の通路の内部に限定されるため、観測対象がＴ字の通路を曲がるといった急激な運動をする場合においても、運動諸元確率分布の予測値と存在検出データとの残差が小さくなり、結果として推定軌道の誤差の増加を抑えることができる。
　実施の形態３に係る移動経路推定装置４Ｂは、実施の形態１に係る移動経路推定装置４よりも少ない存在検出センサの個数で同じ精度の移動経路の推定を行うことができる。

実施の形態４．
　実施の形態１では、観測対象候補が互いに接近しかつ予測モデルから外れた動きである場合、一方の観測対象候補の移動経路と他方の観測対象候補の移動経路が入れ替わる場合がある。例えば、２つの観測対象候補のそれぞれの移動経路が接近し、互いの移動経路が交差せず、かつ観測対象候補が旋回しながら離れる場合、等速直進運動の予測モデルが採用されていると、２つの観測対象候補の移動経路が交差して推定される可能性がある。
　そこで、実施の形態４に係る移動経路推定装置は、観測対象候補の属性を基準に含めた移動経路の推定を行う。

　図１３は、この発明の実施の形態４に係る移動経路推定システム１Ｃの構成を示すブロック図である。図１３において、移動経路推定システム１Ｃは、Ｎ（Ｎは２以上の整数）個の入退検出センサ２Ｂ－１～２Ｂ－Ｎ、Ｍ（Ｍは２以上の整数）個の存在検出センサ３Ａ－１～３Ａ－Ｍ、移動経路推定装置４Ｃおよび移動経路保存部５Ａを備える。移動経路推定装置４Ｃは、消失予測部４０Ｂ、出現予測部４１Ｂ、移動予測部４２Ｂ、データ更新部４３Ａおよび推定部４４Ａを備える。

　入退検出センサ２Ｂ－１～２Ｂ－Ｎは、観測領域に設けられて、入退検出データを検出するセンサであり、光学カメラ、赤外線センサ、および測距センサなどのセンサによって実現される。入退検出データには、観測領域内の検出範囲を通過した観測対象候補が通過方向ごとに検出されたデータの他に、属性検出データも含まれる。属性検出データは、入退検出センサによって検出された観測対象候補の属性を示すデータである。

　観測対象候補の“属性”は、観測対象候補の物体の態様を示す情報であって、例えば、観測対象候補の大きさ、形状、光の反射強度、色彩、周期的な運動、音声、臭気、および付属する物体といった情報が挙げられる。なお、これらの情報は、個々の観測対象候補に固有な情報でなくてもよい。観測対象候補の属性は、数値化された情報（以下、属性値と記載する）で扱われる。入退検出センサ２Ｂ－１～２Ｂ－Ｎは、属性検出データを含む入退検出データおよびセンサ位置を観測時刻ごとに消失予測部４０Ｂおよび出現予測部４１Ｂに送信する。

　存在検出センサ３Ａ－１～３Ａ－Ｍは、観測領域に設けられて、存在検出データを検出するセンサであって、光学カメラ、赤外線センサ、測距センサ、振動センサ、マイクロフォンおよび熱検知センサなどのセンサにより実現される。存在検出データには、検出範囲に存在する観測対象候補が検出されたデータの他に、上記属性検出データも含まれる。
　存在検出センサ３Ａ－１～３Ａ－Ｍは、上記属性検出データを含む存在検出データおよびセンサ位置をデータ更新部４３Ａに送信する。

　消失予測部４０Ｂは、入退検出センサ２Ｂ－１～２Ｂ－Ｎのそれぞれから入力した属性検出データを含む入退検出データおよびセンサ位置、データ更新部４３Ａから入力した前時刻における観測対象候補の存在確率および運動諸元確率分布に基づいて、前時刻から現時刻までの間の観測対象候補の消失が反映された観測対象候補の存在確率の予測値および運動諸元確率分布の予測値と、予測値に対応する観測対象候補の属性を示す属性値を算出する。消失予測部４０Ｂによって算出された観測対象候補の存在確率の予測値、運動諸元確率分布の予測値および属性値は、移動予測部４２Ｂに出力される。

　出現予測部４１Ｂは、入退検出センサ２Ｂ－１～２Ｂ－Ｎのそれぞれから入力した属性検出データを含む入退検出データおよびセンサ位置に基づいて、前記時刻から現時刻までの間に出現した観測対象候補の存在確率の予測値、運動諸元確率分布の予測値および予測値に対応する観測対象候補の属性を示す属性値を算出する。出現予測部４１Ｂによって算出された観測対象候補の存在確率の予測値、運動諸元確率分布の予測値および属性値は、移動予測部４２Ｂに出力される。

　移動予測部４２Ｂは、消失予測部４０Ｂおよび出現予測部４１Ｂによってそれぞれ算出された観測対象候補の存在確率の予測値、運動諸元確率分布の予測値および予測値に対応する観測対象候補の属性を示す属性値に基づいて、前時刻から現時刻までの間の観測対象候補の移動が反映された、観測対象候補の存在確率の予測値、運動諸元確率分布の予測値および予測値に対応する観測対象候補の属性を示す属性値を算出する。移動予測部４２Ｂによって算出された観測対象候補の存在確率の予測値、運動諸元確率分布の予測値および属性値は、データ更新部４３Ａに出力される。

　データ更新部４３Ａは、存在検出センサ３Ａ－１～３Ａ－Ｍのそれぞれから入力した属性検出データを含む存在検出データおよびセンサ位置、および、消失予測部４０Ｂと出現予測部４１Ｂとから入力した観測対象候補の存在確率の予測値、運動諸元確率分布の予測値および予測値に対応する観測対象候補の属性を示す属性値に基づいて、前時刻における観測対象候補の存在確率、運動諸元確率分布および属性値を算出する。
　データ更新部４３Ａは、現時刻における、１つまたは複数の観測対象候補の存在確率、運動諸元確率分布および属性値を算出して、算出した存在確率、運動諸元確率分布および属性値を、推定部４４Ａおよび次の時刻フレームの予測を行う消失予測部４０Ｂに出力する。

　推定部４４Ａは、データ更新部４３Ａから入力した観測対象候補の存在確率、運動諸元確率分布および属性値に基づいて、現時刻に観測対象が観測領域に存在するか否かを推定し、現時刻に存在すると推定した観測対象の運動諸元推定値および属性値を算出する。
　観測対象の運動諸元推定値および属性値は、推定部４４Ａによって移動経路保存部５Ａに出力されて保存される。

　移動経路保存部５Ａに観測時刻ごとに保存された観測対象の運動諸元推定値が時系列に並べられたデータが、観測対象の移動経路を示すデータとなる。移動経路保存部５Ａは、移動経路推定装置４Ｃの利用者からの要求に応じて、観測対象の移動経路を示すデータまたは観測対象の属性を出力する。

　入退検出センサ２Ｂ－１～２Ｂ－Ｎは、図４Ａおよび図４Ｂに示したセンサ群１００によって実現される。存在検出センサ３Ａ－１～３Ａ－Ｍは、センサ群１０１によって実現される。移動経路保存部５Ａは、図４Ａおよび図４Ｂに示した記憶装置１０２によって実現される。
　移動経路推定装置４Ｃにおける、消失予測部４０Ｂ、出現予測部４１Ｂ、移動予測部４２Ｂ、データ更新部４３Ａおよび推定部４４Ａのそれぞれの機能は、処理回路により実現される。すなわち、移動経路推定装置４Ｃは、図１４を用いて後述するステップＳＴ１ｄからステップＳＴ１１ｄまでの処理を実行するための処理回路を備える。
　処理回路は、専用のハードウェアの処理回路１０３であってもよいが、メモリ１０５に記憶されたプログラムを実行するプロセッサ１０４であってもよい。

　次に動作について説明する。
　図１４は、移動経路推定装置４Ｃの動作を示すフローチャートである。図１４のステップＳＴ１ｄ、ステップＳＴ２ｄ、ステップＳＴ５ｄ、ステップＳＴ８ｄ、ステップＳＴ９ｄおよびステップＳＴ１１ｄの処理は、図６のステップＳＴ１ａ、ステップＳＴ２ａ、ステップＳＴ５ａ、ステップＳＴ８ａ、ステップＳＴ９ａおよびステップＳＴ１１ａと同様であるので説明を省略する。
　以降の説明においても、図６と同様に、運動諸元確率分布ｆ_ｋ｜ｋ ^（ｉ）（ｘ）を、重み付けされた粒子の集合で近似することを前提とする。また、記号および添え字のうち、説明がないものは、実施の形態１と同一の意味である。時刻ｋ（現時刻）にｎ（ｎ＝１，２，・・・，Ｎ）番目の入退検出センサ２Ｂ－ｎによって検出された属性検出データを、ｃ_ｋ，ｎ ^{（ＩｎＯｕｔ）}とし、時刻ｋにｍ（ｍ＝１，２，・・・，Ｍ）番目の存在検出センサ３Ａ－ｍによって検出された属性検出データをｃ_ｋ，ｍ ^{（Ｅｘｉｓｔ）}とする。時刻ｋにおける識別番号ｉの観測対象候補に対応する属性値をｃ_ｋ｜ｋ ^（ｉ）とする。なお、属性値をスカラー値としたが、複数の属性値を組み合わせて観測対象候補の属性を表現する場合は、属性値のそれぞれをベクトル値としてもよい。

　ステップＳＴ３ｄにおいて、消失予測部４０Ｂは、ｎ番目の入退検出センサ２Ｂ－ｎによって検出された属性検出データｃ_ｋ，ｎ ^{（ＩｎＯｕｔ）}を含む入退検出データ、センサ位置、および、時刻ｋ－１における観測対象候補の存在確率および運動諸元確率分布に基づいて、前時刻から現時刻までの間の観測対象候補の消失が反映された観測対象候補の存在確率の予測値、運動諸元確率分布の予測値および予測値に対応する観測対象候補の属性値を算出する。

　なお、消失予測部４０Ｂが、ステップＳＴ３ｄの処理を１回目に行う場合、時刻ｋ－１における観測対象候補の存在確率および運動諸元確率分布を用いる。２回目以降にステップＳＴ３ｄの処理を行う場合、消失予測部４０Ａは、直前のステップＳＴ３ｄで算出した存在確率の予測値、運動諸元確率分布の予測値および属性値を用いる。

　例えば、消失予測部４０Ｂは、上記式（７）、上記式（８）、上記式（１１）、下記式（４１）から下記式（４４）までを用いて、時刻ｋ－１から時刻ｋまでの観測対象候補の消失が反映された、観測対象候補の存在確率の予測値ｒ_{ｋ｜ｋ－１，Ｓ} ^（ｉ）、運動諸元確率分布の予測値を示す粒子である｛ｗ_{ｋ｜ｋ－１，Ｓ} ^{（ｉ，Ｊ）}，ｘ_{ｋ｜ｋ－１，Ｓ} ^{（ｉ，Ｊ）}｝_{ｊ＝１，・・・，Ｊ}、および、これらの予測値に対応する観測対象候補の属性を示す属性値ｃ_{ｋ｜ｋ－１，Ｓ} ^（ｉ）を算出する。識別番号ｉは、時刻ｋ－１の観測対象候補と同一である。｜｜ｘ｜｜_２はベクトルｘのユークリッドノルムである。β_{ｋ，ＦＬＳ} ^{（ＩｎＯｕｔ）}およびｄ_Ｔｈは、事前設定パラメータであり、特に、β_{ｋ，ＦＬＳ} ^{（ＩｎＯｕｔ）}は、入退検出センサが誤検出する確率である。

　上記式（７）、上記式（８）、上記式（１１）、上記式（４１）から上記式（４４）を用いて、観測対象候補の存在確率の予測値、運動諸元確率分布の予測値、およびこれらの予測値に対応する観測対象候補の属性値を算出する処理は、識別番号ｉの観測対象候補の運動諸元確率分布の平均位置Ｈ_Ｐｏｓｍ_{ｋ－１｜ｋ－１} ^（ｉ）が、観測対象候補の退出を検出したｎ番目の入退検出センサ２Ｂ－ｎが設置されたセンサ位置ｚ_ｎ ^{（ＩｎＯｕｔ）}に近く、かつ、識別番号ｉの観測対象候補の属性値が示す属性と入退検出センサ２Ｂ－ｎによって検出された属性検出データが示す属性とが同一であれば、識別番号ｉの観測対象候補の存在確率をβ_{ｋ，ＦＬＳ} ^{（ＩｎＯｕｔ）}倍に減少させる処理に相当する。この処理により、入退検出センサ２Ｂ－ｎによって検出された観測対象候補の退出を、観測対象候補の予測値に反映させることができる。

　なお、上記説明では、観測対象候補の存在確率をβ_{ｋ，ＦＬＳ} ^{（ＩｎＯｕｔ）}倍して観測対象候補の退出を反映させる場合を示したが、これに限定されるものではない。
　例えば、β_{ｋ，ＦＬＳ} ^{（ＩｎＯｕｔ）}＝０として、上記式（４３）に示すｄ^{（ｉ，ｎ）}が閾値よりも小さい場合に、消失予測部４０Ｂが、識別番号ｉの観測値候補の存在確率および運動諸元確率分布（粒子）を削除してもよい。

　また、存在確率を減らす観測対象候補を、運動諸元確率分布の平均位置のみに基づいて選択する場合を示したが、これに限定されるものではない。
　例えば、消失予測部４０Ｂは、運動諸元確率分布の平均速度ベクトルの方位角度と、入退検出センサが、観測領域から退出する観測対象候補を検出する方向の方位角度と、の差分を用いて、存在確率を減らす観測対象候補を選択してもよい。

　ステップＳＴ４ｄにおいて、出現予測部４１Ｂは、ｎ番目の入退検出センサ２Ｂ－ｎによって検出された属性検出データｃ_ｋ，ｎ ^{（ＩｎＯｕｔ）}を含む入退検出データ、および上記式（２）で示した入退検出センサ２Ｂ－ｎのセンサ位置に基づいて、前時刻から現時刻までの間に出現した観測対象候補の存在確率の予測値ｒ_{ｋ｜ｋ－１，Ｂ} ^（ｉｎ）、運動諸元確率分布の予測値｛ｗ_{ｋ｜ｋ－１，Ｂ} ^{（ｉｎ，Ｊ）}，ｘ_{ｋ｜ｋ－１，Ｂ} ^{（ｉｎ，Ｊ）}｝_{ｊ＝１，・・・，Ｊ}、および予測値に対応する観測対象候補の属性値ｃ_{ｋ｜ｋ－１，Ｂ} ^（ｉｎ）を算出する。識別番号ｉ_ｎは、ｎ番目の入退検出センサ２Ｂ－ｎの周辺に出現したと予測された観測対象候補の識別番号であり、既存の観測対象候補の識別番号と異なる値とする。

　例えば、出現予測部４１Ｂは、上記式（１２）から上記式（１６）と下記式（４５）を用いて、観測対象候補の存在確率の予測値ｒ_{ｋ｜ｋ－１，Ｂ} ^（ｉｎ）、運動諸元確率分布の予測値｛ｗ_{ｋ｜ｋ－１，Ｂ} ^{（ｉｎ，Ｊ）}，ｘ_{ｋ｜ｋ－１，Ｂ} ^{（ｉｎ，Ｊ）}｝_{ｊ＝１，・・・，Ｊ}および予測値に対応する観測対象候補の属性値ｃ_{ｋ｜ｋ－１，Ｂ} ^（ｉｎ）を算出する。ここで、β_{ｋ，ＦＬＳ} ^{（ＩｎＯｕｔ）}、Ｖ_Ｂ ^（３）およびＶ_Ｂ ^（４）は事前設定パラメータであり、特に、β_{ｋ，ＦＬＳ} ^{（ＩｎＯｕｔ）}は、入退検出センサ２Ｂ－ｎが誤検出する確率である。

　上記式（１２）から上記式（１６）と上記式（４５）を用いて観測対象候補の存在確率の予測値、運動諸元確率分布の予測値および属性値を算出する処理は、ｎ番目の入退検出センサ２Ｂ－ｎの設置位置を初期位置として新規の観測対象候補の粒子を生成し、生成した粒子のそれぞれの速度を一様分布の乱数として、存在確率と粒子の重み係数とを固有値とする識別番号ｉ_ｎの観測対象候補を追加する処理に相当する。この処理により、時刻ｋに入退検出センサ２Ｂ－ｎによって観測領域への進入が検出された観測対象候補の予測値が得られる。

　新たに出現した観測対象候補に関する粒子の初期速度を、上記式（１５）および上記式（１６）に示した０を中心とする一様分布の乱数としたが、ｎ番目の入退検出センサ２Ｂ－ｎが、観測領域に進入する観測対象候補を検出する方向に偏った分布の乱数を用いてもよい。

　ステップＳＴ６ｄにおいて、移動予測部４２Ｂは、時刻ｋ－１から時刻ｋまでの間の観測対象候補の消失が反映された存在確率の予測値および運動諸元確率分布の予測値と予測値に対応する観測対象候補の属性値、および、時刻ｋ－１から時刻ｋまでの間に出現した観測対象候補の存在確率の予測値および運動諸元確率分布の予測値と予測値に対応する観測対象候補の属性値に基づいて、時刻ｋ－１から時刻ｋまでの間の観測対象候補の移動が反映された観測対象候補の存在確率の予測値、運動諸元確率分布の予測値、および予測値に対応する観測対象候補の属性値を算出する。

　移動予測部４２Ｂは、時刻ｋ－１から時刻ｋまで間の観測対象候補の消失が反映された観測対象候補の存在確率の予測値ｒ_{ｋ｜ｋ－１，Ｓ} ^（ｉｓ）、時刻ｋ－１から時刻ｋまで間の観測対象候補の消失が反映された運動諸元確率分布の予測値｛ｗ_{ｋ｜ｋ－１，Ｓ} ^{（ｉｓ，ｊ）}，ｘ_{ｋ｜ｋ－１，Ｓ} ^{（ｉｓ，ｊ）}｝_{ｊ＝１，・・・，Ｊ}、これらの予測値に対応する観測対象候補の属性値ｃ_{ｋ｜ｋ－１，Ｓ} ^（ｉｓ）、時刻ｋ－１から時刻ｋまでに出現した観測対象候補の存在確率の予測値ｒ_{ｋ｜ｋ－１，Ｂ} ^（ｉＢ）、時刻ｋ－１から時刻ｋまでに出現した観測対象候補の運動諸元確率分布の予測値｛ｗ_{ｋ｜ｋ－１，Ｂ} ^{（ｉＢ，ｊ）}，ｘ_{ｋ｜ｋ－１，Ｂ} ^{（ｉＢ，ｊ）}｝_{ｊ＝１，・・・，Ｊ}、および、これらの予測値に対応する観測対象候補の属性値ｃ_{ｋ｜ｋ－１，Ｂ} ^（ｉＢ）に基づいて、時刻ｋ－１から時刻ｋまで間の観測対象候補の移動が反映された、観測対象候補の存在確率の予測値ｒ_{ｋ｜ｋ－１} ^（ｉ）および運動諸元確率分布の予測値｛ｗ_{ｋ｜ｋ－１} ^{（ｉ，ｊ）}，ｘ_{ｋ｜ｋ－１} ^{（ｉ，ｊ）}｝_{ｊ＝１，・・・，Ｊ}とこれらの予測値に対応する観測対象候補の属性値ｃ_{ｋ｜ｋ－１} ^（ｉ）を算出する。
　ここで、ｉ_Ｓは、ステップＳＴ３ｄで予測値が算出された観測対象候補の識別番号である。ｉ_Ｂは、ステップＳＴ４ｄで予測値が算出された観測対象候補の識別番号である。識別番号ｉは、ｉ_Ｓとｉ_Ｂの両方を含む識別番号である。ｉ_Ｓの集合をＬ_{ｋ｜ｋ－１，Ｓ}とし、ｉ_Ｂの集合をＬ_{ｋ｜ｋ－１，Ｂ}とすることで、ｉは、上記式（１７）で定義される。

　移動予測部４２Ｂは、上記式（１８）から上記式（２１）までおよび下記式（４６）を用いて観測対象候補の移動が反映された観測対象候補の存在確率の予測値ｒ_{ｋ｜ｋ－１} ^（ｉ）および観測対象候補の運動諸元確率分布の予測値｛ｗ_{ｋ｜ｋ－１} ^{（ｉ，ｊ）}，ｘ_{ｋ｜ｋ－１} ^{（ｉ，ｊ）}｝_{ｊ＝１，・・・，Ｊ}と属性値ｃ_{ｋ｜ｋ－１} ^（ｉ）を算出する。
　ただし、Φ_ｋは、等速直進運動の予測モデルに基づいて状態ベクトルを遷移させる行列である。τ_ｋは、時刻ｋ－１から時刻ｋまでの経過時間である。Ｑ_ｋは、等速直進運動の予測モデルの誤差を表すシステム雑音を示す行列であって、事前設定パラメータである。

　ステップＳＴ７ｄにおいて、データ更新部４３Ａは、存在検出センサ３Ａ－１～３Ａ－Ｍから受信した属性検出データｃ_ｋ，ｍ ^{（Ｅｘｉｓｔ）}を含む存在検出データおよびセンサ位置、時刻ｋ－１から時刻ｋまでの間の観測対象候補の移動が反映された観測対象候補の存在確率の予測値および運動諸元確率分布の予測値と属性値に基づいて、時刻ｋにおける観測対象候補の存在確率および運動諸元確率分布を算出（更新）する。

　図１５は、図１４のステップＳＴ７ｄの詳細な処理を示すフローチャートである。
　ステップＳＴ７ｄ－１において、データ更新部４３Ａは、上記式（２３）を用いて、識別番号ｉの観測対象候補の存在確率の予測値ｒ_{ｋ｜ｋ－１} ^（ｉ）のそれぞれの値に基づいて、識別番号の部分集合ごとの信頼度ｂ_{ｋ｜ｋ－１}（Ｉ^（λ））を算出する。λ＝１，２，・・・，Λである。

　ステップＳＴ７ｄ－２において、データ更新部４３Ａは、存在検出センサ３Ａ－１～３Ａ－Ｍのうち、時刻ｋにおいて観測対象候補を検出した存在検出センサが設置されたセンサ位置｛ｚ_ｋ，ｍ ^{（Ｅｘｉｓｔ）}｝_ｍ∈Ｍ＊、部分集合信頼度ｂ_{ｋ｜ｋ－１}（Ｉ^（λ））、および識別番号ｉの運動諸元確率分布の予測値｛ｗ_{ｋ｜ｋ－１} ^{（ｉ，ｊ）}，ｘ_{ｋ｜ｋ－１} ^{（ｉ，ｊ）}｝_{ｊ＝１，・・・，Ｊ}に基づいて、写像ごとの信頼度ｂ’_ｋ｜ｋ（θ_ｋ ^{（λ，ｈ）}）および写像ごとの識別番号ｉの運動諸元確率分布｛ｗ_ｋ｜ｋ ^{（ｉ，ｊ）}（θ_ｋ ^{（λ，ｈ）}），ｘ_ｋ｜ｋ ^{（ｉ，ｊ）}（θ_ｋ ^{（λ，ｈ）}）｝_{ｊ＝１，・・・，Ｊ}を算出する。写像θ_ｋ ^{（λ，ｈ）}は、実施の形態１で示したように、部分集合Ｉ^（λ）の要素を、時刻ｋにおいて観測対象候補を検出した存在検出センサの存在検出センサ番号の少なくとも１つに対応付ける写像と定義する。これを数式で表現すると、上記式（２３）のようになる。

　ここで、データ更新部４３Ａは、上記式（２４）、上記式（２４ａ）および上記式（２５）を用いて、写像ごとの信頼度ｂ’_ｋ｜ｋ（θ_ｋ ^{（λ，ｈ）}）および写像ごとの識別番号ｉの運動諸元確率分布｛ｗ_ｋ｜ｋ ^{（ｉ，ｊ）}（θ_ｋ ^{（λ，ｈ）}），ｘ_ｋ｜ｋ ^{（ｉ，ｊ）}（θ_ｋ ^{（λ，ｈ）}）｝_{ｊ＝１，・・・，Ｊ}を算出する。
　上記式（２４）、上記式（２４ａ）および上記式（２５）における変数は、上記式（２６）、上記式（２８）、上記式（２９）および下記式（４７）に従って算出される。
　下記式（４７）におけるｐ_ｋ，ＤＴは、時刻ｋにおいて存在検出センサが検出に失敗せずに観測対象を検出する確率である。πは円周率である。

　ステップＳＴ７ｄ－３において、データ更新部４３Ａは、ステップＳＴ７ｄ－２で算出した、写像ごとの信頼度ｂ’_ｋ｜ｋ（θ_ｋ ^{（λ，ｈ）}）を規格化し、写像ごとの識別番号ｉの運動諸元確率分布における重み係数｛ｗ’_ｋ｜ｋ ^{（ｉ，ｊ）}（θ_ｋ ^{（λ，ｈ）}）｝_{ｊ＝１，・・・，Ｊ}を規格化する。規格化には、上記式（３０）および上記式（３１）が用いられる。

　ステップＳＴ７ｄ－４において、データ更新部４３Ａは、規格化された写像ごとの信頼度ｂ_ｋ｜ｋ（θ_ｋ ^{（λ，ｈ）}）および重み係数が規格化された写像ごとの識別番号ｉの運動諸元確率分布｛ｗ_ｋ｜ｋ ^{（ｉ，ｊ）}（θ_ｋ ^{（λ，ｈ）}），ｘ_ｋ｜ｋ ^{（ｉ，ｊ）}（θ_ｋ ^{（λ，ｈ）}）｝_{ｊ＝１，・・・，Ｊ}に基づいて、時刻ｋにおける識別番号ｉの観測対象候補の存在確率ｒ_ｋ｜ｋ ^（ｉ）、識別番号ｉの観測対象候補の運動諸元確率分布｛ｗ_ｋ｜ｋ ^{（ｉ，ｊ）}，ｘ_ｋ｜ｋ ^{（ｉ，ｊ）}｝_{ｊ＝１，・・・，Ｊ}、および時刻ｋにおける識別番号ｉの観測対象候補の属性値ｃ_ｋ｜ｋ ^（ｉ）を算出する。

　データ更新部４３Ａは、上記式（３２）、上記式（３３）および上記式（３４）を用いて、時刻ｋにおける識別番号ｉの観測対象候補の存在確率ｒ_ｋ｜ｋ ^（ｉ）および識別番号ｉの観測対象候補の運動諸元確率分布｛ｗ_ｋ｜ｋ ^{（ｉ，ｊ）}，ｘ_ｋ｜ｋ ^{（ｉ，ｊ）}｝_{ｊ＝１，・・・，Ｊ}を算出する。データ更新部４３Ａは、下記式（４８）を用いて、時刻ｋにおける識別番号ｉの観測対象候補の属性値ｃ_ｋ｜ｋ ^（ｉ）を算出する。

　図１４の説明に戻る。
　ステップＳＴ１０ｄにおいて、推定部４４Ａは、時刻ｋにおける識別番号ｉの観測対象候補の運動諸元確率分布に基づいて、当該観測対象候補に対応する観測対象の運動諸元と属性値を推定する。例えば、推定部４４は、識別番号ｉの観測対象候補に対応する観測対象の運動諸元推定値μ_ｋ ^（ｉ）を、上記式（３６）に従って算出し、識別番号ｉの観測対象候補に対応する観測対象の属性値の推定値ｃ_ｋ ^（ｉ）を、下記式（４９）に従って算出する。

　移動経路保存部５Ａには、異なる時刻フレームに保存された識別番号が同一の運動諸元推定値同士が時系列順に保存される。これにより、運動諸元推定値は、観測対象の推定された移動経路として移動経路保存部５Ａに保存される。
　なお、移動経路保存部５Ａには、観測対象の位置以外の運動諸元および観測対象の属性値も保存され、移動経路保存部５に保存されたデータは、利用者からの要求に応じて抽出される。

　以上のように、実施の形態４に係る移動経路推定装置４Ｃにおいて、消失予測部４０Ｂおよび出現予測部４１Ｂが、属性検出データを含む入退検出データを用いて、観測対象候補の存在確率の予測値、運動諸元確率分布の予測値および予測値に対応する観測対象候補の属性を示す属性値を算出する。移動予測部４２Ｂが、消失予測部４０Ｂおよび出現予測部４１Ｂからそれぞれ入力した予測値および属性値に基づいて、観測対象候補の存在確率の予測値、運動諸元確率分布の予測値および予測値に対応する観測対象候補の属性値を算出する。データ更新部４３Ａが、属性検出データを含む存在検出データ、存在検出センサ位置、移動予測部４２Ｂによって算出された観測対象候補の存在確率の予測値、運動諸元確率分布の予測値および属性値に基づいて、現時刻における観測対象候補の存在確率、運動諸元確率分布および属性値を算出する。推定部４４Ａが、データ更新部４３Ａによって算出された観測対象候補の存在確率、運動諸元確率分布および属性値に基づいて、観測対象の運動諸元推定値および属性値を算出する。
　このように構成することで、データ更新部４３Ａは、存在検出センサによって検出された属性検出データと観測対象候補の属性値が一致しない場合、観測対象候補の存在確率を減少させるため、異なる属性の観測対象候補間で推定された移動経路が入れ替わる頻度を減らすことができる。例えば、２つの観測対象候補が接近しつつ旋回する場合、接近した観測対象が互いに異なる属性であれば、誤って入れ替わった移動経路が推定されることがない。

実施の形態５．
　図１６は、この発明の実施の形態５に係る移動経路推定システム１Ｄの構成を示すブロック図である。移動経路推定システム１Ｄは、観測領域に設定した複数のノード（拠点）のうちのいずれかのノードに観測対象候補が出現し、観測対象候補がノードから消失することを予測して、観測対象候補が存在したノード間を繋いだ移動経路を推定する。
　図１６に示すように、移動経路推定システム１Ｄは、Ｎ（Ｎは２以上の整数）個の入退検出センサ２Ｃ－１～２Ｃ－Ｎ、Ｍ（Ｍは２以上の整数）個の存在検出センサ３Ｂ－１～３Ｂ－Ｍ、移動経路推定装置４Ｄおよび移動経路保存部５Ｂを備えて構成される。
　移動経路推定装置４Ｄは、消失予測部４０Ｃ、出現予測部４１Ｃ、移動予測部４２Ｃ、データ更新部４３Ｂ、推定部４４Ｂおよびデータ保存部４６を備える。

　入退検出センサ２Ｃ－１～２Ｃ－Ｎは、観測領域への進入または観測領域からの退出が可能な位置に設定された全てのノードに設けられて、入退検出データを検出するセンサである。入退検出センサ２Ｃ－１～２Ｃ－Ｎは、光学カメラ、赤外線センサ、および測距センサなどのセンサによって実現される。入退検出データは、ノードから観測領域外へ退出する観測対象候補および観測領域外からノードへ進入する観測対象候補が検出されたデータである。

　入退検出センサ２Ｃ－１～２Ｃ－Ｎは、入退検出データおよびノード番号を含む検出信号を、観測時刻ごとに消失予測部４０Ｃおよび出現予測部４１Ｃへ送信する。
　ノード番号は、観測領域に設定された全てのノードに付与された通し番号である。
　また、入退検出センサ２Ｃ－１～２Ｃ－Ｎから消失予測部４０Ｃおよび出現予測部４１Ｃに送信されるノード番号は、入退検出データが得られたノードを示す入退検出ノード情報に相当する。

　存在検出センサ３Ｂ－１～３Ｂ－Ｍは、観測領域の内部に設定された全てまたは一部のノードに設けられて、存在検出データを検出するセンサである。存在検出センサ３Ｂ－１～３Ｂ－Ｍは、光学カメラ、赤外線センサ、測距センサ、振動センサ、マイクロフォンおよび熱検知センサなどのセンサによって実現される。また、存在検出データは、センサの検出範囲に存在する観測対象候補が検出されたデータである。

　存在検出センサ３Ｂ－１～３Ｂ－Ｍは、存在検出データおよびノード番号を含む検出信号を、データ更新部４３Ｂに送信する。なお、存在検出センサ３Ｂ－１～３Ｂ－Ｍからデータ更新部４３Ｂに送信されるノード番号は、存在検出データが得られたノードを示す存在検出ノード情報に相当する。入退検出センサ２Ｃ－１～２Ｃ－Ｎおよび存在検出センサ３Ｂ－１～３Ｂ－Ｍによって検出された複数の観測対象候補には、それぞれの観測対象候補に固有な識別番号が付与される。

　入退検出センサ２Ｃ－１～２Ｃ－Ｎおよび存在検出センサ３Ｂ－１～３Ｂ－Ｍが、検出データを移動経路推定装置４Ｄに出力する時間の区切りを“時刻フレーム”と呼び、センサそれぞれの時刻フレームは同期しているものとする。実施の形態１と同様に、センサによって観測対象候補が検出された観測時刻に相当する時刻フレームを“現時刻”と呼び、現時刻の１時刻フレーム前の過去の時刻に相当する時刻フレームを“前時刻”と呼ぶ。
　前時刻を現時刻よりも１時刻フレーム前の時刻と定義したが、これは例示であり、実施の形態５では、現時刻よりも一定時間だけ過去の時刻であれば、２以上の時刻フレームを遡った時刻であってもよい。

　図１７は、実施の形態５における観測領域１０Ｂの例を示す図であり、５つのノード２００～２０５およびノード間経路３００を示している。ノード間経路３００は、ノード間を繋いだ経路であり、観測対象候補は、ノード間経路３００に沿って移動する。ノード間経路３００は、観測領域１０Ｂ内の全てのノード間に存在するとは限らない。図１７の例では、ノード２０３とノード２０５との間にはノード間経路が存在しない。
　なお、実施の形態５においても観測対象候補の個数は未知であり、複数の観測対象候補が同時に移動することがあり、また観測領域から観測対象候補が消失する場合があるものとする。

　観測対象候補となる人Ａ３は、観測領域１０Ｂの外からノード２００に進入し、矢印で示すように、ノード２００とノード２０１との間のノード間経路３００を通って、ノード２０１に到達する。このとき、人Ａ３は、“ノード２０１から移動する”および“ノード２０１に留まる”のうちのいずれかを選択する。図１７の例では、人Ａ３は、矢印で示すように、ノード２０１から、ノード２０３、ノード２０４およびノード２０５の順で移動して、ノード２０５から観測領域１０Ｂの外へ退出する。

　観測対象候補が人である場合を示したが、ノード間を移動する観測対象候補としては、駅間を移動する電車または人、空港領域を移動する航空機、自動車専用道において駐車場間を移動する自動車、港間を移動する船舶、部屋間を移動する人または携帯機器、市街地に点在する拠点（例えば、コンビニエンスストア）の近傍を通過する人または自動車、国土内の都市を移動する人または携帯機器などが挙げられる。

　図１６の説明に戻る。
　移動経路推定装置４Ｄは、入退検出センサ２Ｃ－１～２Ｃ－Ｎから入力した入退検出データ、および存在検出センサ３Ｂ－１～３Ｂ－Ｍから入力した存在検出データを用いて、観測対象候補の運動諸元推定値を算出する。運動諸元推定値は、観測対象候補の移動経路を構成するデータであって、観測対象候補ごとに移動経路保存部５Ｂに保存される。移動経路保存部５Ｂに保存された移動経路を示すデータは、移動経路推定システム１Ｄの利用者からの要求に応じて出力される。

　消失予測部４０Ｃは、入退検出センサ２Ｃ－１～２Ｃ－Ｎから受信した入退検出データとノード番号とを含む検出信号、前時刻における観測対象候補の存在確率および前時刻における観測対象候補のノード間確率分布に基づいて、前時刻から現時刻までの間の観測対象候補の消失が反映された、観測対象候補の存在確率の予測値およびノード間確率分布の予測値を算出する。なお、上記ノード番号は、入退検出データが得られたノードのノード番号である。消失予測部４０Ｃによって算出された存在確率の予測値およびノード間確率分布の予測値は、移動予測部４２Ｃに出力される。

　ノード間確率分布は、観測対象候補が存在する位置の確率分布であって、例えば、観測対象候補がノードに存在する確率が離散的な空間内で定義された確率分布である。
　すなわち、ノード間確率分布は、Ｎ個のスカラー値から構成されており、ノード番号ｎにおけるノード間確率分布は、観測領域のｎ番目のノードに観測対象候補が存在する確率に相当する。

　出現予測部４１Ｃは、入退検出センサ２Ｃ－１～２Ｃ－Ｎから受信した入退検出データおよびノード番号を含む検出信号に基づいて、前時刻から現時刻までの間に出現した観測対象候補の存在確率の予測値およびノード間確率分布の予測値を算出する。
　なお、上記ノード番号は、入退検出データが得られたノードのノード番号である。出現予測部４１Ｃによって算出された存在確率の予測値およびノード間確率分布の予測値は、移動予測部４２Ｃに出力される。

　移動予測部４２Ｃは、消失予測部４０Ｃおよび出現予測部４１Ｃから入力した観測対象候補の存在確率の予測値およびノード間確率分布の予測値、およびデータ保存部４６から取得したノード間経路データに基づいて、前時刻から現時刻までの間の観測対象候補の移動が反映された観測対象候補の存在確率の予測値およびノード間確率分布の予測値を算出する。移動予測部４２Ｃによって算出された存在確率の予測値およびノード間確率分布の予測値は、データ更新部４３Ｂに出力される。

　ノード間経路データは、観測対象候補が移動可能なノード間経路を示すデータである。データ保存部４６には、観測領域１０Ｂで観測対象候補が移動可能なノード間経路を示すノード間経路データが保存される。

　データ更新部４３Ｂは、存在検出センサ３Ｂ－１～３Ｂ－Ｍから受信した存在検出データおよびノード番号を含む検出信号、移動予測部４２Ｃによって算出された観測対象候補の存在確率の予測値および運動諸元確率分布の予測値に基づいて、現時刻における観測対象候補の存在確率および運動諸元確率分布を算出する。上記ノード番号は、存在検出データが得られたノードのノード番号である。データ更新部４３Ｂは、現時刻における１つまたは複数の観測対象候補の存在確率およびノード間確率分布を算出し、算出した存在確率およびノード間確率分布を、推定部４４Ｂおよび次の時刻フレームの予測を行う消失予測部４０Ｃに出力する。

　推定部４４Ｂは、データ更新部４３Ｂから入力した現時刻における観測対象候補の存在確率およびノード間確率分布に基づいて、現時刻に観測対象が観測領域に存在するか否かを推定し、現時刻に存在すると推定された観測対象の運動諸元推定値を算出する。
　推定部４４Ｂによって算出された運動諸元推定値は、移動経路保存部５Ｂに出力されて保存される。移動経路保存部５Ｂに保存された観測対象の運動諸元推定値が時系列に並べられたデータが、観測対象の移動経路を示すデータとなる。
　なお、移動経路保存部５Ｂに保存された観測対象の移動経路を示すデータは、移動経路推定装置４Ｄの利用者の要求に応じて表示装置などに出力される。

　入退検出センサ２Ｃ－１～２Ｃ－Ｎは、図４Ａおよび図４Ｂに示したセンサ群１００によって実現される。存在検出センサ３Ｂ－１～３Ｂ－Ｍは、センサ群１０１によって実現される。さらに、移動経路保存部５Ｂおよびデータ保存部４６は、図４Ａおよび図４Ｂに示した記憶装置１０２によって実現される。
　移動経路推定装置４Ｄにおける、消失予測部４０Ｃ、出現予測部４１Ｃ、移動予測部４２Ｃ、データ更新部４３Ｂおよび推定部４４Ｂのそれぞれの機能は、処理回路により実現される。すなわち、移動経路推定装置４Ｄは、図１８を用いて後述するステップＳＴ１ｅからステップＳＴ５ｅまでの処理を実行するための処理回路を備える。
　処理回路は、専用のハードウェアの処理回路１０３であってもよいが、メモリ１０５に記憶されたプログラムを実行するプロセッサ１０４であってもよい。
　なお、データ保存部４６は、移動経路推定装置４Ｄとは別に設けられた記憶装置に構築してもよい。この場合、移動予測部４２Ｃは、移動経路推定装置４Ｄとは別に設けられた記憶装置と通信して、当該記憶装置からノード間経路データを取得する。

　次に動作について説明する。
　図１８は、実施の形態５に係る移動経路推定方法を示すフローチャートである。
　消失予測部４０Ｃが、入退検出センサ２Ｃ－１～２Ｃ－Ｎから受信した入退検出データおよびノード番号、データ更新部４３Ｂから入力した前時刻における観測対象候補の存在確率およびノード間確率分布に基づいて、前時刻から現時刻までの間の観測対象候補の消失が反映された観測対象候補の存在確率の予測値およびノード間確率分布の予測値を算出する（ステップＳＴ１ｅ）。

　次に、出現予測部４１Ｃが、入退検出センサ２Ｃ－１～２Ｃ－Ｎから受信した入退検出データおよびノード番号に基づいて、前時刻から現時刻までの間に出現した観測対象候補の存在確率の予測値およびノード間確率分布の予測値を算出する（ステップＳＴ２ｅ）。

　移動予測部４２Ｃが、消失予測部４０Ｃおよび出現予測部４１Ｃから入力した観測対象候補の存在確率の予測値およびノード間確率分布の予測値、およびデータ保存部４６から取得したノード間経路データに基づいて、前時刻から現時刻までの間の観測対象候補の移動が反映された観測対象候補の存在確率の予測値およびノード間確率分布の予測値を算出する（ステップＳＴ３ｅ）。

　データ更新部４３Ｂが、存在検出センサ３Ｂ－１～３Ｂ－Ｍから受信した存在検出データおよびノード番号、および、移動予測部４２Ｃから入力した観測対象候補の存在確率の予測値およびノード間確率分布の予測値に基づいて、現時刻における観測対象候補の存在確率およびノード間確率分布を算出する（ステップＳＴ４ｅ）。

　推定部４４Ｂが、データ更新部４３Ｂから入力した観測対象候補の存在確率およびノード間確率分布に基づいて、現時刻に観測対象が観測領域に存在するか否かを推定し、現時刻に存在すると推定した観測対象の移動経路を推定する（ステップＳＴ５ｅ）。例えば、推定部４４Ｂは、観測対象の移動経路を示す運動諸元推定値を算出し、算出した運動諸元推定値を移動経路保存部５Ｂに保存する。

　次に、観測対象の移動経路の推定処理について詳細に説明する。
　以降では、観測対象の運動諸元として、観測対象候補が存在する位置、すなわち、観測対象候補が存在するノードを推定する場合について説明する。なお、ノード間を移動する観測対象の運動諸元は、観測対象候補が存在する位置（ノード）に限定されるものではなく、推定対象の運動諸元には、観測対象候補の移動速度または加速度を含めてもよい。

　移動経路推定装置４Ｄは、観測対象候補の運動が時刻フレームごとに独立な確率過程によって選択される運動モデルを採用している。すなわち、観測対象候補が特定のノードへ移動すること、または、観測対象候補が特定のノードに留まることは、時刻フレームごとに独立な確率過程により選択される。ただし、これは例示であり、前時刻における運動が現時刻における運動に影響する場合、例えば速度または加速度の運動諸元が運動に影響を及ぼす場合には、等速モデルまたは等加速度モデルを用いてもよい。

　以降では、時刻フレーム間の経過時間が、観測対象候補がノード間の移動に要する時間に比べて十分に短いことを前提として、観測対象候補は、時刻フレーム間に複数のノードを経由したノード間の移動ができないものとする。
　ただし、時刻フレーム間に観測対象候補が経由し得るノードの個数の最大値が既知であれば、ノードを経由した移動に関する予測値を算出してもよい。観測対象候補のノード間の移動に要する時間の長さに応じて、観測対象候補がノード間を移動している状態を設定するモデルを採用してもよい。

　また、観測対象候補は、観測領域内の全てのノードに進入することが可能であり、退出も可能である。入退検出センサは、観測領域に設定された全てのノードに設置されているものとする。なお、観測対象候補が進入または退出できないノードを観測領域に設定した場合、このノードには、入退検出センサ２Ｃ－ｎを設置しなくてもよく、このノードにおける観測対象候補の出現および消失を予測しなくてもよい。

　入退検出センサ２Ｃ－１～２Ｃ－Ｎは、入退検出データの誤検出があるものと仮定し、存在検出センサ３Ｂ－１～３Ｂ－Ｍについても、存在検出データの誤検出または存在検出データの検出失敗があるものと仮定する。

　さらに、ノード間確率分布を、重み付けされた有限個数の“粒子”の集合で近似する。ノード間確率分布の値が大きい領域には、重み係数が大きい粒子が密に分布し、ノード間確率分布の値が小さい領域には、重み係数が小さな粒子が疎に分布する。なお、近似する方法は、前述した方法に限定されるものではなく、例えば、離散値化したガウス分布のそれぞれの線形結合によってノード間確率分布を近似してもよい。

　前述したノード番号とは別に、観測領域に設定された複数のノードのそれぞれには識別番号が付与される。以降では、ノードの総数をＮとする。Ｎは、入退検出センサの個数と同数である。観測対象候補が存在するノードの識別番号を、観測対象候補の運動諸元であるスカラー値ｘで表す。
　また、ｎ番目（ｎ＝１，２，・・・，Ｎ）の入退検出センサ２Ｃ－ｎが存在するノードのノード番号をｚ_ｋ，ｎ ^{（ＩｎＯｕｔ）}とする。ｍ番目（ｍ＝１，２，・・・，Ｍ）の存在検出センサ３Ｂ－ｍが存在するノードのノード番号をｚ_ｋ，ｍ ^{（Ｅｘｉｓｔ）}とする。

　実施の形態１と同様に、時刻ｋにおける観測対象候補の識別番号の集合をＬ_ｋで表す。時刻ｋにおける識別番号ｉの観測対象候補の存在確率をｒ_ｋ｜ｋ ^（ｉ）とし、時刻ｋ－１における識別番号ｉの観測対象候補の存在確率をｒ_{ｋ－１｜ｋ－１} ^（ｉ）とし、時刻ｋ－１における識別番号ｉの観測対象候補の存在確率を用いて算出された、時刻ｋにおける識別番号ｉの観測対象候補の存在確率の予測値をｒ_{ｋ｜ｋ－１} ^（ｉ）とする。

　時刻ｋにおける、識別番号ｉの観測対象候補のノード間確率分布をｆ_ｋ｜ｋ ^（ｉ）（ｘ）と表す。ｆ_ｋ｜ｋ ^（ｉ）（ｘ）は、時刻ｋにおいて識別番号ｉの観測対象候補が観測領域に存在する場合に、この観測対象候補がノード番号ｘに存在する確率を示している。離散化された確率分布の定義より、任意のｉについて、下記式（５０）が成立する。

　重み付けされた粒子の集合でノード間確率分布ｆ_ｋ｜ｋ ^（ｉ）（ｘ）を近似することで、識別番号ｉの観測対象候補のノード間確率分布ｆ_ｋ｜ｋ ^（ｉ）（ｘ）は上記式（４）で表すことができる。また、ｆ_{ｋ｜ｋ－１} ^（ｉ）（ｘ）は、時刻ｋ－１における識別番号ｉの観測対象候補のノード間確率分布を用いて算出された、時刻ｋおける識別番号ｉの観測対象候補のノード間確率分布の予測値である。重み付けされた粒子の集合でｆ_{ｋ｜ｋ－１} ^（ｉ）（ｘ）を近似することで、識別番号ｉの観測対象候補のノード間確率分布の予測値ｆ_{ｋ｜ｋ－１} ^（ｉ）（ｘ）は上記式（５）で表すことができる。

　前述したノード間確率分布の近似を前提として、以降では、時刻ｋにおける識別番号ｉの観測対象候補のノード間確率分布ｆ_ｋ｜ｋ ^（ｉ）（ｘ）を、Ｊ個の重み係数ｗ_ｋ｜ｋ ^{（ｉ、ｊ）}と粒子位置ｘ_ｋ｜ｋ ^{（ｉ、ｊ）}から算出するものとして説明する。これは、観測対象候補のノード間確率分布の予測値ｆ_{ｋ｜ｋ－１} ^（ｉ）（ｘ）においても同様である。

　データ保存部４６に保存されているノード間経路データは、Ｎ行Ｎ列の行列Ψとして、下記式（５１）で表される。ｎ１行ｎ２列（ｎ１とｎ２はともに１以上Ｎ以下の整数）の要素ψ_{ｎ１，ｎ２}は、識別番号ｎ１のノードから識別番号ｎ２のノードへのノード間経路が存在するか否かを示す要素であり、識別番号ｎ１のノードから、識別番号ｎ２のノードへのノード間経路が存在すれば“１”、存在しなければ“０”となる。Ψの対角成分は全て１とする。

　図１９は、移動経路推定装置４Ｄの動作を示すフローチャートである。
　消失予測部４０Ｃおよび出現予測部４１Ｃは、入退検出センサ２Ｃ－１～２Ｃ－Ｎから受信した時刻ｋにおけるＮ個の入退検出データのうち、観測領域に対する観測対象候補の進入または退出が検出されたデータを、任意の順番で選択して入力する（ステップＳＴ１ｆ）。これにより入退検出データを選択して処理するループが開始される。ステップＳＴ２ｆからステップＳＴ４ｆでは、ｎ番目（ｎ＝１，２，・・・，Ｎのいずれか）の入退検出センサ２Ｃ－ｎによって検出された入退検出データが選択されたものと仮定する。

　消失予測部４０Ｃおよび出現予測部４１Ｃは、ｎ番目の入退検出センサ２Ｃ－ｎによって検出された入退検出データに基づいて、観測対象候補の移動方向が観測領域への進入方向であるか否かを判定する（ステップＳＴ２ｆ）。
　消失予測部４０Ｃは、観測対象候補の移動方向が観測領域からの退出方向であると判定した場合（ステップＳＴ２ｆ；ＮＯ）、ステップＳＴ３ｆへ移行する。
　出現予測部４１Ｃは、観測対象候補の移動方向が観測領域への進入方向であると判定した場合（ステップＳＴ２ｆ；ＹＥＳ）、ステップＳＴ４ｆへ移行する。

　ステップＳＴ３ｆにおいて、消失予測部４０Ｃは、入退検出センサ２－ｎによって検出された入退検出データ、入退検出センサ２－ｎが設置されたノードのノード番号、時刻ｋ－１の観測対象候補の存在確率およびノード間確率分布に基づいて、時刻ｋ－１から時刻ｋまでの間の観測対象候補の消失が反映された、観測対象候補の存在確率の予測値およびノード間確率分布の予測値を算出する。

　なお、消失予測部４０Ｃが、ステップＳＴ３ｆの処理を１回目に行う場合、時刻ｋ－１における観測対象候補の存在確率およびノード間確率分布を用いる。２回目以降にステップＳＴ３ｆの処理を行う場合、消失予測部４０Ｃは、直前のステップＳＴ３ｆで算出した存在確率の予測値およびノード間確率分布の予測値を用いる。

　消失予測部４０Ｃは、下記式（５２）、上記式（７）、上記式（８）および下記式（５３）を用いて、時刻ｋ－１から時刻ｋまでの間の観測対象候補の消失が反映された識別番号ｉの観測対象候補の存在確率の予測値ｒ_{ｋ｜ｋ－１，Ｓ} ^（ｉ）を算出し、時刻ｋ－１から時刻ｋまでの間の観測対象候補の消失が反映された識別番号ｉの観測対象候補のノード間確率分布の予測値を表す粒子である｛ｗ_{ｋ｜ｋ－１，Ｓ} ^{（ｉ，ｊ）}，ｘ_{ｋ｜ｋ－１，Ｓ} ^{（ｉ，ｊ）}｝_{ｊ＝１，・・・，Ｊ}を算出する。ただし、識別番号ｉは前時刻の観測対象候補と同一である。下記式（５２）におけるβ_{ｋ，ＦＬＳ} ^{（ＩｎＯｕｔ）}は、入退検出センサが誤検出する確率である。

　上記式（７）、上記式（８）、上記式（５２）および上記式（５３）を用いて観測対象候補の存在確率の予測値およびノード間確率分布の予測値を算出する処理は、識別番号ｉの観測対象候補のノード間確率分布を示す粒子のうち、合計の重み係数が最も多いノードのノード番号μ_{ｋ－１｜ｋ－１} ^（ｉ）が、観測対象候補の退出を検出したｎ番目の入退検出センサ２Ｃ－ｎが設置されたノードのノード番号ｚ_ｎ ^{（ＩｎＯｕｔ）}と同一であれば、識別番号ｉの観測対象候補の存在確率をβ_{ｋ，ＦＬＳ} ^{（ＩｎＯｕｔ）}倍に減少させる処理に相当する。この処理により、入退検出センサ２Ｃ－ｎによって検出された観測対象候補の退出を、観測対象候補の予測値に反映させることができる。

　識別番号ｉの観測対象候補の存在確率をβ_{ｋ，ＦＬＳ} ^{（ＩｎＯｕｔ）}倍することで、識別番号ｉの観測対象候補の存在確率に対して観測対象候補の退出（消失）を反映させたが、これに限定されるものではない。例えば、β_{ｋ，ＦＬＳ} ^{（ＩｎＯｕｔ）}＝０として、上記式（９）に示すｄ^{（ｉ，ｎ）}が閾値よりも小さい場合に、消失予測部４０Ｃが、識別番号ｉの観測値候補の存在確率およびノード間確率分布（粒子）を削除してもよい。

　ステップＳＴ４ｆにおいて、出現予測部４１Ｃは、入退検出センサ２Ｃ－ｎによって検出された入退検出データ、および入退検出センサ２－ｎが設置されたノードのノード番号に基づいて、時刻ｋ－１から時刻ｋまでの間に出現した観測対象候補の存在確率の予測値およびノード間確率分布の予測値を算出する。

　出現予測部４１Ｃは、上記式（１２）、上記式（１３）および下記式（５４）に用いて、時刻ｋ－１から時刻ｋまでの間に出現した観測対象候補の存在確率の予測値ｒ_{ｋ｜ｋ－１，Ｂ} ^（ｉｎ）を算出し、時刻ｋ－１から時刻ｋまでに出現した観測対象候補のノード間確率分布の予測値を表す粒子である｛ｗ_{ｋ｜ｋ－１，Ｂ} ^{（ｉｎ，ｊ）}，ｘ_{ｋ｜ｋ－１，Ｂ} ^{（ｉｎ，ｊ）}｝_{ｊ＝１，・・・，Ｊ}を算出する。識別番号ｉ_ｎは、ｎ番目の入退検出センサ２Ｃ－ｎの周辺に出現したと予測された観測対象候補の識別番号であり、既存の観測対象候補の識別番号と異なる値である。

　上記式（１２）、上記式（１３）および上記式（５４）に用いて、存在確率の予測値およびノード間確率分布の予測値を算出する処理は、ｎ番目の入退検出センサ２－ｎの設置位置を初期位置として新規の観測対象候補の粒子を生成する処理に相当する。この処理により、時刻ｋに入退検出センサ２Ｃ－ｎによって観測領域への進入が検出された観測対象候補の予測値が得られる。

　ステップＳＴ５ｆにおいて、消失予測部４０Ｃおよび出現予測部４１Ｃは、入退検出センサ２Ｃ－１～２Ｃ－Ｎから受信した時刻ｋにおけるＮ個の入退検出データのうち、観測領域に対する観測対象候補の進入または退出が検出された全てのデータを選択したか否かを判定する。このとき、未選択のデータがあると判定した場合（ステップＳＴ５ｆ；ＮＯ）、消失予測部４０Ｃおよび出現予測部４１Ｃは、未選択の入退検出データを選択して、ステップＳＴ２ｆの処理に戻る。全てのデータを選択していれば（ステップＳＴ５ｆ；ＹＥＳ）、ステップＳＴ６ｆの処理に移行する。

　ステップＳＴ６ｆにおいて、移動予測部４２Ｃは、時刻ｋ－１から時刻ｋまでの間の観測対象候補の消失が反映された存在確率の予測値およびノード間確率分布の予測値、時刻ｋ－１から時刻ｋまでの間に出現した観測対象候補の存在確率の予測値およびノード間確率分布の予測値、およびノード間経路データに基づいて、時刻ｋ－１から時刻ｋまでの間の観測対象候補の移動が反映された、観測対象候補の存在確率の予測値およびノード間確率分布の予測値を算出する。

　移動予測部４２Ｃは、時刻ｋ－１から時刻ｋまで間の観測対象候補の消失が反映された観測対象候補の存在確率の予測値ｒ_{ｋ｜ｋ－１，Ｓ} ^（ｉｓ）、時刻ｋ－１から時刻ｋまで間の観測対象候補の消失が反映されたノード間確率分布の予測値｛ｗ_{ｋ｜ｋ－１，Ｓ} ^{（ｉｓ，ｊ）}，ｘ_{ｋ｜ｋ－１，Ｓ} ^{（ｉｓ，ｊ）}｝_{ｊ＝１，・・・，Ｊ}、時刻ｋ－１から時刻ｋまでの間に出現した観測対象候補の存在確率の予測値ｒ_{ｋ｜ｋ－１，Ｂ} ^（ｉＢ）、時刻ｋ－１から時刻ｋまでの間に出現した観測対象候補のノード間確率分布の予測値｛ｗ_{ｋ｜ｋ－１，Ｂ} ^{（ｉＢ，ｊ）}，ｘ_{ｋ｜ｋ－１，Ｂ} ^{（ｉＢ，ｊ）}｝_{ｊ＝１，・・・，Ｊ}、およびノード間経路データ（上記式（５１）に示した行列Ψ）に基づいて、時刻ｋ－１から時刻ｋまで間の観測対象候補の移動が反映された観測対象候補の存在確率の予測値ｒ_{ｋ｜ｋ－１} ^（ｉ）およびノード間確率分布の予測値｛ｗ_{ｋ｜ｋ－１} ^{（ｉ，ｊ）}，ｘ_{ｋ｜ｋ－１} ^{（ｉ，ｊ）}｝_{ｊ＝１，・・・，Ｊ}を算出する。ｉ_Ｓは、ステップＳＴ３ｆで予測値が算出された観測対象候補の識別番号である。ｉ_Ｂは、ステップＳＴ４ｆで予測値が算出された観測対象候補の識別番号である。識別番号ｉは、ｉ_Ｓとｉ_Ｂの両方を含む識別番号である。ｉ_Ｓの集合をＬ_{ｋ｜ｋ－１，Ｓ}とし、ｉ_Ｂの集合をＬ_{ｋ｜ｋ－１，Ｂ}とすることで、ｉは上記式（１７）で定義される。

　移動予測部４２Ｃは、上記式（１８）、上記式（１９）および下記式（５５）を用いて、観測対象候補の移動が反映された観測対象候補の存在確率の予測値ｒ_{ｋ｜ｋ－１} ^（ｉ）および観測対象候補のノード間確率分布の予測値｛ｗ_{ｋ｜ｋ－１} ^{（ｉ，ｊ）}，ｘ_{ｋ｜ｋ－１} ^{（ｉ，ｊ）}｝_{ｊ＝１，・・・，Ｊ}を算出する。ただし、｛ｎ｜Ψ_ｎ’，ｎ＝１｝は、ψ_ｎ’，ｎ＝１を満たす値ｎの集合であり、Ｕ（ｘ；｛ａ_１，・・・，ａ_Ａ｝）は、離散値の集合｛ａ_１，・・・，ａ_Ａ｝のいずれかの要素を同じ確率で選択したときに離散値ｘが従う一様分布の確率密度関数である。すなわち、下記式（５５）は、前時刻の粒子位置から移動可能なノードのうち、いずれか１つのノードを一様な乱数で選択することを示している。

　なお、ノード間確率分布を構成する粒子においてノード間経路が選択される確率が一様であるものとしたが、これに限定されるものではない。
　例えば、ノード間経路ごとの起こりやすさに応じたパラメータによって、観測対象候補同士で移動先のノードが選択される確率を互いに異なる確率にしてもよい。

　ステップＳＴ７ｆにおいて、データ更新部４３Ｂは、存在検出センサ３Ｂ－１～３Ｂ－Ｍから受信した存在検出データおよびノード番号、時刻ｋ－１から時刻ｋまでの間の観測対象候補の移動が反映された観測対象候補の存在確率の予測値およびノード間確率分布の予測値に基づいて、時刻ｋにおける観測対象候補の存在確率およびノード間確率分布を算出（更新）する。

　図２０は、図１９のステップＳＴ７ｆの詳細な処理を示すフローチャートである。
　ステップＳＴ７ｆ－１において、データ更新部４３Ｂは、識別番号ｉの観測対象候補の存在確率の予測値ｒ_{ｋ｜ｋ－１} ^（ｉ）のそれぞれの値に基づいて、識別番号の部分集合ごとの信頼度ｂ_{ｋ｜ｋ－１}（Ｉ^（λ））を算出する。λ＝１，２，・・・，Λである。
　データ更新部４３Ｂは、上記式（２２）を用いて、部分集合ごとの信頼度ｂ_{ｋ｜ｋ－１}（Ｉ^（λ））を算出する。識別番号の部分集合Ｉ^（λ）とは、０個以上の識別番号ｉから構成される集合である。例えば、ｉ＝１，２である場合、識別番号の部分集合Ｉ^（λ）は、λ＝１，２，３，４として、Ｉ^（１）＝φ、Ｉ^（２）＝｛１｝、Ｉ^（３）＝｛２｝、Ｉ^（４）＝｛１，２｝の４通りとなる。ここで、φは空集合を表す。

　ステップＳＴ７ｆ－２において、データ更新部４３Ｂは、存在検出センサ３Ｂ－１～３Ｂ－Ｍのうち、時刻ｋにおいて観測対象候補を検出した存在検出センサが設置されたノードのノード番号｛ｚ_ｋ，ｍ ^{（Ｅｘｉｓｔ）}｝_ｍ∈Ｍ＊、部分集合信頼度ｂ_{ｋ｜ｋ－１}（Ｉ^（λ））および識別番号ｉのノード間確率分布の予測値｛ｗ_{ｋ｜ｋ－１} ^{（ｉ，ｊ）}，ｘ_{ｋ｜ｋ－１} ^{（ｉ，ｊ）}｝_{ｊ＝１，・・・，Ｊ}に基づいて、写像ごとの信頼度ｂ’_ｋ｜ｋ（θ_ｋ ^{（λ，ｈ）}）と写像ごとの識別番号ｉのノード間確率分布｛ｗ_ｋ｜ｋ ^{（ｉ，ｊ）}（θ_ｋ ^{（λ，ｈ）}），ｘ_ｋ｜ｋ ^{（ｉ，ｊ）}（θ_ｋ ^{（λ，ｈ）}）｝_{ｊ＝１，・・・，Ｊ}を算出する。写像θ_ｋ ^{（λ，ｈ）}は、上記実施の形態１で示したように、部分集合Ｉ^（λ）の要素を、時刻ｋにおいて観測対象候補を検出した存在検出センサの存在検出センサ番号の少なくとも１つに対応付ける写像であると定義する。これを数式で表現すると、上記式（２３）のようになる。

　データ更新部４３Ｂは、上記式（２４）、上記式（２４ａ）および上記式（２５）を用いて、写像ごとの信頼度ｂ’_ｋ｜ｋ（θ_ｋ ^{（λ，ｈ）}）および写像ごとの識別番号ｉのノード間確率分布｛ｗ_ｋ｜ｋ ^{（ｉ，ｊ）}（θ_ｋ ^{（λ，ｈ）}），ｘ_ｋ｜ｋ ^{（ｉ，ｊ）}（θ_ｋ ^{（λ，ｈ）}）｝_{ｊ＝１，・・・，Ｊ}を算出する。
　ただし、上記式（２４）、上記式（２４ａ）および上記式（２５）における変数は、上記式（２６）および下記式（５６）に従って算出される。
　なお、下記式（５６）におけるｐ_ｋ，ＤＴは、時刻ｋにおいて存在検出センサが検出に失敗せずに観測対象候補を検出する確率である。β_{ｋ，ＦＬＳ} ^{（Ｅｘｉｓｔ）}は、時刻ｋにおいて存在検出センサが誤検出する確率である。

　ステップＳＴ７ｆ－３において、データ更新部４３Ｂは、写像ごとの信頼度ｂ’_ｋ｜ｋ（θ_ｋ ^{（λ，ｈ）}）を規格化し、写像ごとの識別番号ｉのノード間確率分布における重み係数｛ｗ’_ｋ｜ｋ ^{（ｉ，ｊ）}（θ_ｋ ^{（λ，ｈ）}）｝_{ｊ＝１，・・・，Ｊ}を規格化する。データ更新部４３Ｂは、上記式（３０）および上記式（３１）を用いて、信頼度ｂ’_ｋ｜ｋ（θ_ｋ ^{（λ，ｈ）}）を規格化し、写像ごとの識別番号ｉのノード間確率分布における重み係数｛ｗ’_ｋ｜ｋ ^{（ｉ，ｊ）}（θ_ｋ ^{（λ，ｈ）}）｝_{ｊ＝１，・・・，Ｊ}を規格化して、信頼度ｂ_ｋ｜ｋ（θ_ｋ ^{（λ，ｈ）}）および重み係数｛ｗ_ｋ｜ｋ ^{（ｉ，ｊ）}（θ_ｋ ^{（λ，ｈ）}）｝_{ｊ＝１，・・・，Ｊ}を得る。

　ステップＳＴ７ｆ－４において、データ更新部４３Ｂは、上記式（３２）、上記式（３３）および上記式（３４）を用いて、規格化された写像ごとの信頼度ｂ_ｋ｜ｋ（θ_ｋ ^{（λ，ｈ）}）、および、重み係数が規格化された写像ごとの識別番号ｉのノード間確率分布｛ｗ_ｋ｜ｋ ^{（ｉ，ｊ）}（θ_ｋ ^{（λ，ｈ）}），ｘ_ｋ｜ｋ ^{（ｉ，ｊ）}（θ_ｋ ^{（λ，ｈ）}）｝_{ｊ＝１，・・・，Ｊ}に基づき、時刻ｋにおける識別番号ｉの観測対象候補の存在確率ｒ_ｋ｜ｋ ^（ｉ）と、識別番号ｉの観測対象候補のノード間確率分布｛ｗ_ｋ｜ｋ ^{（ｉ，ｊ）}，ｘ_ｋ｜ｋ ^{（ｉ，ｊ）}｝_{ｊ＝１，・・・，Ｊ}を算出する。

　図１９の説明に戻る。
　ステップＳＴ８ｆにおいて、推定部４４Ｂは、時刻ｋ（現時刻）における観測対象候補のうちのいずれか１つを選択して入力する。これにより、観測対象候補を選択して処理するループが開始される。以降では、識別番号ｉの観測対象候補が選択された場合を例に挙げて説明する。

　次に、推定部４４Ｂは、上記式（３５）に従い、ステップＳＴ８ｆで選択した識別番号ｉの観測対象候補の存在確率ｒ_ｋ｜ｋ ^（ｉ）が閾値ｒ_Ｔｈよりも大きいか否かを判定する（ステップＳＴ９ｆ）。識別番号ｉの観測対象候補の存在確率ｒ_ｋ｜ｋ ^（ｉ）が閾値ｒ_Ｔｈよりも大きいと判定された場合（ステップＳＴ９ｆ；ＹＥＳ）、ステップＳＴ１０ｆの処理に移行する。一方、識別番号ｉの観測対象候補の存在確率ｒ_ｋ｜ｋ ^（ｉ）が閾値ｒ_Ｔｈ以下であれば（ステップＳＴ９ｆ；ＮＯ）、ステップＳＴ１１ｆの処理に移行する。

　ステップＳＴ１０ｆにおいて、推定部４４Ｂは、時刻ｋにおける識別番号ｉの観測対象候補のノード間確率分布に基づいて、当該観測対象候補に対応する観測対象の運動諸元を推定する。例えば、推定部４４Ｂは、識別番号ｉの観測対象候補に対応する観測対象の運動諸元推定値、すなわち、観測対象候補が存在すると推定されるノード番号μ_ｋ ^（ｉ）を、下記式（５７）に従って算出する。

　移動経路保存部５Ｂには、異なる時刻フレームに保存された識別番号が同一の運動諸元推定値同士が時系列順に保存される。これにより、運動諸元推定値は、観測対象の推定された移動経路として移動経路保存部５Ｂに保存される。

　以上のように、実施の形態５に係る移動経路推定装置４Ｄは、観測対象候補の消失および出現が反映された観測対象候補の存在確率の予測値およびノード間確率分布の予測値に基づいて、観測時刻における観測対象候補の存在確率およびノード間確率分布を更新する。移動経路推定装置４Ｄは、更新した観測対象候補の存在確率およびノード間確率分布に基づいて、観測対象の移動経路を示す運動諸元推定値を算出する。
　このように、移動経路推定装置４Ｄは、観測対象候補の消失および出現を予測して、個数が未知の観測対象の移動経路を推定することができる。

　また、移動予測部４２Ｃは、ノード間経路に制限された観測対象候補の移動を予測するため、観測対象候補のノード間確率分布が実現可能な範囲で絞り込まれる。これにより、観測対象候補がノード間経路のみを移動するという条件において、観測対象候補の移動可能なノード間経路を制限しない場合に比べて、観測対象候補の移動の予測値と実際の観測対象の移動との乖離が小さくなる。このため、より誤差の小さい移動経路推定値を得ることができる。従って、個数が未知の観測対象候補のノード間の移動経路を推定することができる。

　なお、本発明は上記実施の形態に限定されるものではなく、本発明の範囲内において、実施の形態のそれぞれの自由な組み合わせまたは実施の形態のそれぞれの任意の構成要素の変形もしくは実施の形態のそれぞれにおいて任意の構成要素の省略が可能である。

　この発明に係る移動経路推定装置は、個数が未知の観測対象の移動経路を推定することができるので、例えば、監視領域の人の入退を監視する監視システムに利用可能である。

　１，１Ａ～１Ｄ　移動経路推定システム、２－１～２－Ｎ，２－ｎ，２Ａ－１～２Ａ－Ｎ，２Ａ－ｎ，２Ｂ－１～２Ｂ－Ｎ，２Ｂ－ｎ，２Ｃ－１～２Ｃ－Ｎ，２Ｃ－ｎ　入退検出センサ、３－１～３－Ｍ，３－ｍ，３Ａ－１～３Ａ－Ｍ，３Ａ－ｍ，３Ｂ－１～３Ｂ－Ｍ，３Ｂ－ｍ　存在検出センサ、４，４Ａ～４Ｄ　移動経路推定装置、５，５Ａ，５Ｂ　移動経路保存部、１０，１０Ａ，１０Ｂ　観測領域、１１－ｎ，１２－ｍ　検出範囲、１３ａ，１３ｂ　移動経路、１４　進入不可領域、４０，４０Ａ～４０Ｃ　消失予測部、４１，４１Ａ～４１Ｃ　出現予測部、４２，４２Ａ～４２Ｃ　移動予測部、４３，４３Ａ，４３Ｂ　データ更新部、４４，４４Ａ，４４Ｂ　推定部、４５，４６　データ保存部、１００，１０１　センサ群、１０２　記憶装置、１０３　処理回路、１０４　プロセッサ、１０５　メモリ、２００～２０５　ノード、３００　ノード間経路。

Claims

　観測領域内の検出範囲を通過した観測対象候補が通過方向ごとに検出された入退検出データ、前記入退検出データが得られた入退検出センサ位置、観測時刻よりも一定時間だけ過去の時刻において観測対象候補が前記観測領域に存在する可能性を示す存在確率および前記過去の時刻における観測対象候補の運動諸元確率分布に基づいて、前記過去の時刻から前記観測時刻までの間の観測対象候補の消失が反映された、観測対象候補の存在確率の予測値および運動諸元確率分布の予測値を算出する消失予測部と、
　前記入退検出データおよび前記入退検出センサ位置に基づいて、前記過去の時刻から前記観測時刻までの間に出現した観測対象候補の存在確率の予測値および運動諸元確率分布の予測値を算出する出現予測部と、
　前記消失予測部および前記出現予測部によってそれぞれ算出された観測対象候補の存在確率の予測値および運動諸元確率分布の予測値に基づいて、前記過去の時刻から前記観測時刻までの間の観測対象候補の移動が反映された、観測対象候補の存在確率の予測値および運動諸元確率分布の予測値を算出する移動予測部と、
　検出範囲に存在する観測対象候補が検出された存在検出データ、前記存在検出データが得られた存在検出センサ位置、および前記移動予測部によって算出された観測対象候補の存在確率の予測値および運動諸元確率分布の予測値に基づいて、前記観測時刻における観測対象候補の存在確率および運動諸元確率分布を算出するデータ更新部と、
　前記データ更新部によって算出された観測対象候補の存在確率および運動諸元確率分布に基づいて、前記観測時刻に観測対象が前記観測領域に存在するか否かを推定し、前記観測時刻に存在すると推定した観測対象の運動諸元推定値を算出する推定部とを備えたこと
　を特徴とする移動経路推定装置。
　前記消失予測部は、前記入退検出データ、前記入退検出センサ位置、前記過去の時刻における観測対象候補の存在確率および運動諸元確率分布、および検出範囲からの観測対象候補の退出が誤検出される確率を用いて、観測対象候補の消失が反映された観測対象候補の存在確率の予測値および運動諸元確率分布の予測値を算出し、
　前記出現予測部は、前記入退検出データ、前記入退検出センサ位置および検出範囲への観測対象候補の進入が誤検出される確率を用いて、観測対象候補の出現が反映された観測対象候補の存在確率の予測値および運動諸元確率分布の予測値を算出すること
　を特徴とする請求項１記載の移動経路推定装置。
　前記データ更新部は、前記存在検出データ、前記存在検出センサ位置、前記移動予測部によって算出された観測対象候補の存在確率の予測値および運動諸元確率分布の予測値、検出範囲に存在しない観測対象候補が誤検出される確率および検出範囲に存在する観測対象候補が検出される確率を用いて、前記観測時刻における観測対象候補の存在確率および運動諸元確率分布を算出すること
　を特徴とする請求項１または請求項２記載の移動経路推定装置。
　前記データ更新部は、複数の観測対象候補と複数の前記存在検出データとが対応付けられる可能性が反映された観測対象候補の運動諸元確率分布を算出すること
　を特徴とする請求項１記載の移動経路推定装置。
　前記入退検出データは、検出範囲を通過した観測対象候補の数を含み、
　前記消失予測部は、検出範囲を通過した観測対象候補の数を含む前記入退検出データ、前記入退検出センサ位置、および前記過去の時刻における観測対象候補の存在確率および運動諸元確率分布に基づいて、前記過去の時刻から前記観測時刻までの間の観測対象候補の消失が反映された観測対象候補の存在確率の予測値および運動諸元確率分布の予測値を算出し、
　前記出現予測部は、検出範囲を通過した観測対象候補の数を含む前記入退検出データおよび前記入退検出センサ位置に基づいて、前記過去の時刻から前記観測時刻までの間に出現した観測対象候補の存在確率の予測値および運動諸元確率分布の予測値を算出すること
　を特徴とする請求項１記載の移動経路推定装置。
　前記移動予測部は、観測対象候補が存在し得ない領域を示す領域データ、および、前記消失予測部および前記出現予測部によってそれぞれ算出された観測対象候補の存在確率の予測値および運動諸元確率分布の予測値に基づいて、前記過去の時刻から前記観測時刻までの間の観測対象候補の移動が反映された観測対象候補の存在確率の予測値および運動諸元確率分布の予測値を算出すること
　を特徴とする請求項１記載の移動経路推定装置。
　前記入退検出データは、観測対象候補の属性を含み、
　前記存在検出データは、観測対象候補の属性を含み、
　前記消失予測部は、観測対象候補の属性を含む前記入退検出データ、前記入退検出センサ位置、および、前記過去の時刻における観測対象候補の存在確率および運動諸元確率分布に基づいて、前記過去の時刻から前記観測時刻までの間の観測対象候補の消失が反映された、観測対象候補の存在確率の予測値、運動諸元確率分布の予測値および予測値に対応する観測対象候補の属性を示す属性値を算出し、
　前記出現予測部は、観測対象候補の属性を含む前記入退検出データおよび前記入退検出センサ位置に基づいて、前記過去の時刻から前記観測時刻までの間に出現した観測対象候補の存在確率の予測値、運動諸元確率分布の予測値および予測値に対応する観測対象候補の属性を示す属性値を算出し、
　前記移動予測部は、前記消失予測部および前記出現予測部によってそれぞれ算出された観測対象候補の存在確率の予測値、運動諸元確率分布の予測値、および予測値に対応する観測対象候補の属性を示す属性値に基づいて、前記過去の時刻から前記観測時刻までの間の観測対象候補の移動が反映された、観測対象候補の存在確率の予測値、運動諸元確率分布の予測値および予測値に対応する観測対象候補の属性を示す属性値を算出し、
　前記データ更新部は、観測対象候補の属性を含む前記存在検出データ、前記存在検出センサ位置、前記移動予測部によって算出された観測対象候補の存在確率の予測値、運動諸元確率分布の予測値および属性値に基づいて、前記観測時刻における観測対象候補の存在確率、運動諸元確率分布および属性値を算出し、
　前記推定部は、前記データ更新部によって算出された観測対象候補の存在確率、運動諸元確率分布および属性値に基づいて、前記観測時刻に観測対象が前記観測領域に存在するか否かを推定し、前記観測時刻に存在すると推定した観測対象の運動諸元推定値および属性値を算出すること
　を特徴とする請求項１記載の移動経路推定装置。
　観測領域に設定されたノードから前記観測領域外へ移動する観測対象候補および前記観測領域外からノードへ移動する観測対象候補が検出された入退検出データ、前記入退検出データが得られたノードを示す入退検出ノード情報、観測時刻よりも一定時間だけ過去の時刻において観測対象候補が前記観測領域に存在する可能性を示す存在確率、および前記過去の時刻において観測対象候補が存在する位置の確率分布であるノード間確率分布に基づいて、前記過去の時刻から前記観測時刻までの間の観測対象候補の消失が反映された、観測対象候補の存在確率の予測値およびノード間確率分布の予測値を算出する消失予測部と、
　前記入退検出データおよび前記入退検出ノード情報に基づいて、前記過去の時刻から前記観測時刻までの間に出現した観測対象候補の存在確率の予測値およびノード間確率分布の予測値を算出する出現予測部と、
　前記消失予測部および前記出現予測部によってそれぞれ算出された観測対象候補の存在確率の予測値およびノード間確率分布の予測値、および観測対象候補が移動可能なノード間の経路を示すノード間経路データに基づいて、前記過去の時刻から前記観測時刻までの間の観測対象候補の移動が反映された、観測対象候補の存在確率の予測値およびノード間確率分布の予測値を算出する移動予測部と、
　検出範囲に存在する観測対象候補が検出された存在検出データ、前記存在検出データが得られたノードを示す存在検出ノード情報および前記移動予測部によって算出された観測対象候補の存在確率の予測値およびノード間確率分布の予測値に基づいて、前記観測時刻における観測対象候補の存在確率およびノード間確率分布を算出するデータ更新部と、
　前記データ更新部によって算出された観測対象候補の存在確率およびノード間確率分布に基づいて、前記観測時刻に観測対象が前記観測領域に存在するか否かを推定し、前記観測時刻に存在すると推定した観測対象の運動諸元推定値を算出する推定部とを備えたこと
　を特徴とする移動経路推定装置。
　前記入退検出データを検出する複数の入退検出センサと、
　前記存在検出データを検出する複数の存在検出センサと、
　請求項１から請求項８のうちのいずれか１項記載の移動経路推定装置とを備えたこと
　を特徴とする移動経路推定システム。
　消失予測部が、観測領域内の検出範囲を通過した観測対象候補が通過方向ごとに検出された入退検出データ、前記入退検出データが得られた入退検出センサ位置、観測時刻よりも一定時間だけ過去の時刻において観測対象候補が前記観測領域に存在する可能性を示す存在確率および前記過去の時刻における観測対象候補の運動諸元確率分布に基づいて、前記過去の時刻から前記観測時刻までの間の観測対象候補の消失が反映された、観測対象候補の存在確率の予測値および運動諸元確率分布の予測値を算出するステップと、
　出現予測部が、前記入退検出データおよび前記入退検出センサ位置に基づいて、前記過去の時刻から前記観測時刻までの間に出現した観測対象候補の存在確率の予測値および運動諸元確率分布の予測値を算出するステップと、
　移動予測部が、前記消失予測部および前記出現予測部によってそれぞれ算出された観測対象候補の存在確率の予測値および運動諸元確率分布の予測値に基づいて、前記過去の時刻から前記観測時刻までの間の観測対象候補の移動が反映された、観測対象候補の存在確率の予測値および運動諸元確率分布の予測値を算出するステップと、
　データ更新部が、検出範囲に存在する観測対象候補が検出された存在検出データ、前記存在検出データが得られた存在検出センサ位置、および前記移動予測部によって算出された観測対象候補の存在確率の予測値および運動諸元確率分布の予測値に基づいて、前記観測時刻における観測対象候補の存在確率および運動諸元確率分布を算出するステップと、
　推定部が、前記データ更新部によって算出された観測対象候補の存在確率および運動諸元確率分布に基づいて、前記観測時刻に観測対象が前記観測領域に存在するか否かを推定し、前記観測時刻に存在すると推定した観測対象の運動諸元推定値を算出するステップとを備えたこと
　を特徴とする移動経路推定方法。