(57)【特許請求の範囲】
【請求項１】
複数の粒子状対象物の時系列画像を撮像する撮像手段と、
撮像された前記時系列画像内の各画素に対応する各小領域内の前記粒子状対象物に関する特徴情報を算出する特徴情報算出手段と、
前記特徴情報の時間変化を検出する時間変化検出手段と、
前記特徴情報の前記時間変化に基づいて、前記各小領域の中で前記粒子状対象物を含む小領域を選択し、
選択された前記小領域の周辺部分に存在する境界部における前記時間変化に関する所定の演算を行うことにより、異常領域を特定する因子特定手段と、を有することを特徴とする粒子群運動解析システム。
【請求項２】
前記因子特定手段によってなされる前記所定の演算は、
前記小領域内において他の領域よりも大きな前記時間変化の原因となっている領域を抽出することによって前記異常領域を特定することを特徴とする請求項１に記載の粒子群運動解析システム。
【請求項３】
前記境界部における前記時間変化についての所定の演算は、
前記小領域内において前記境界部の内側に含まれる特定領域の前記時間変化を、前記境界部の前記時間変化へ加算すること、または、前記境界部の前記時間変化を前記小領域内における前記特定領域の前記時間変化へ加算することを特徴とする請求項１または２に記載の粒子群運動解析システム。
【請求項４】
前記特徴情報は、
前記小領域における前記粒子状対象物の数量と、
前記小領域における前記粒子状対象物の領域の占める割合と、
のうちの少なくとも一つであることを特徴とする請求項１から３のいずれか一項に記載の粒子群運動解析システム。
【請求項５】
前記粒子状対象物に応じて、前記所定の小領域の少なくとも大きさと、形状と、異方性とのいずれか一つを可変とすることを特徴とする請求項１から４のいずれか一項に記載の粒子群運動解析システム。
【請求項６】
前記粒子状対象物は、生物の細胞内における粒子であることを特徴とする請求項１から５のいずれか一項に記載の粒子群運動解析システム。
【請求項７】
前記粒子状対象物である前記生物の細胞内における前記粒子は、複数の異なる前記時系列画像において、少なくとも、消失することと、出現することと、大きさを変えることとのいずれか一つを行うことを特徴とする請求項６に記載の粒子群運動解析システム。
【請求項８】
前記粒子状対象物は、不特定多数の人間が集合した集団における人間であることを特徴とする請求項１から５のいずれか一項に記載の粒子群運動解析システム。
【請求項９】
複数の粒子状対象物の画像を時系列順序で撮像する撮像ステップと、
撮像された前記画像内の各画素に対応する各小領域内の粒子状対象物に関する特徴情報を算出する特徴情報算出ステップと、
前記特徴情報の時間変化を検出する時間変化検出ステップと、
前記特徴情報の前記時間変化に基づいて、前記各小領域の中で前記粒子状対象物を含む小領域を選択し、
選択された前記小領域の周辺部分に存在する境界部における前記時間変化に関する所定の演算を行うことにより、異常領域を特定する因子特定ステップと、を有することを特徴とする粒子群運動解析方法。
【請求項１０】
粒子状対象物の運動を解析する粒子群運動解析システム用のプログラムであり、コンピュータにより読み取り可能なプログラムであって、
複数の前記粒子状対象物の時系列画像を撮像する撮像手段と、
撮像された前記時系列画像内の各画素に対応する各小領域内の粒子状対象物に関する特徴情報を算出する特徴情報算出手段と、
前記特徴情報の時間変化を検出する時間変化検出手段と、
前記特徴情報の前記時間変化に基づいて、前記各小領域の中で前記粒子状対象物を含む小領域を選択し、
選択された前記小領域の周辺部分に存在する境界部における前記時間変化に関する所定の演算を行うことにより、異常領域を特定する因子特定手段としてコンピュータを機能させることを特徴とするプログラム。
【請求項１１】
複数の粒子状対象物の時系列画像を撮像する撮像手段と、
撮像された前記時系列画像内の各画素に対応する各小領域内の前記粒子状対象物に関する特徴情報を算出する特徴情報算出手段と、
前記特徴情報の時間変化を検出する時間変化検出手段と、
前記時間変化の原因となった前記粒子状対象物を含む領域または前記時間変化の原因となった前記粒子状対象物の近傍の領域である異常領域を特定する因子特定手段とを有し、
前記特徴情報は、
前記小領域における前記粒子状対象物の領域の占める割合であることを特徴とする粒子群運動解析システム。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、粒子群運動解析システム、特に細胞内の粒子の動きや群集内の人間の動き等を解析するシステム、粒子群運動解析方法及びプログラムに関するものである。
【背景技術】
【０００２】
従来、顕微鏡観察下における細胞の画像において、細胞をトラッキングする技術や細胞にマーキングする技術が提案されている（例えば、日本国特開２００３−２０７７１９号公報、日本国特許第２９３０６１８号公報参照）。また、顕微鏡を用いた生物細胞の蛍光観察により、細胞内のＤＮＡの観察や染色体異常の解析等が行われるようになっている。生物の細胞内には多くの粒子が存在している。粒子は、基本的にランダムな動きをする。そして、粒子の中には、他の粒子とは異なる特別な働きをするものがある。特別な働きをする粒子は、他の粒子とは異なる特殊な動きをする。このような特別な働きをする粒子を検出できることが望まれている。
しかしながら、日本国特開２００３−２０７７１９号公報に開示された構成では、観察対象である細胞が常に観察視野内に存在するように細胞をトラッキングするものである。また、日本国特許第２９３０６１８号公報に開示された構成では、分裂を繰りかえす細胞に次々とマーキングを行うものである。いずれの構成によっても、細胞内の複数の粒子のうちから、特殊な働きをする粒子、即ち特殊な動きをする粒子を検出することはできない。
本発明は、上記に鑑みてなされたものであって、複数の粒子状対象物の中から、他の粒子状対象物とは異なる動きをする粒子状対象物を検出できる粒子群運動解析システム、粒子群運動解析方法及びプログラムを提供することを目的とする。
【発明の開示】
【０００３】
上述した課題を解決し、目的を達成するために、本発明によれば、複数の粒子状対象物の時系列画像を撮像する撮像手段と、撮像された時系列画像内の各画素に対応する各小領域内の粒子状対象物に関する特徴情報を算出する特徴情報算出手段と、特徴情報の時間変化を検出する時間変化検出手段と、特徴情報の時間変化に基づいて、各小領域の中で粒子状対象物を含む小領域を選択し、選択された小領域の周辺部分に存在する境界部における時間変化に関する所定の演算を行うことにより、異常領域を特定する因子特定手段と、を有することを特徴とする粒子群運動解析システムを提供できる。
また、本発明の好ましい態様によれば、因子特定手段によってなされる所定の演算は、小領域内において他の領域よりも大きな時間変化の原因となっている領域を抽出することによって異常領域を特定することが望ましい。
また、本発明の好ましい態様によれば、境界部における時間変化についての所定の演算は、小領域内において境界部の内側に含まれる特定領域の時間変化を、境界部の時間変化へ加算すること、または、境界部の時間変化を小領域内における特定領域の時間変化へ加算することが望ましい。
また、本発明の好ましい態様によれば、特徴情報は、小領域における粒子状対象物の数量と、小領域における粒子状対象物の領域の占める割合と、のうちの少なくとも一つであることが望ましい。
また、本発明の好ましい態様によれば、粒子状対象物に応じて、所定の小領域の少なくとも大きさと、形状と、異方性とのいずれか一つを可変とすることが望ましい。
また、本発明の好ましい態様によれば、粒子状対象物は、生物の細胞内における粒子であることが望ましい。
また、本発明の好ましい態様によれば、粒子状対象物である生物の細胞内における粒子は、複数の異なる時系列画像において、少なくとも、消失することと、出現することと、
大きさを変えることとのいずれか一つを行うことが望ましい。
また、本発明の好ましい態様によれば、粒子状対象物は、不特定多数の人間が集合した集団における人間であることが望ましい。
また、第２の本発明によれば、複数の粒子状対象物の画像を時系列順序で撮像する撮像ステップと、撮像された画像内の各画素に対応する各小領域内の粒子状対象物に関する特徴情報を算出する特徴情報算出ステップと、特徴情報の時間変化を検出する時間変化検出ステップと、特徴情報の時間変化に基づいて、各小領域の中で粒子状対象物を含む小領域を選択し、選択された小領域の周辺部分に存在する境界部における時間変化に関する所定の演算を行うことにより、異常領域を特定する因子特定ステップと、を有することを特徴とする粒子群運動解析方法を提供できる。
また、第３の本発明によれば、粒子状対象物の運動を解析する粒子群運動解析システム用のプログラムであり、コンピュータにより読み取り可能なプログラムであって、複数の粒子状対象物の時系列画像を撮像する撮像手段と、撮像された時系列画像内の各画素に対応する各小領域内の粒子状対象物に関する特徴情報を算出する特徴情報算出手段と、特徴情報の時間変化を検出する時間変化検出手段と、特徴情報の時間変化に基づいて、各小領域の中で粒子状対象物を含む小領域を選択し、選択された小領域の周辺部分に存在する境界部における時間変化に関する所定の演算を行うことにより、異常領域を特定する因子特定手段としてコンピュータを機能させることを特徴とするプログラムを提供できる。
また、別の本発明によれば、複数の粒子状対象物の時系列画像を撮像する撮像手段と、撮像された時系列画像内の各画素に対応する各小領域内の前記粒子状対象物に関する特徴情報を算出する特徴情報算出手段と、特徴情報の時間変化を検出する時間変化検出手段と、時間変化の原因となった粒子状対象物を含む領域または時間変化の原因となった粒子状対象物の近傍の領域である異常領域を特定する因子特定手段とを有し、特徴情報は、小領域における粒子状対象物の領域の占める割合であることを特徴とする粒子群運動解析システムを提供できる。
本発明に係る粒子群運動解析システムでは、撮像手段は、複数の粒子状対象物の時系列画像を撮像する。特徴情報算出手段は、撮像された時系列画像内の各画素に対応する各小領域内の粒子状対象物に関する特徴情報を算出する。時間変化検出手段は、特徴情報の時間変化を検出する。そして、因子特定手段は、時間変化の原因となった粒子状対象物を含む領域、または時間変化の原因となった粒子状対象物の近傍の領域である異常領域を特定する。かかる構成において、一定方向に移動する複数の粒子状対象物のうち、他と異なる動きをする粒子状対象物が存在しているときを考える。例えば、特徴情報として粒子の数量を用いると、時間変化検出手段は、粒子の数量のダイバージェンスを算出する。そして、因子特定手段は、所定の演算処理に従って、ダイバージェンスの変化がより大きい異常領域を特定する。異常領域内またはその近傍には、他と異なる動きをする粒子状対象物が存在している。これにより、複数の粒子状対象物の中から、他の粒子状対象物とは異なる動きをする粒子状対象物を検出できるという効果を奏する。
【図面の簡単な説明】
【０００４】
図１は、実施例１に係る粒子群運動解析システムの機能ブロックを示す図である。図２は、顕微鏡ユニットで観察される画像を示す図である。図３は、粒子群運動解析方法の手順を示すフローチャートである。図４は、小領域Ｒ（ｉ，ｊ）の構成を示す図である。図５Ａは、粒子が移動する前の状態を示す図である。図５Ｂは、粒子が移動した後の状態を示す図である。図６Ａは、ダイバージェンスを示す図である。図６Ｂは、投票後の値を示す図である。図７は、異常粒子が存在する状態を示す図である。図８Ａは、異常粒子が存在するときのダイバージェンスを示す図である。図８Ｂは、異常粒子が存在するときの投票後の値を示す図である。図９Ａ〜図９Ｃは、特徴情報の例を示す図である。図１０は、小領域Ｒ（ｉ，ｊ）の他の構成を示す図である。図１１Ａ〜図１１Ｅは、リンキング処理を説明する図である。図１２は、実施例２に係る粒子群運動解析システムの機能ブロックを示す図である。図１３は、監視カメラで観察される時系列画像を示す図である。
【発明を実施するための最良の形態】
【０００５】
以下に、本発明に係る粒子群運動解析システム、粒子群運動解析方法及びプログラムの実施例を図面に基づいて詳細に説明する。なお、この実施例により、この発明が限定されるものではない。
【実施例１】
【０００６】
図１は、本発明の実施例１に係る粒子群運動解析システム１００の機能ブロックを示す。顕微鏡ユニット１０１は、例えば、蛍光観察用の顕微鏡である。顕微鏡ユニット１０１により、生物細胞の蛍光観察を行うことができる。撮像部１０２は、細胞内の複数の粒子の時系列画像を撮像する。時系列画像とは、運動している粒子状対象物を時間的順序に従って撮像した画像をいう。
図２は、顕微鏡ユニット１０１で観察される時系列画像２００を示す。時系列画像２００において、細胞ＣＥ内には複数の粒子Ｐ１、Ｐ２・・・Ｐｎ（以下、適宜「粒子Ｐｎ」という。）が存在している。粒子Ｐｎは、粒子状対象物に対応し、様々な大きさの球形状の粒子である。そして、各粒子Ｐｎは、ランダムに動いている。
図１に戻って説明を続ける。記憶部１０３は、撮像部１０２で撮像された複数の時系列画像２００を記憶する。特徴情報算出部１０４は、撮像された時系列画像２００内の各画素に対応する各小領域Ｒ（ｉ，ｊ）における粒子Ｐｎに関する特徴情報を算出する。以下、特徴情報の例である粒子Ｐｎの数量を用いて説明を行う。なお、「特徴情報」の他の例の詳細は後述する。
時間変化検出部１０５は、粒子Ｐｎの数量の時間変化、例えばダイバージェンス（発散）を算出する。因子特定部１０６は、ダイバージェンスの変化の原因となった粒子Ｐｎを含む領域またはその近傍の領域である異常領域を特定する。制御部１０７は、システム全体の制御を行う。また、コンピュータ１０８は、入出力制御、演算処理、時系列画像や解析結果の表示処理等を行う。
粒子群運動解析システム１００のハードウエアのためのプログラムとして、粒子状対象物の運動を解析する粒子群運動解析システム用のプログラムであり、コンピュータにより読み取り可能なプログラムであって、複数の粒子状対象物の時系列画像を撮像する撮像手段と、撮像された時系列画像内の各画素に対応した各小領域内の粒子状対象物に関する特徴情報を算出する特徴情報算出手段と、特徴情報の時間変化を検出する時間変化検出手段と、時間変化の原因となった粒子状対象物を含む領域または粒子状対象物の近傍の領域である異常領域を特定する因子特定手段としてコンピュータを機能させることを特徴とするプログラムを用いることが望ましい。
特徴情報算出部１０４、時間変化検出部１０５、因子特定部１０６の機能は、例えば、ＣＰＵを用いて実現できる。なお、特徴情報算出部１０４、時間変化検出部１０５、因子特定部１０６の機能を情報記憶媒体１０９からコンピュータ１０８にプログラムを読み取らせて実現することも可能である。
また、情報記憶媒体１０９としては、例えば、フレキシブルディスク、ＣＤ−ＲＯＭ、光磁気ディスク、ＩＣカード、ＲＯＭカートリッジ、パンチカードやバーコード等の符号が印刷された印刷物、コンピュータの内部記憶装置（ＲＡＭやＲＯＭ等のメモリ）及び外部記憶装置等、コンピュータによって読み取ることができる種々の媒体を用いることができる。そのプログラムの読み取り方式は接触方式であっても、非接触方式であっても良い。
さらに、情報記憶媒体１０９に代えて、上述した各機能を実現するためにプログラムを、伝送路を介してホスト装置等からダウンロードすることによって上述した各機能を実現することも可能である。
図３は、粒子群運動解析の手順を示すフローチャートである。ステップＳ３０１において、撮像部１０２は、複数の粒子Ｐｎを含む細胞ＣＥの時系列画像を撮像する。記憶部１０３は、撮像された複数の時系列画像を記憶する。ステップＳ３０２において、特徴情報算出部１０４は、複数の時系列画像２００に基づいて特徴情報を算出する。ステップＳ３０３において、時間変化検出部１０５は、特徴情報の時間変化、例えばダイバージェンスを検出する。ステップＳ３０４において、因子特定部１０６は、ダイバージェンスの変化の原因となった粒子Ｐｎを含む領域またはその近傍の領域である異常領域を特定する。
次に、粒子群運動解析システム１００による解析手順を具体例に基づいてさらに説明する。まず、図２で示した時系列画像２００において、検出対象領域を設定する。検出対象領域は、粒子Ｐｎの運動を解析する領域である。次に、検出対象領域を直交格子状の画素Ｉ（ｉ，ｊ）に分割する。そして、画像中の検出対象となっている領域の各画素Ｉ（ｉ，ｊ）について、所定の小領域Ｒ（ｉ，ｊ）を定める。また、検出対象領域の画素数と同じ数の投票箱Ｖ（ｉ，ｊ）を準備する。なお、任意の一つの画素Ｉ（ｉ，ｊ）に着目し、この画素Ｉ（ｉ，ｊ）を適宜「注目画素Ｉ（ｉ，ｊ）」という。注目画素Ｉ（ｉ，ｊ）は、特定領域に対応する。
図４は、小領域Ｒ（ｉ，ｊ）を示す。小領域Ｒ（ｉ，ｊ）は、注目画素Ｉ（ｉ，ｊ）を中心とする３画素×３画素の大きさの正方形領域である。また、図４において、注目画素Ｉ（ｉ，ｊ）に隣接する周囲の８個の画素を斜線を付して示す境界部Ｂ（ｉ，ｊ）という。
（ダイバージェンスΔＳ（ｉ，ｊ）の説明）
図５Ａは、検出対象領域の画素Ｉ（ｉ，ｊ）と粒子Ｐｎとの関係を示す。なお、以下、簡単のため、粒子Ｐｎは全て同一の大きさの円形として記載する。図５Ａにおいて、斜線を付した注目画素Ｉ（ｉ，ｊ）内に記載されている数字「３」は、小領域Ｒ（ｉ，ｊ）内に存在する粒子Ｐｎの数量の総和を示している。
図５Ｂは、全ての粒子Ｐｎが、図５Ｂにおいて矢印で示すように右方向へ１画素分だけ移動している定常状態における、検出対象領域の画素Ｉ（ｉ，ｊ）と粒子Ｐｎとの関係を示す。「定常状態」とは、検出対象領域において、粒子Ｐｎの数量に増減が無い状態という。換言すると、例えば、定常状態とは、複数の粒子Ｐｎがランダムに移動している状態、または全ての粒子Ｐｎが一定の方向へ常に移動している状態である。そして、「定常状態」とは異なる状態、例えば他の粒子と異なる動きをしている粒子を「異常粒子」という。なお、理解の容易のため、複数の粒子Ｐｎのうち粒子Ｐ１〜Ｐ４について付し矢印を付し、他の粒子については矢印を省略する。
図５Ｂにおいて、図５Ａと同様に斜線を付して示す注目画素Ｉ（ｉ，ｊ）に着目する。そして、注目画素Ｉ（ｉ，ｊ）に対応する小領域Ｒ（ｉ，ｊ）内に存在する粒子Ｐｎの数量の総和値を算出する。例えば、注目画素Ｉ（ｉ，ｊ）について、小領域Ｒ（ｉ，ｊ）内に存在する粒子Ｐｎの数量の総和値は「３」である。
図６Ａは、注目画素Ｉ（ｉ，ｊ）についてのダイバージェンス（発散）ΔＳ（ｉ，ｊ）を示す。ダイバージェンスΔＳ（ｉ，ｊ）は、粒子Ｐｎが移動した後の小領域Ｒ（ｉ，ｊ）内の粒子数の総和値から、粒子Ｐｎが移動する前の小領域Ｒ（ｉ，ｊ）内の粒子数の総和値を減算した値である。ダイバージェンスΔＳ（ｉ，ｊ）は、時間変化に対応する。
例えば、注目画素Ｉ（ｉ，ｊ）のダイバージェンスΔＳ（ｉ，ｊ）は、次式で示される。
ΔＳ（ｉ，ｊ）＝３（移動後）−３（移動前）
＝０
なお、ダイバージェンスΔＳ（ｉ，ｊ）のマイナスの符号は、移動後に総和値が減少したことを示す。また、ダイバージェンスΔＳ（ｉ，ｊ）のプラスの符号は、移動後に総和値が増加したことを示す。
さらに好ましくは、複数のダイバージェンスΔＳ（ｉ，ｊ）を時間平均することが望ましい。これにより、ノイズの影響を低減できる。
（境界投票の説明）
次に、以下の式（１）〜（８）に示す手順により、境界部Ｂ（ｉ，ｊ）の各画素に対して、注目画素Ｉ（ｉ，ｊ）のダイバージェンスΔＳ（ｉ，ｊ）を投票（ｖｏｔｉｎｇ）する。但し、ダイバージェンスΔＳ（ｉ，ｊ）を投票する前の投票箱をＶ’（ｉ，ｊ）とする。「投票」は、例えば、特徴情報に関する所定の２次元ヒストグラムを算出する手順に対応する。
Ｖ（ｉ−１，ｊ−１）＝Ｖ’（ｉ−１，ｊ−１）＋ΔＳ（ｉ，ｊ）・・・（１）
Ｖ（ｉ，ｊ−１）    ＝Ｖ’（ｉ，ｊ−１）    ＋ΔＳ（ｉ，ｊ）・・・（２）
Ｖ（ｉ＋１，ｊ−１）＝Ｖ’（ｉ＋１，ｊ−１）＋ΔＳ（ｉ，ｊ）・・・（３）
Ｖ（ｉ＋１，ｊ）    ＝Ｖ’（ｉ＋１，ｊ）    ＋ΔＳ（ｉ，ｊ）・・・（４）
Ｖ（ｉ＋１，ｊ＋１）＝Ｖ’（ｉ＋１，ｊ＋１）＋ΔＳ（ｉ，ｊ）・・・（５）
Ｖ（ｉ，ｊ＋１）    ＝Ｖ’（ｉ，ｊ＋１）    ＋ΔＳ（ｉ，ｊ）・・・（６）
Ｖ（ｉ−１，ｊ＋１）＝Ｖ’（ｉ−１，ｊ＋１）＋ΔＳ（ｉ，ｊ）・・・（７）
Ｖ（ｉ−１，ｊ）    ＝Ｖ’（ｉ−１，ｊ）    ＋ΔＳ（ｉ，ｊ）・・・（８）
図４に斜線を付して示す境界部Ｂ（ｉ，ｊ）の各画素に対して、上記式（１）〜（８）に示す加算を行う。これにより、図４において、境界部Ｂ（ｉ，ｊ）の各画素に対して右回りに投票を行うことになる。そして、上記式（１）〜（８）に示す投票を検出対象領域の全ての画素について行う。
全ての画素について上記式（１）〜（８）に示す投票をした後の新たな投票箱Ｖ（ｉ，ｊ）の値を次式（９）に示す。
Ｖ（ｉ，ｊ）＝ΔＳ（ｉ−１，ｊ−１）＋ΔＳ（ｉ，ｊ−１）
＋ΔＳ（ｉ＋１，ｊ−１）
＋ΔＳ（ｉ＋１，ｊ）
＋ΔＳ（ｉ＋１，ｊ＋１）＋ΔＳ（ｉ，ｊ＋１）
＋ΔＳ（ｉ−１，ｊ＋１）
＋ΔＳ（ｉ−１，ｊ）・・・（９）
最終的に投票箱Ｖ（ｉ，ｊ）を得るためには、以下の（Ａ）、（Ｂ）のいずれの演算方法でも良い。
（Ａ）上記式（１）〜（８）に従ってダイバージェンスΔＳ（ｉ，ｊ）を境界部Ｂ（ｉ，ｊ）に投票する方法（以下、「順投票」という。）
（Ｂ）上記式（９）に従って、全ての画素についてダイバージェンスΔＳ（ｉ，ｊ）を計算した後に投票箱Ｖ（ｉ，ｊ）を算出する方法（以下、「逆投票」という。）
順投票では、演算に必要なメモリ容量が小さくて良いという利点を有する。また、逆投票では、計算速度が速いという利点を有する。
図６Ｂは、例えば、逆投票より得られた投票箱Ｖ（ｉ，ｊ）の結果を示す。図５Ｂに示したように、全ての粒子Ｐｎは定常状態にあり、一定の方向、本例では右方向に移動している。このため、投票箱Ｖ（ｉ，ｊ）の値は、一定の範囲０〜２内の大きさである。
次に、他の粒子とは異なる動きをする異常粒子が存在する場合について説明する。図７は、異常粒子Ｐａｂが存在するときの、検出対象領域の画素Ｉ（ｉ，ｊ）と粒子Ｐｎとの関係を示す。図７において、異常粒子Ｐａｂを除く他の粒子Ｐｎは矢印で示すように右方向へ１画素分だけ移動している。これに対して、三角形で示した異常粒子Ｐａｂは、図７の左方向へ１画素分だけ移動する。
図８Ａは、異常粒子Ｐａｂが存在している画素におけるダイバージェンスΔＳ（ｉ，ｊ）を算出した結果を示す。また、図８Ｂは、例えば、逆投影により得られた投票箱Ｖ（ｉ，ｊ）の値を示す。図８Ｂから明らかなように、異常粒子Ｐａｂが存在している画素では、投票箱Ｖ（ｉ，ｊ）の値（＝６）が、他の画素における投票箱Ｖ（ｉ，ｊ）の値よりも顕著に大きい。この結果、投票箱Ｖ（ｉ，ｊ）の値に基づいて、異常粒子Ｐａｂを含む小領域Ｒ（ｉ，ｊ）を検出できる。
顕微鏡観察、とくに共焦点（コンフォーカル）顕微鏡による観察では、粒子Ｐｎが焦点深度方向に上下運動することで、観察視野内において消失することと、出現することと、大きさを変えることとのいずれか一つを行うことが多い。本実施例では、このような動きをする粒子を高い感度で検出できる。
特に、複数の時系列画像のいずれかにおいて、着目している粒子Ｐｎが消失したときは、いわゆるリンキング（Ｌｉｎｋｉｎｇ）処理を行うことが望ましい。図１１Ａ〜図１１Ｅは、リンキング処理を説明する図である。図１１Ａは、時刻ｔ０で撮像された粒子Ｐｔｏを黒丸印で示す。図１１Ｂは、粒子Ｐｔ０が移動して、時刻ｔ１で撮像された粒子Ｐｔ１を黒丸印で示す。そして、図１１Ｃにおいて、時刻ｔ２では粒子が消失している。そして、図１１Ｄでは、時刻ｔ３で再び粒子Ｐｔ３が黒丸印で示すように現れている。このとき、図１１Ｅに示すように、粒子Ｐｔ２（破線の三角形で示す）の検出が不確定である時刻ｔ２における位置を時間的に補間して推定する。このように、時間補間により同じ粒子であると同定（推定）できる。これにより、粒子が突然消失するようなときでも、高感度で異常粒子を検出できる
また、小領域Ｒ（ｉ，ｊ）の大きさは、３画素×３画素の正方形形状に限られない。例えば、粒子Ｐｎの大きさ、数量、存在密度、移動速度ベクトル分布、単位時間当たりの移動量の最大値等に応じて、小領域Ｒ（ｉ，ｊ）の少なくとも大きさと、形状と、異方性とのいずれか一つを可変とすることができる。これにより、異常粒子を最適化した条件で高感度で検出できる。また、リアルタイム処理においては、前フレームの処理結果をもとに小領域Ｒ（ｉ，ｊ）の大きさ、形状、異方性を変更しても良い。なお、「異方性」とは、注目画素Ｉ（ｉ，ｊ）が小領域Ｒ（ｉ，ｊ）の中心位置からシフトした位置に設定されていることをいう。例えば、長方形の小領域Ｒ（ｉ，ｊ）に対して、注目画素Ｉ（ｉ，ｊ）が長方形の短辺近傍に存在する場合である。
上述の説明では、一個ずつの粒子Ｐｎの数量をカウントする例を用いている。粒子Ｐｎのカウントについて、具体的に説明する。まず、蛍光された細胞を明るさで２値化（閾値処理）する。２値化された画像についてラベリング処理を行う。これにより、孤立領域を検出して、孤立領域の数量をカウントする。例えば、図９Ａに示すように、粒子Ｐｎが孤立しているときは、孤立領域の数量は画素の数量に対応する。
また、図９Ｂに示すように、２つの粒子が近接、重複している場合もある。このときは、斜線を付して示す２つの画素が孤立領域に対応する。また、粒子Ｐｎのサンプルの画像をテンプレートとしてマッチング処理を行うこともできる。このときは、マッチング処理結果の類似度のピークをカウントする。さらに、例えば、図９Ｃに示すように、１個の粒子Ｐｎが複数の画素にわたる大きさの場合もある。この場合は、まず、蛍光された細胞を明るさで２値化する。次に、小領域Ｒ（ｉ，ｊ）内で閾値以上の領域をカウントする。そして、小領域Ｒ（ｉ，ｊ）に対する粒子Ｐｎの面積の割合を特徴情報とする。
また、特徴情報の他の例を説明する。蛍光された細胞のとき、各画素における粒子Ｐｎの明るさのそのままの値（輝度値）、または上述のテンプレートによるマッチング処理結果の類似度を算出する。次に、算出された明るさ、または類似度を小領域Ｒ（ｉ，ｊ）内において積分する。そして、積分された値を粒子Ｐｎの平均的な大きさで除算する。この除算された結果を特徴情報として用いることもできる。
本実施例によれば、小領域Ｒ（ｉ，ｊ）（＝検出フィルタ）を設定し、ダイバージェンスΔＳ（ｉ，ｊ）と、順投影または逆投影とを行うことにより、他の粒子Ｐｎとは異なる動きの異常粒子Ｐａｂを高感度に検出できる。なお、本実施例では、境界部Ｂ（ｉ，ｊ）を小領域Ｒ（ｉ，ｊ）の内側に設定している。しかしながら、これに限られず、境界部Ｂ（ｉ，ｊ）を小領域Ｒ（ｉ，ｊ）の外側に設定すること、及び境界部Ｂ（ｉ，ｊ）を小領域Ｒ（ｉ，ｊ）の外側と内側とに半分ずつの画素にまたがるように設定することのいずれの構成とすることもできる。
（実施例１の変形例）
次に、上記実施例１の変形例について説明する。上記実施例１では、注目画素Ｉ（ｉ，ｊ）に隣接する８個の画素を境界部Ｂ（ｉ，ｊ）としている。これに対して、本変形例では、注目画素Ｉ（ｉ，ｊ）に離間して存在する１６個の画素を境界部Ｂ（ｉ，ｊ）とする。
図１０は、本変形例における小領域Ｒ（ｉ，ｊ）を示す。小領域Ｒ（ｉ，ｊ）は、注目画素Ｉ（ｉ，ｊ）を中心とする５画素×５画素の大きさの正方形形状の領域である。そして、境界部Ｂ（ｉ，ｊ）は、図１０で斜線を付して示すような、注目画素Ｉ（ｉ，ｊ）に対して離間して存在する１６個の画素である。
そして、実施例１と同様に、注目画素Ｉ（ｉ，ｊ）についてのダイバージェンスΔＳ（ｉ，ｊ）を算出する。ダイバージェンスΔＳ（ｉ，ｊ）は、粒子Ｐｎが移動した後の小領域Ｒ（ｉ，ｊ）内の粒子数の総和値から、粒子Ｐｎが移動する前の小領域Ｒ（ｉ，ｊ）内の粒子数の総和値を減算した値である。
次に、以下の式（１０）〜（２５）に示す手順により、境界部Ｂ（ｉ，ｊ）の各画素に対して、注目画素Ｉ（ｉ，ｊ）のダイバージェンスΔＳ（ｉ，ｊ）を投票（ｖｏｔｉｎｇ）する。但し、ダイバージェンスΔＳ（ｉ，ｊ）を投票する前の投票箱をＶ’（ｉ，ｊ）とする。
Ｖ（ｉ−２，ｊ−２）＝Ｖ’（ｉ−２，ｊ−２）＋ΔＳ（ｉ，ｊ）
・・・（１０）
Ｖ（ｉ−１，ｊ−２）＝Ｖ’（ｉ−１，ｊ−２）＋ΔＳ（ｉ，ｊ）
・・・（１１）
Ｖ（ｉ，ｊ−２）    ＝Ｖ’（ｉ，ｊ−２）    ＋ΔＳ（ｉ，ｊ）
・・・（１２）
Ｖ（ｉ＋１，ｊ−２）＝Ｖ’（ｉ＋１，ｊ−２）＋ΔＳ（ｉ，ｊ）
・・・（１３）
Ｖ（ｉ＋２，ｊ−２）＝Ｖ’（ｉ＋２，ｊ−２）＋ΔＳ（ｉ，ｊ）
・・・（１４）
Ｖ（ｉ＋２，ｊ−１）＝Ｖ’（ｉ＋２，ｊ−１）＋ΔＳ（ｉ，ｊ）
・・・（１５）
Ｖ（ｉ＋２，ｊ）    ＝Ｖ’（ｉ＋２，ｊ）    ＋ΔＳ（ｉ，ｊ）
・・・（１６）
Ｖ（ｉ＋２，ｊ＋１）＝Ｖ’（ｉ＋２，ｊ＋１）＋ΔＳ（ｉ，ｊ）
・・・（１７）
Ｖ（ｉ＋２，ｊ＋２）＝Ｖ’（ｉ＋２，ｊ＋２）＋ΔＳ（ｉ，ｊ）
・・・（１８）
Ｖ（ｉ＋１，ｊ＋２）＝Ｖ’（ｉ＋１，ｊ＋２）＋ΔＳ（ｉ，ｊ）
・・・（１９）
Ｖ（ｉ，ｊ＋２）    ＝Ｖ’（ｉ，ｊ＋２）    ＋ΔＳ（ｉ，ｊ）
・・・（２０）
Ｖ（ｉ−１，ｊ＋２）＝Ｖ’（ｉ−１，ｊ＋２）＋ΔＳ（ｉ，ｊ）
・・・（２１）
Ｖ（ｉ−２，ｊ＋２）＝Ｖ’（ｉ−２，ｊ＋２）＋ΔＳ（ｉ，ｊ）
・・・（２２）
Ｖ（ｉ−２，ｊ＋１）＝Ｖ’（ｉ−２，ｊ＋１）＋ΔＳ（ｉ，ｊ）
・・・（２３）
Ｖ（ｉ−２，ｊ）    ＝Ｖ’（ｉ−２，ｊ）    ＋ΔＳ（ｉ，ｊ）
・・・（２４）
Ｖ（ｉ−２，ｊ−１）＝Ｖ’（ｉ−２，ｊ−１）＋ΔＳ（ｉ，ｊ）
・・・（２５）
図１０に示す境界部Ｂ（ｉ，ｊ）の各画素に対して、上記式（１０）〜（２５）に示す加算を行う。これにより、図１０において、境界部Ｂ（ｉ，ｊ）の各画素に対して右回りに投票を行うことになる。そして、上記式（１０）〜（２５）に示す投票を検出対象領域の全ての画素について行う。
全ての画素について上記式（１０）〜（２５）に示す投票後の新たな投票箱Ｖ（ｉ，ｊ）の値を次式（２６）に示す。
Ｖ（ｉ，ｊ）＝ΔＳ（ｉ−２，ｊ−２）＋ΔＳ（ｉ−１，ｊ−２）
＋ΔＳ（ｉ，ｊ−２）
＋ΔＳ（ｉ＋１，ｊ−２）＋ΔＳ（ｉ＋２，ｊ−２）
＋ΔＳ（ｉ＋２，ｊ−１）＋ΔＳ（ｉ＋２，ｊ）
＋ΔＳ（ｉ＋２，ｊ＋１）
＋ΔＳ（ｉ＋２，ｊ＋２）＋ΔＳ（ｉ＋１，ｊ＋２）
＋ΔＳ（ｉ，ｊ＋２）
＋ΔＳ（ｉ−１，ｊ＋２）＋ΔＳ（ｉ−２，ｊ＋２）
＋ΔＳ（ｉ−２，ｊ＋１）＋ΔＳ（ｉ−２，ｊ）
＋ΔＳ（ｉ−２，ｊ−１）
・・・（２６）
このように、本変形例においても上記実施例１と同様に、最終的に投票箱Ｖ（ｉ，ｊ）を得るためには、順投票と逆投票とのいずれの手順でも良い。
【実施例２】
【０００７】
図１２は、本発明の実施例２に係る粒子群運動解析システム３００の機能ブロックを示す。上記実施例１では、顕微鏡ユニットにより、生物細胞の蛍光観察を行っている。これに対して、本実施例では、粒子状対象物として群集における人間を対象とする点が異なる。上記実施例１と同一の部分には同一の符号を付し、重複する説明は省略する。
監視カメラ３０１は、不特定多数の人間が集合する場所、例えば、駅、空港、公共施設等に設置されている。そして、監視カメラ３０１は、駅構内の群集における人間をモニタする。撮像部１０２は、群集の時系列画像を撮像する。なお、図１２では、監視カメラ３０１と撮像部１０２とを別々に示している。しかしながら、監視カメラ３０１が撮像部１０２を含む構成とすることもできる。
図１３は、時系列画像４００を示す。監視カメラ３０１が撮像する視野内には、複数の人間ＰＳ１、ＰＳ２・・・ＰＳｎ（以下、適宜「人間ＰＳｎ」という。）が存在している。人間ＰＳ１、ＰＳ２・・・ＰＳｎは、矢印で示す方向へ歩行している。
粒子群運動解析システム３００は、上記実施例１で述べた運動解析手順と同様の手順により、他の人間と異なる動きをする人間、即ち異常行動をする人間を検出することができる。
また、小領域Ｒ（ｉ，ｊ）の面積を大きくすることにより、群集内の人間の大まかな流れの時間変化を検出できる。これにより、駅における混雑を解消するために有用な群集行動解析を行うことができる。
また、本実施例において、人間ＰＳｎの頭部や頭髪は、形状や大きさ等が予め略決まっている。従って、人間ＰＳｎの頭部や頭髪のサンプル画像をテンプレートにしておくこともできる。そして、テンプレートを用いて、マッチング処理を行う。これにより、人間ＰＳｎを同定することができる。特定の人間ＰＳｎを同定した後は、同定された人間ＰＳｎをトラッキング（追跡）することも可能である。
本実施例によれば、群集内において他の人間と異なる動きをする人間を検出できる。また、上述のように、大まかな人間の流れを検出することで、例えば駅構内における警備員の効率的な配置に有用な情報を提供できる。
【実施例３】
【０００８】
次に、本発明の実施例３について説明する。上記実施例においては、異常粒子Ｐａｂが存在している画素では、投票箱Ｖ（ｉ，ｊ）の値が、他の画素における投票箱Ｖ（ｉ，ｊ）の値よりも顕著に大きい。これにより、投票箱Ｖ（ｉ，ｊ）の値に基づいて、異常粒子Ｐａｂを含む小領域Ｒ（ｉ，ｊ）を検出できる。
ここで、投票箱Ｖ（ｉ，ｊ）の値を絶対値として扱うこと、例えば、最大値が１に近いときは異常粒子が存在している状態、最大値が１に近くなければ、単なるノイズである状態という判断ができると、さらに利便性が向上する。
このため、本実施例では、後述する正規化するステップをさらに有している。また、本実施例は、小領域Ｒ（ｉ，ｊ）の形状、時間的要素には全く限定されず、様々な形状、大きさの小領域Ｒ（ｉ，ｊ）を設定することができる。
次式（２７）は、小領域Ｒ（ｉ，ｊ）に粒子が存在する割合を示している。なお、上述したように、小領域Ｒ（ｉ，ｊ）の形状は任意の形状とすることができる。このため、座標（ｉ，ｊ）の代わりに任意の位置（ｘ，ｙ）を用いて説明を行う。
時間ｔ、位置（ｘ，ｙ）において小領域Ｒ（ｘ，ｙ）に粒子が存在している割合ρ（ｘ，ｙ，ｔ）は、次式（２７）のように表される。
0004934018.tif 000002
ここで、Ｐ（ｘ，ｙ，ｔ）は、確率関数である（０≦Ｐ≦１）。Ｐ（ｘ，ｙ，ｔ）は、イメージデータに対応している。
また、小領域Ｒ（ｘ，ｙ）に存在している粒子の数量ｎ（ｘ，ｙ，ｔ）は、次式（２８）で示すことができる。
0004934018.tif 000003
ここで、ρｏは、小領域Ｒ（ｘ，ｙ）において一つの粒子が占める割合の平均値である。例えば、一つの粒子が占める面積が２画素×２画素の粒子を探すときは、
ρｏ＝２×２＝４
となる。
そして、小領域Ｒ（ｘ，ｙ）における粒子数の変化は、次式（２９）で表すことができる。なお、Ｋ（ｘ，ｙ，ｔ）は、上記実施例におけるダイバージェンスΔＳ（ｘ，ｙ，ｔ）に対応する。
0004934018.tif 000004
式（２９）において、Ｋ＝０のとき、定常状態となっている。これに対して、Ｋ≠０のとき、小領域Ｒ（ｘ、ｙ）にＫ個の異常粒子が存在している。そして、ダイバージェンス理論によれば、Ｋ≠０のとき小領域Ｒ（ｘ、ｙ）の境界において何かが起きているということが言える。しかしながら、Ｋ≠０に関連している境界画素の正確な位置まではわからない。
そこで、本実施例では、次式（３０）で示すような正規化された関数Ｍ（ｘ，ｙ）を用いる。
0004934018.tif 000005
関数Ｍは、投票される位置または投票する位置においては、所定の値を有している。これに対して、投票しない位置においては、ゼロである。そして、例えば注目画素位置・時間を（ｘ’，ｙ’，ｔ’）とすると、関数Ｍ（ｘ，ｙ）の値が大きいほど、特定の位置・時間（ｘ’＋ｘ，ｙ’＋ｙ，ｔ’）に、Ｋ（ｘ’，ｙ’，ｔ’）個の異常粒子のうち、いずれかが存在している傾向がある。また、小領域Ｒ（ｉ、ｊ）が大きくなると、関数Ｍの値は小さくなる。
そして、時刻ｔ、位置（ｘ，ｙ）に存在する異常粒子の平均個数Ｕは、次式（３１）で表すことができる。
0004934018.tif 000006
式（３１）は、上述した順投票の手順に相当する。本願の発明者らは、本実施例で説明したアルゴリズムにより画像処理する実験を、コンフォーカル顕微鏡により得られたビデオ画像に対して行っている。一つの異常粒子を含む時間ｔ＝１〜ｔ＝７の７つの画像を用いる。粒子を５１２画素×５１２画素のオリジナルイメージの中から、１２０画素×１２０画素のイメージ領域を抽出する。小領域Ｒ（ｉ、ｊ）の大きさは、５画素×５画素とする。また、小領域Ｒ（ｉ、ｊ）における重み（ｗｅｉｇｈｔ）は、各画素に対して１／２５である。
単にＫを演算するときのみに比較して、さらにＵまで演算すると、異常粒子の軌跡がより明確に認識でき、また、ノイズも軽減できる。なお、この画像処理には、１．１ＧＨｚのＣＰＵを用いたとき、約０．５秒ほどかかった。
本実施例によれば、絶対値により異常粒子の存在を判断できる、また、ノイズを低減することもできる。
また、本発明は、その趣旨を逸脱しない範囲で様々な変形例をとることができる。
【産業上の利用可能性】
【０００９】
以上のように、本発明に係る粒子群運動解析システムは、細胞内の粒子の動きや群集内の人間の動き、特に他とは異なる動きをする粒子等を検出するときに有用である。
【図１】
0004934018.tif 000007
【図２】
0004934018.tif 000008
【図３】
0004934018.tif 000009
【図４】
0004934018.tif 000010
【図５Ａ】
0004934018.tif 000011
【図５Ｂ】
0004934018.tif 000012
【図６Ａ】
0004934018.tif 000013
【図６Ｂ】
0004934018.tif 000014
【図７】
0004934018.tif 000015
【図８Ａ】
0004934018.tif 000016
【図８Ｂ】
0004934018.tif 000017
【図９Ａ】
0004934018.tif 000018
【図９Ｂ】
0004934018.tif 000019
【図９Ｃ】
0004934018.tif 000020
【図１０】
0004934018.tif 000021
【図１１Ａ】
0004934018.tif 000022
【図１１Ｂ】
0004934018.tif 000023
【図１１Ｃ】
0004934018.tif 000024
【図１１Ｄ】
0004934018.tif 000025
【図１１Ｅ】
0004934018.tif 000026
【図１２】
0004934018.tif 000027
【図１３】
0004934018.tif 000028