JP2024052595A - 核医学診断装置、データ処理方法及びプログラム - Google Patents

Abstract

【課題】画質を向上させること。【解決手段】 実施形態に係る核医学診断装置は、取得部と、算出部と、特定部とを備える。取得部は、第１の検出器で検出される第１の光子数情報と、前記第１の検出器とは異なる第２の検出器で検出される第２の光子数情報とを取得する。算出部は、前記第１の光子数情報に基づいて前記第１の検出器に対応する第１の発光確率モデルと、前記第２の光子数情報に基づいて前記第２の検出器に対応する第２の発光確率モデルとを算出する。特定部は、前記第１の発光確率モデルと前記第２の発光確率モデルとに基づいて、前記第１の検出器と前記第２の検出器とによって規定されるＬＯＲ（Ｌｉｎｅ Ｏｆ Ｒｅｓｐｏｎｓｅ）の検出時間差に関する確率分布モデルを特定する。【選択図】図１

Description

本明細書及び図面に開示の実施形態は、核医学診断装置、データ処理方法及びプログラムに関する。

ＰＥＴ画像の再構成を行う際において、ＴｏＦ（Ｔｉｍｅ Ｏｆ Ｆｌｉｇｈｔ）カーネルを用いて画像再構成が行われる。ここで、画像再構成に用いられるＴｏＦカーネルの時間幅が、検出器のＴｏＦ性能に整合するものでない場合、再構成後の画像の画質が低下する場合がある。従って、検出器のＴｏＦ性能と、画像再構成に用いるＴｏＦカーネルとが、できるだけ一致することが望ましい。

また、再構成画像の画質を高画質にするための一つの方法としては、各検出器の性能を均一にし、検出器ごとの性能のばらつきを小さくすることが考えられるが、コスト面での負担が大きい場合がある。

Eric Berg and Simon R Cherry, "Using convolutional neural networks to estimate time-of-flight from PET detector waveforms", vol. 63, no 2, Phys. Med. Biol, 2018年

本明細書及び図面の開示の実施形態が解決しようとする課題の一つは、画質を向上させることである。ただし、本明細書及び図面に開示の実施形態により解決しようとする課題は上記課題に限られない。後述する実施形態に示す各構成による各効果に対応する課題を他の課題として位置づけることもできる。

実施形態に係る核医学診断装置は、取得部と、算出部と、特定部とを備える。取得部は、第１の検出器で検出される第１の光子数情報と、前記第１の検出器とは異なる第２の検出器で検出される第２の光子数情報とを取得する。算出部は、前記第１の光子数情報に基づいて前記第１の検出器に対応する第１の発光確率モデルと、前記第２の光子数情報に基づいて前記第２の検出器に対応する第２の発光確率モデルとを算出する。特定部は、前記第１の発光確率モデルと前記第２の発光確率モデルとに基づいて、前記第１の検出器と前記第２の検出器とによって規定されるＬＯＲ（Ｌｉｎｅ Ｏｆ Ｒｅｓｐｏｎｓｅ）の検出時間差に関する確率分布モデルを特定する。

図１は、第１の実施形態及び第２の実施形態に係る核医学診断装置の構成の一例を示した図である。 図２は、第１の実施形態に係る核医学診断装置が行う処理の流れを示したフローチャートである。 図３は、第１の実施形態に係る核医学診断装置が行う処理について説明した図である。 図４は、図２のステップＳ３００の処理の一例をより詳細に説明したフローチャートである。 図５は、第１の実施形態に係る核医学診断装置が行う処理について説明した図である。 図６は、第１の実施形態に係る核医学診断装置が行う処理について説明した図である。 図７は、第１の実施形態に係る核医学診断装置が行う処理について説明した図である。 図８は、図２のステップＳ３００の処理の別の一例をより詳細に説明したフローチャートである。 図９は、第１の実施形態に係る核医学診断装置が行う処理について説明した図である。 図１０は、第１の実施形態に係る核医学診断装置が行う処理について説明した図である。 図１１は、第１の実施形態に係る計算結果の一例について説明した図である。 図１２は、第１の実施形態に係る計算結果の一例について説明した図である。 図１３は、第２の実施形態に係る核医学診断装置が行う処理の流れを示したフローチャートである。 図１４は、図１３のステップＳ３０の処理の一例をより詳細に説明したフローチャートである。 図１５は、図１４のステップＳ３００の処理の一例をより詳細に説明したフローチャートである。 図１６は、図１４のステップＳ３００の処理の別の一例をより詳細に説明したフローチャートである。 図１７は、第２の実施形態に係る核医学診断装置が行う処理について説明した図である。 図１８は、第２の実施形態に係る方法の検証結果の一例について示した図である。 図１９は、第２の実施形態に係る方法の検証結果の一例について示した図である。 図２０は、第２の実施形態に係る方法の検証結果の一例について示した図である。

（第１の実施形態）
以下、図面を参照しながら、核医学診断装置、データ処理方法及びプログラムの実施形態について詳細に説明する。

図１は、第１の実施形態に係る核医学診断装置としてのＰＥＴ装置１００の構成を示す図である。図１に示すように、第１の実施形態に係るＰＥＴ装置１００は、架台装置１と、コンソール装置２とを備える。架台装置１は、検出器３と、フロントエンド回路１０２と、天板１０３と、寝台１０４と、寝台駆動部１０６とを備える。

検出器３は、被検体に投与され集積した薬剤から放出される陽電子が周囲の組織の電子と対消滅を起こすことで発生したガンマ線が発光体と相互作用することにより励起状態となった物質が再び基底状態に遷移する際に再放出される光であるシンチレーション光（蛍光）を検出することにより、放射線を検出する検出器である。また、実施形態においては、検出器３は、チェレンコフ光も検出することができる。検出器３は、被検体に投与され集積した薬剤から放出される陽電子が周囲の組織の電子と対消滅を起こすことで発生したガンマ線の放射線のエネルギー情報を検出する。検出器３は、被写体Ｐの周囲をリング状に取り囲むように複数配置され、例えば複数の検出器ブロックからなる。

検出器３は、典型的には、シンチレータ結晶と、光検出素子からなる光検出面とからなる。

シンチレータ結晶の材質としては、例えばチェレンコフ光を発生させるのに適した材質、例えば、ゲルマニウム酸ビスマス（ＢＧＯ、Ｂｉｓｍｕｔｈ Ｇｅｒｍａｎｉｕｍ Ｏｘｉｄｅ）や、鉛ガラス（ＳｉＯ₂＋ＰｂＯ）、フッ化鉛（ＰｂＦ₂）、ＰＷＯ（ＰｂＷＯ₄）等の鉛化合物を用いることができる。また、別の例として、例えば、ＬＹＳＯ（Ｌｕｔｅｔｉｕｍ Ｙｔｔｒｉｕｍ Ｏｘｙｏｒｔｈｏｓｉｌｉｃａｔｅ）、ＬＳＯ（Ｌｕｔｅｔｉｕｍ Ｏｘｙｏｒｔｈｏｓｉｌｉｃａｔｅ）、ＬＧＳＯ（Ｌｕｔｅｔｉｕｍ Ｇａｄｏｌｉｎｉｕｍ Ｏｘｙｏｒｔｈｏｓｉｌｉｃａｔｅ）等やＢＧＯ等の、シンチレータ結晶を用いてもよい。光検出面３ｂを構成する光検出素子は、例えば複数の画素からなり、これらの画素それぞれは、例えばＳＰＡＤ（Ｓｉｎｇｌｅ Ｐｈｏｔｏｎ Ａｖａｌａｎｃｈｅ Ｄｉｏｄｅ）で構成される。また、検出器３の構成は、上述の例に限られず、一例として、光検出素子は、例えばＳｉＰＭ（Ｓｉｌｉｃｏｎ ｐｈｏｔｏｍｕｌｔｉｐｌｉｅｒ）や、光電子増倍管を用いてもよい。

なお、シンチレータ結晶はモノリシック結晶であってもよく、また光検出素子からなる光検出面は、シンチレータ結晶の６面に例えば配置されてもよい。なお、以下の第１の実施形態では、検出器３のシンチレータ結晶は、モノリシック結晶でない場合をまずは例にとり説明する。

また、架台装置１は、フロントエンド回路１０２により、検出器３の出力信号から計数情報を生成し、生成した計数情報を、コンソール装置２の記憶部１３０に格納する。なお、検出器３は、複数のブロックに区分けされ、フロントエンド回路１０２を備える。

フロントエンド回路１０２は、検出器３の出力信号をデジタルデータに変換し、計数情報を生成する。この計数情報には、消滅ガンマ線の検出位置、エネルギー値、及び検出時間が含まれる。例えば、フロントエンド回路１０２は、シンチレーション光を同じタイミングで電気信号に変換した複数の光検出素子を特定する。そして、フロントエンド回路１０２は、消滅ガンマ線が入射したシンチレータの位置を示すシンチレータ番号（Ｐ）を特定する。消滅ガンマ線が入射したシンチレータ位置を特定する手段は、各光検出素子の位置及び電気信号の強度に基づいて重心演算を行うことによって特定してもよい。また、シンチレータと光検出素子の各々の素子サイズが対応している場合には、例えば最大出力が得られた光検出素子に対応するシンチレータを消滅ガンマ線が入射したシンチレータ位置と仮定し、さらにシンチレータ間散乱を考慮して最終的に特定するなどする。

また、フロントエンド回路１０２は、各光検出素子から出力された電気信号の強度を積分計算するあるいは電気信号強度が閾値を超えた時間（Ｔｉｍｅ ｏｖｅｒ Ｔｈｒｅｓｈｏｌｄ）を計測し検出器３に入射した消滅ガンマ線のエネルギー値（Ｅ）を特定する。また、フロントエンド回路１０２は、検出器３によって消滅ガンマ線によるシンチレーション光が検出された検出時間（Ｔ）を特定する。なお、検出時間（Ｔ）は、絶対時刻であってもよいし、撮影開始時点からの経過時間であってもよい。このように、フロントエンド回路１０２は、シンチレータ番号（Ｐ）、エネルギー値（Ｅ）、及び検出時間（Ｔ）を含む計数情報を生成する。

なお、フロントエンド回路１０２は、例えば、ＣＰＵ（Ｃｅｎｔｒａｌ Ｐｒｏｃｅｓｓｉｎｇ Ｕｎｉｔ）、ＧＰＵ（Ｇｒａｐｈｉｃａｌ Ｐｒｏｃｅｓｓｉｎｇ Ｕｎｉｔ）或いは、特定用途向け集積回路（Ａｐｐｌｉｃａｔｉｏｎ Ｓｐｅｃｉｆｉｃ Ｉｎｔｅｇｒａｔｅｄ Ｃｉｒｃｕｉｔ：ＡＳＩＣ）、プログラマブル論理デバイス（例えば、単純プログラマブル論理デバイス（Ｓｉｍｐｌｅ Ｐｒｏｇｒａｍｍａｂｌｅ Ｌｏｇｉｃ Ｄｅｖｉｃｅ：ＳＰＬＤ）、複合プログラマブル論理デバイス（Ｃｏｍｐｌｅｘ Ｐｒｏｇｒａｍｍａｂｌｅ Ｌｏｇｉｃ Ｄｅｖｉｃｅ：ＣＰＬＤ）、及びフィールドプログラマブルゲートアレイ（Ｆｉｅｌｄ Ｐｒｏｇｒａｍｍａｂｌｅ Ｇａｔｅ Ａｒｒａｙ：ＦＰＧＡ））等の回路により実現される。フロントエンド回路１０２は、フロントエンド部の一例である。

天板１０３は、被写体Ｐが載置されるベッドであり、寝台１０４の上に配置される。寝台駆動部１０６は、処理回路１５０の制御機能１０５ｅによる制御の下、天板１０３を移動させる。例えば、寝台駆動部１０６は、天板１０３を移動させることで、被写体Ｐを架台装置１の撮影口内に移動させる。

コンソール装置２は、操作者によるＰＥＴ装置１００の操作を受け付け、ＰＥＴ画像の撮影を制御するとともに、架台装置１によって収集された計数情報を用いてＰＥＴ画像を再構成する。図１に示すように、コンソール装置２は、処理回路１５０と、入力装置１１０と、ディスプレイ１２０と、記憶部１３０とを備える。なお、コンソール装置２が備える各部は、バスを介して接続される。処理回路１５０の詳細については後述する。

入力装置１１０は、ＰＥＴ装置１００の操作者によって各種指示や各種設定の入力に用いられるマウスやキーボード等であり、入力された各種指示や各種設定を、処理回路１５０に転送する。例えば、入力装置１１０は、撮影開始指示の入力に用いられる。

ディスプレイ１２０は、操作者によって参照されるモニター等であり、処理回路１５０による制御の下、被写体の呼吸波形やＰＥＴ画像を表示したり、操作者から各種指示や各種設定を受け付けるためのＧＵＩ（Ｇｒａｐｈｉｃａｌ Ｕｓｅｒ Ｉｎｔｅｒｆａｃｅ）を表示したりする。

記憶部１３０は、ＰＥＴ装置１００において用いられる各種データを記憶する。記憶部１３０は、例えば、メモリで構成され、一例として、ＲＡＭ（Ｒａｎｄｏｍ Ａｃｃｅｓｓ Ｍｅｍｏｒｙ）、フラッシュメモリ（ｆｌａｓｈ ｍｅｍｏｒｙ）等の半導体メモリ素子や、ハードディスク、光ディスク等によって実現される。記憶部１３０は、シンチレータ番号（Ｐ）、エネルギー値（Ｅ）、及び検出時間（Ｔ）が対応づけられた情報である計数情報、同時計数情報の通し番号であるコインシデンスＮｏ．に計数情報の組が対応づけられた同時計数情報、または同時計数情報を集計して得られる射影データ、再構成されたＰＥＴ画像等を記憶する。

処理回路１５０は、取得機能１５０ａ、算出機能１５０ｂ、特定機能１５０ｃ、再構成機能１５０ｄ、制御機能１５０ｅ、受付機能１５０ｆ、画像生成機能１５０ｇ、表示制御機能１５０ｈを有する。なお、取得機能１５０ａ、算出機能１５０ｂ、特定機能１５０ｃの各機能については、後ほど詳しく説明する。

第１の実施形態及び第２の実施形態では、取得機能１５０ａ、算出機能１５０ｂ、特定機能１５０ｃ、再構成機能１５０ｄ、制御機能１５０ｅ、受付機能１５０ｆ、画像生成機能１５０ｇ、表示制御機能１５０ｈにて行われる各処理機能は、コンピュータによって実行可能なプログラムの形態で記憶部１３０へ記憶されている。処理回路１５０はプログラムを記憶部１３０から読み出し、実行することで各プログラムに対応する機能を実現するプロセッサである。換言すると、各プログラムを読み出した状態の処理回路１５０は、図１の処理回路１５０内に示された各機能を有することになる。

なお、図１においては単一の処理回路１５０にて、取得機能１５０ａ、算出機能１５０ｂ、特定機能１５０ｃ、再構成機能１５０ｄ、制御機能１５０ｅ、受付機能１５０ｆ、画像生成機能１５０ｇ、表示制御機能１５０ｈにて行われる処理機能が実現されるものとして説明するが、複数の独立したプロセッサを組み合わせて処理回路１５０を構成し、各プロセッサがプログラムを実行することにより機能を実現するものとしても構わない。換言すると、上述のそれぞれの機能がプログラムとして構成され、１つの処理回路１５０が各プログラムを実行する場合であってもよい。別の例として、特定の機能が専用の独立したプログラム実行回路に実装される場合であってもよい。

上記説明において用いた「プロセッサ」という文言は、例えば、ＣＰＵ（Ｃｅｎｔｒａｌ Ｐｒｏｃｅｓｓｉｎｇ Ｕｎｉｔ）、ＧＰＵ（Ｇｒａｐｈｉｃａｌ Ｐｒｏｃｅｓｓｉｎｇ Ｕｎｉｔ）或いは、特定用途向け集積回路（Ａｐｐｌｉｃａｔｉｏｎ Ｓｐｅｃｉｆｉｃ Ｉｎｔｅｇｒａｔｅｄ Ｃｉｒｃｕｉｔ：ＡＳＩＣ）、プログラマブル論理デバイス（例えば、単純プログラマブル論理デバイス（Ｓｉｍｐｌｅ Ｐｒｏｇｒａｍｍａｂｌｅ Ｌｏｇｉｃ Ｄｅｖｉｃｅ：ＳＰＬＤ）、複合プログラマブル論理デバイス（Ｃｏｍｐｌｅｘ Ｐｒｏｇｒａｍｍａｂｌｅ Ｌｏｇｉｃ Ｄｅｖｉｃｅ：ＣＰＬＤ）、及びフィールドプログラマブルゲートアレイ（Ｆｉｅｌｄ Ｐｒｏｇｒａｍｍａｂｌｅ Ｇａｔｅ Ａｒｒａｙ：ＦＰＧＡ））等の回路を意味する。プロセッサは記憶部１３０に保存されたプログラムを読み出し実行することで機能を実現する。

なお、図１において、取得機能１５０ａ、算出機能１５０ｂ、特定機能１５０ｃ、再構成機能１５０ｄ、制御機能１５０ｅ、受付機能１５０ｆ、画像生成機能１５０ｇ、表示制御機能１５０ｈは、それぞれ、取得部、算出部、特定部、再構成部、制御部、受付部、生成部、表示制御部の一例である。なお、処理回路１５０に代えて、フロントエンド回路１０２が、算出部、取得部、特定部等の処理を担ってもよい。

処理回路１５０は、再構成機能１５０ｄにより、フロントエンド回路１０２から取得したデータに基づいて、ＰＥＴ画像の再構成を行い、画像生成機能１５０ｇにより、画像の生成を行う。処理回路１５０は、制御機能１５０ｅにより、架台装置１及びコンソール装置２を制御することによって、ＰＥＴ装置１００の全体制御を行う。例えば、処理回路１５０は、制御機能１５０ｅにより、ＰＥＴ装置１００における撮影を制御する。また、処理回路１０５は、制御機能１５０ｅにより、寝台駆動部１０６を制御する。処理回路１５０は、受付機能１５０ｆにより、入力装置１１０を通じて、ユーザから情報の入力を受け付ける。また、処理回路１５０は、表示制御機能１５０ｈにより、ＰＥＴ画像その他のデータをディスプレイ１２０に表示させる。また、処理回路１５０は、画像生成機能１５０ｇにより、種々の画像の生成を行う。

続いて、第１の実施形態に係る背景について簡単に説明する。

ＰＥＴ画像の再構成を行う際において、ＴｏＦ（Ｔｉｍｅ Ｏｆ Ｆｌｉｇｈｔ）カーネルを用いて画像再構成が行われる。ＴｏＦカーネルとは、２つの検出器で検出された信号の検出時間差の関数として表された、検出イベントの確率密度関数である。ここで、画像再構成に用いられるＴｏＦカーネルの時間幅が、検出器のＴｏＦ性能に整合するものでない場合、再構成後の画像の画質が低下する場合がある。例えば、検出器のＴｏＦ性能と比較して、ＴｏＦカーネルの時間幅が小さい場合には、画像再構成後の放射性物質の位置の推定位置が、本来の分布よりも局所的に集積してしまう。一方、検出器のＴｏＦ性能と比較して、ＴｏＦカーネルの時間幅が大きい場合には、画像再構成後の放射線物質の位置の推定位置が、本来の分布より周辺部分に分布してしまう。従って、検出器のＴｏＦ性能と、画像再構成に用いるＴｏＦカーネルとが、できるだけ一致することが望ましい。

ここで、再構成画像の画質を高画質にするための一つの方法としては、各検出器の性能を均一にし、検出器ごとの性能のばらつきを小さくすることが考えられる。しかしながら、各検出器の性能を均一にするのは、コスト面での負担が大きい場合がある。そこで、第１の実施形態では、発光モデルに基づいて、検出器毎にＴｏＦカーネルの算出を行う。このことで、検出器ごとの性能のばらつきを、ＴｏＦカーネルを検出器ごとに行うことにより吸収することができ、画質を向上させることができる。

具体的には、第１の実施形態に係る核医学診断装置は、取得部と、算出部と、特定部とを備える。取得部は、第１の検出器で検出される第１の光子数情報と、第１の検出器とは異なる第２の検出器で検出される第２の光子数情報とを取得する。算出部は、第１の光子数情報に基づいて第１の検出器に対応する第１の発光確率モデルと、第２の光子数情報に基づいて第２の検出器に対応する第２の発光確率モデルとを算出する。特定部は、第１の発光確率モデルと第２の発光確率モデルとに基づいて、第１の検出器と第２の検出器とによって規定されるＬＯＲ（Ｌｉｎｅ Ｏｆ Ｒｅｓｐｏｎｓｅ）の検出時間差に関する確率分布モデルを特定する。

また、第１の実施形態に係るデータ処理方法は、第１の検出器で検出される第１の光子数情報と、前記第１の検出器とは異なる第２の検出器で検出される第２の光子数情報とを取得し、第１の光子数情報に基づいて第１の検出器に対応する第１の発光確率モデルと、第２の光子数情報に基づいて第２の検出器に対応する第２の発光確率モデルとを算出し、第１の発光確率モデルと第２の発光確率モデルとに基づいて、第１の検出器と第２の検出器とによって規定されるＬＯＲの検出時間差に関する確率分布モデルを特定する。

また、第１の実施形態に係るプログラムは、第１の検出器で検出される第１の光子数情報と、第１の検出器とは異なる第２の検出器で検出される第２の光子数情報とを取得し、第１の光子数情報に基づいて第１の検出器に対応する第１の発光確率モデルと、第２の光子数情報に基づいて第２の検出器に対応する第２の発光確率モデルとを算出し、第１の発光確率モデルと第２の発光確率モデルとに基づいて、第１の検出器と第２の検出器とによって規定されるＬＯＲの検出時間差に関する確率分布モデルを特定する処理をコンピュータに実行させる。

以下、図２～図１０を用いて、第１の実施形態に係る核医学診断装置が行う処理について、詳しく説明する。

図２に、第１の実施形態に係るＰＥＴ装置１００が行う処理の流れが示されている。以下のフローチャートでは、簡単のため、第１の検出器ｌと第２の検出器ｍに係る確率分布モデル（ＴｏＦカーネル）を生成する場合についてのみ記載している。しかしながら、典型的には、画像再構成においては、すべてのｌとｍの組について処理回路１５０は確率分布モデルを生成した上でステップＳ６００の画像再構成処理を行う。この場合、処理回路１５０は、すべてのｌとｍの組についてステップＳ１００からステップＳ５００の処理を繰り返した上で、ステップＳ６００の処理を行う。

ステップＳ１００において、処理回路１５０は、取得機能１５０ａにより、第１の検出器ｌで検出される第１の光子数情報を取得する。

ここで、第１の光子数情報とは、例えば、当該検出器で検出される標準的な一イベントについて、平均的に検出されるシンチレーション光子数Ｎ_Ｓ;lである。また、別の例として、第１の光子数情報とは、例えば、当該検出器で検出される標準的な一イベントについて、平均的に検出されるチェレンコフ光子数Ｎ_Ｃ;lである。また、別の例として、第１の光子数情報とは、当該検出器で検出される標準的な一イベントについて、平均的に検出されるシンチレーション光子数Ｎ_Ｓ;l及びシンチレーション光子数Ｎ_Ｓ;lである。すなわち、第１の光子数情報とは、当該検出における検出事象に対応するシンチレーション光または／及びチェレンコフ光の光子数の期待値である。

ＰＥＴ装置１００は、所定のキャリブレーションスキャンを行うことで、第１の光子数情報を特定してもよい。具体的には、ＰＥＴ装置１００は、所定のキャリブレーションスキャンを行うことで、第１の検出器ｌで取得されるエネルギースペクトルに基づいて、第１の光子数情報及び第２の光子数情報を取得してもよい。

図３に、検出器３内での発光モデルと、処理回路１５０が取得機能１５０ａにより取得する光子スペクトルの概略が示されている。検出器３のシンチレータ１１にガンマ線１０が入射すると、ガンマ線１０はシンチレータ１１にエネルギーを光電吸収またはコンプトン散乱の形で付与するとともに、高エネルギー電子１２が発生する。曲線１６は、高エネルギー電子１２のエネルギーＥ_ｅを時間ｔの関数として示したものである。高エネルギー電子１２は、シンチレータ１１中を走行するとともにエネルギーを失っていくが、エネルギーが閾値Ｅ_ｔｈより大きい間は、チェレンコフ光１３を周囲に放出する。

曲線１７は、チェレンコフ光の光子密度ｎ_ｃ（ｔ）を、時間ｔのガウス関数に近似したものである。具体的には、チェレンコフ光の光子密度ｎ_ｃ（ｔ）は、Ｎ_ｃをチェレンコフ光の光子数、σ_Ｃをチェレンコフ光スペクトルの半値幅、ｔ_０をチェレンコフ光のピーク時刻として、以下の式（１）で与えられる。

一方、高エネルギー電子１２は、走行に伴い、シンチレータ１１の中で価電子１４を励起し、励起された価電子１４が基底状態に戻るとき、シンチレーション光１５が放出される。曲線１８は、価電子の励起密度ｎ_ｅｘ（ｔ）を、時刻ｔのガウス関数で近似したものである。価電子の励起密度ｎ_ｅｘ（ｔ）は、Ｎ_ｓをシンチレーション光の光子数、σ_Ｓをシンチレーション光スペクトルの半値幅、ｔ_１をシンチレーション光のピーク時刻として、以下の式（２）で与えられる。

曲線２１は、時刻２０において励起された価電子がもとの状態に遷移していく様子を示している。励起された価電子がもとの状態に遷移するのに伴ってシンチレーション光が発生する。曲線２２は、シンチレーション光の光子密度ｎ_ｓ（ｔ）を、時刻ｔの関数として示したものである。シンチレーション光の光子密度ｎ_ｓ（ｔ）は、緩和定数をτ_Sとして、価電子の励起密度ｎ_ｅｘ（ｔ）を用いて、以下の式（３）で与えられる。

ここで、価電子の励起密度ｎ_ｅｘ（ｔ）の関数形を式（２）で仮定した場合、式（３）に式（２）を代入すると、以下の式（４）が得られる。

すなわち、ＰＥＴ装置１００は、所定のキャリブレーションスキャンを行うことで、第１の検出器ｌ及び第２の検出器ｍの各々で取得されるエネルギースペクトルを取得し、当該エネルギースペクトルを基に曲線１７や曲線２２に示されるスペクトルを生成し、生成したスペクトルに基づいて、第１の光子数情報及び第２の光子数情報を取得してもよい。

なお、発光の分布関数ｐ_ｐｈ（ｔ）は、チェレンコフ光の光子密度ｎ_ｃ（ｔ）とシンチレーション光の光子密度ｎ_ｓ（ｔ）との和に比例するから、以下の式（５）が成り立つ。

ここで、発光の分布関数ｐ_ｐｈ（ｔ）は、以下の式（６）のように規格化されている。

なお、ステップＳ１００において、処理回路１５０は、取得機能１５０ａにより、光子数情報以外の情報を、第１の検出器ｌに係るパラメータセットとして、取得してもよい。一例として、処理回路１５０は、取得機能１５０ａにより、チェレンコフ光に対する光センサの時間分解能σ_Ｃを、第１の光子数情報に代えて、または第１の光子数情報に加えて、第１の検出器ｌにおいて取得するパラメータセットθ_ｌの一つとして、取得してもよい。また、別の例として、処理回路１５０は、取得機能１５０ａにより、シンチレーション光に対する光センサの時間分解能σ_Ｓを、第１の光子数情報に代えて、または第１の光子数情報に加えて、第１の検出器ｌにおいて取得するパラメータセットθ_ｌの一つとして、取得してもよい。また、別の例として、処理回路１５０は、取得機能１５０ａにより、過剰ノイズ発生確率を、第１の検出器ｌにおいて取得するパラメータセットθ_ｌの一つとして、取得してもよい。ここでいう過剰ノイズの発生確率とは、例えばオプティカルクロストークやアフターパルス等のＡＰＤ（Ａｖａｌａｎｃｈｅ ｐｈｏｔｏ ｄｉｏｄｅ）過剰ノイズの発生確率を言う。

ステップＳ２００において、処理回路１５０は、取得機能１５０ａにより、第１の検出器ｌとは異なる第２の検出器ｍで検出される第２の光子数情報を取得する。第２の光子数情報とは、例えば、第２の検出器ｍで検出される標準的な一イベントについて、平均的に検出されるシンチレーション光子数Ｎ_Ｓ;ｍである。また、別の例として、第２の光子数情報とは、例えば、第２の検出器ｍで検出される標準的な一イベントについて、平均的に検出されるチェレンコフ光子数Ｎ_Ｃ;ｍである。また、別の例として、第２の光子数情報とは、第２の検出器で一つのイベントあたり平均的に検出されるシンチレーション光子数Ｎ_Ｓ;ｍ及びシンチレーション光子数Ｎ_Ｓ;ｍである。

また、ＰＥＴ装置１００は、所定のキャリブレーションスキャンを行うことで、第２の光子数情報を特定してもよい。具体的には、ＰＥＴ装置１００は、所定のキャリブレーションスキャンを行うことで、第２の検出器ｍで取得されるエネルギースペクトルに基づいて、第２の光子数情報を取得してもよい。

また、処理回路１５０は、取得機能１５０ａにより、チェレンコフ光に対する光センサの時間分解能σ_Ｃ、シンチレーション光に対する光センサの時間分解能σ_Ｓ、過剰ノイズ発生確率等を、第２の検出器ｍにおいて取得するパラメータセットθ_ｍとして、取得してもよい。

ステップＳ３００において、処理回路１５０は、算出機能１５０ｂにより、ステップＳ１００で取得された第１の光子数情報に基づいて、第１の検出器ｌに対する第１の発光確率モデルp_ｔ；ｌ ^(doi)を算出する。なお、（ｄｏｉ）との記号は、シンチレータ内での発光位置の不定性を考慮する実施形態と、考慮しない実施形態との両方が可能である旨の趣旨の記号であり、処理回路１５０は、算出機能１５０ｂにより、シンチレータ内での発光位置の不定性を考慮しない実施形態においては、第１の発光確率モデルp_ｔ；ｌを算出することを意味し、シンチレータ内での発光位置の不定性を考慮する実施形態においては、後述する、第１の発光確率モデルp_ｔ；ｌ ^ｄｏｉを算出することを意味する。

図４にステップＳ３００の処理の一例がより詳細に説明されている。図４のフローチャートは、ステップＳ３００の処理の一例をより詳細に説明したフローチャートである。

ステップＳ３１０において、処理回路１５０は、算出機能１５０ｂにより、ステップＳ１００で取得された第１の検出器ｌにおける第１の光子数情報に基づいて、発光の分布関数モデルp_ｐｈ，ｌ (t)を生成する。発光の分布関数p_ｐｈ (t)は、時刻ｔに発光が起こる確率を示す分布関数である。発光の分布関数モデルp_ｐｈ,l (t)は、時刻に関して全時間で積分すると１になる。また、発光の分布関数モデルp_ｐｈ，ｌ (t)は、検出器ごとに定められる。

処理回路１５０は、算出機能１５０ｂにより、第１の光子数情報に基づいて、例えば式（１）、（４）、（５）、（６）等を用いて発光の分布関数モデルp_ｐｈ,l(t)を生成する。第１の光子数情報が、検出器ごとに取得されているので、処理回路１５０は、算出機能１５０ｂにより、発光の分布関数モデルp_ｐｈ,l (t)を、検出器ごとに算出する。図５に、発光の分布関数p_ｐｈ(t)の一例が示されている。図５において、曲線３０は、発光の分布関数p_ｐｈ(t)を示している。また、検出時刻を定義するためのフォトンの閾値（Ｎ’）を併記する。

続いて、ステップＳ３２０Ａにおいて、処理回路１５０は、算出機能１５０ｂにより、発光の分布関数モデルp_ｐｈ,l (t)に基づいて、第１の発光確率モデルｐ_ｔ,l(t)を算出する。

ここで、第１の発光確率モデルｐ_ｔ,l(≧Ｎ’，t)は、検出事象に含まれる複数の光子の内の最初のＮ’個の光子を検出する確率密度に関するモデルである。換言すると、第１の発光確率モデルｐ_ｔ,l(≧Ｎ’，t)は、検出事象に含まれる複数の光子の内の所定の数の光子を検出する確率密度に関するモデルである。図６に、第１の発光確率モデルｐ_ｔ,l(≧Ｎ’，t)が示されている。図６のグラフ３１は、第１の発光確率モデルｐ_ｔ,l(≧Ｎ’，t)の概形を示している。なお、第１の発光確率モデルｐ_ｔ,l(≧Ｎ’，t)は、検出器ごとに定められる。以下、必要に応じて、検出器の添え字については適宜省略する。

検出事象に含まれる複数の光子の内の最初のＮ’個の光子を検出する確率密度である第１の発光確率モデルｐ_t（≧Ｎ’，ｔ）は、以下のように算出することができる。

はじめに、時刻ｔまでにＮ’個未満のフォトンしか検出しない確率Ｐ’（≧Ｎ’，ｔ）に関する方程式を導出する。これは、時刻ｔまでに高々Ｎ’－１個のフォトンしか検出していない確率であり、時刻ｔまでにｋ個のフォトンを検出している確率Ｐｄ（ｋ，ｔ）を用いて以下の式（７）で表される。

時刻ｔまでにｋ個のフォトンを検出している確率は、フォトンの検出事象を２項分布で評価することで、以下の式（８）のように微分形式で表現できる。

この微分方程式を順次ｋ＝０、１，２，…に対して解くと、帰納的に、以下の式（９）
が成り立つ。

式（７）より、時刻ｔまでにＮ’個未満のフォトンしか検出しない確率Ｐ’（≧Ｎ’，ｔ）は以下の式（１０）のように書くことが出来る。

これを式（８）に代入することで、検出事象に含まれる複数の光子の内の最初のＮ’個の光子を検出する確率密度である第１の発光確率モデルｐ_t（≧Ｎ’，ｔ）は、以下の式（１１）のように表すことが出来る。

ここで、関数F（ｔ）は発光の分布関数ｐ_ｐｈ（ｔ）の累積分布関数であり、以下の式（１２）のように書くことが出来る。

以下の実施例では簡単のため、初めてフォトンを検出した時刻で評価した確率分布、すなわち、Ｎ’＝０とした場合を例にとる。すなわち、初の光子を検出する確率密度である第１の発光確率モデルをp_ｔ (t)とし、以下の式（１３）のように陽に表す。

すなわち、処理回路１５０は、算出機能１５０ｂにより、式（１２）及び式（１３）を用いることで、発光の分布関数ｐ_ph(t)に基づいて、第１の発光確率モデルp_ｔ(t)を計算することができる。このように、処理回路１５０は、算出機能１５０ｂにより、式（１２）及び式（１３）に第１の検出器における発光の分布関数モデルp_ｐｈ (t)を代入し、第１の光子数情報を式（１３）の右辺のＮに代入して式（１３）の右辺を評価することで、第１の光子数情報及び第１の検出器における発光の分布関数モデルp_ｐｈ (t)に基づいて、第１の発光確率モデルｐ_ｔ(t)を算出する。

なお、上述の発光の分布関数ｐ_ph(t)及び第１の第１の発光確率モデルp_t(t)は、各検出器ごとに定めることができる。従って、検出器の添え字を明示的に記載すると、式（１３）は、以下の式（１４）のように書くこともできる。

ここで、p_ｔ；ｌ(t)は検出器ｌにおける第１の発光確率モデル、ｐ_ｐｈ；ｌ(t)は検出器ｌにおける発光の分布関数、Ｎ_ｓ；ｌは検出器ｌにおけるシンチレーション光の光子数、Ｎ_ｃ；ｌは検出器ｌにおけるチェレンコフ光の光子数を示す。この式よりわかるように、第１の発光確率モデルp_t,l(t)は、検出事象に対応する光子数の期待値Ｎ_ｓ；ｌ、Ｎ_ｃ；ｌに基づく。

なお、第１の発光確率モデルを算出の際、処理回路１５０は、算出機能１５０ｂにより、シンチレータ内での発光位置の不定性を考慮した計算を行うこともできる。例えば、図７に示されるように、ガンマ線４０がシンチレータ１１の深さｘにある点４２で発光し、矢印４３で示されるように、シンチレーション光が物質中の光速ｃ／ｎで屈折率ｎの媒体中を伝播してシンチレータ１１の端４４に到達した場合、シンチレータ内の発光位置によって、ガンマ線４０がシンチレータに入射してからシンチレータ１１の端に到達するまでの時刻が異なるので、信号を検出する時刻が異なる。

すなわち、ガンマ線４０がシンチレータ１１に入射してから深さｘの地点まで到達するのに要する時刻ｔ_１は、以下のように、式（１５）で与えられる。

また、点４２で発光した光が深さｘからシンチレータ１１の端まで到達するのに要する時刻ｔ_２は、以下のように、式（１６）で与えられる。

従って、ガンマ線４０がシンチレータ１１に入射してからシンチレーション光がシンチレータ１１の端まで到達するのに要する時刻ｔ’は、以下のように、式（１７）で与えられる。

式（１７）より、シンチレータ内の発光位置ｘによって、信号を検出する時刻が異なることがわかる。

このように、シンチレータ内の発光位置ｘを考慮した取り扱いを考える場合、以下の式（１８）が成り立つ。

ここで、p_ｔ、ｌ ^ｄｏｉ(t)は、シンチレータ内での発光位置の不定性を考慮した第１の発光確率モデルであり、ｐ_ｐｈ、ｌ ^ｄｏｉ(t)は、シンチレータの位置の発光位置の不定性を考慮した発光の分布関数モデルである。

なお、式（１８）は、処理回路１５０が取得機能１５０ａによりステップＳ１００において取得するパラメータθ_ｌを関数の引数として明示的に表現すると、以下の式（１９）及び式（２０）のようにも書くことができる。

ここで、θ_ｌは、検出器ｌにおいて取得されるパラメータセットを示しており、例えば、シンチレーション光子数Ｎ_ｓ；ｌ、チェレンコフ光子数Ｎ_ｃ；ｌ、チェレンコフ光に対する光センサの時間分解能σ_Ｃ、シンチレーション光に対する光センサの時間分解能σ_Ｓ等のパラメータの組を表す。

式（１８）に戻り、シンチレータの位置の発光位置の不定性を考慮した場合のステップＳ３００の処理について説明する。図８は、シンチレータ内での発光位置の不定性を考慮した場合の、発光確率モデルの算出の流れを示したフローチャートである。

図４で示した場合と同様に、ステップＳ３１０において、処理回路１５０は、算出機能１５０ｂにより、第１の検出器ｌにおける第１の光子数情報等のパラメータθ_ｌに基づいて、発光の分布関数モデルｐ_ph,l(t)を生成する。

続いて、ステップＳ３１１において、処理回路１５０は、算出機能１５０ｂにより、シンチレータ内での位置ｘごとに、発光確率密度ｐ_ｄｏｉ（ｘ）を、減弱係数μの値に基づいて算出する。位置ｘにおける発光確率密度ｐ_ｄｏｉ（ｘ）は、以下の式（２１）で表される。

なお、発光確率密度ｐ_ｄｏｉ（ｘ）は、以下の式（２２）に示されるように、規格化されている。

なお、図７において、グラフ４５には、減弱率μexp(-μx)が示されている。

図８に戻り、ステップＳ３１２において、処理回路１５０は、算出機能１５０ｂにより、ステップＳ３１０で生成された分布関数モデルｐ_ｐｈ、ｌ（ｔ）と発光確率密度ｐ_ｄｏｉ（ｘ）とに基づいて、発光位置不確定性を考慮した発光の分布関数モデルｐ_ｐｈ；ｌ ^ｄｏｉ（ｔ）を算出する。

具体的には、処理回路１５０は、算出機能１５０ｂにより、発光位置不確定性を考慮した発光の分布関数モデルｐ_ｐｈ；ｌ ^ｄｏｉ（ｔ）を、ステップＳ３１０で生成された分布関数モデルｐ_ｐｈ、ｌ（ｔ）と、ステップＳ３１１で算出された発光確率密度ｐ_ｄｏｉ（ｘ）とに基づいて、次の式（２３）に基づいて算出する。

ここで、Ｌはシンチレータの長さであり、ｘは発光位置であり、ｔ’は式（１７）において説明した通り、ガンマ線がシンチレータに入射してからシンチレーション光がシンチレータの端まで到達するまでの遅延時間である。

続いて、ステップＳ３２０Ｂにおいて、処理回路１５０は、算出機能１５０ｂにより、発光位置不確定性を考慮した発光の分布関数モデルｐ_ｐｈ；ｌ ^ｄｏｉ（ｔ）に基づいて、発光位置不確定性を考慮した発光確率モデルｐ_ｔ；ｌ ^ｄｏｉ（θ_l;t）を算出する。具体的には、式（２４）及び式（２５）に示されるように、処理回路１５０は、算出機能１５０ｂにより、式（１３）と同様に、発光位置不確定性を考慮した発光の分布関数モデルｐ_ｐｈ ^ｄｏｉ（ｔ）に基づいて、発光位置不確定性を考慮した発光確率モデルｐ^ｄｏｉ（t）を算出する。

なお、式（２３）において、式（１８）を基に変数変換すると、以下の式（２６）が得られるので、処理回路１５０は、算出機能１５０ｂにより、式（２４）の右辺の評価において、式（２６）を用いて算出を行ってもよい。

このようにして、処理回路１５０は、算出機能１５０ｂにより、シンチレータ内での発光位置の不定性を考慮して、第１の発光確率モデルおよび第２の発光確率モデルを算出する。

図２に戻り、ステップＳ４００において、処理回路１５０は、算出機能１５０ｂにより、ステップＳ３００と同様の処理を行うことにより、ステップＳ２００で取得された第２の光子数情報に基づいて、第２の検出器ｍに対する第２の発光確率モデルp_ｔ；ｍ ^(doi)を算出する。

具体的には、処理回路１５０は、算出機能１５０ｂにより、ステップＳ２００で取得された第２の検出器ｍにおける第２の光子数情報に基づいて、例えば式（１）、（４）、（５）、（６）等に相当する式を用いて、発光の分布関数モデルp_ｐｈ；ｍ (t)を生成する。

続いて、処理回路１５０は、算出機能１５０ｂにより、第２の光子数情報及び第２の検出器ｍにおける発光の分布関数モデルp_ｐｈ；ｍ(t)に基づいて、第２の発光確率モデルｐ_ｔ；ｍ(t)を、例えば式（１３）及び式（１４）等に相当する式を用いて算出する。

また、第２の発光確率モデルp_t,ｍ(t)は、検出事象に対応する光子数の期待値Ｎ_ｓ；ｍ、Ｎ_ｃ；ｍに基づく。

なお、第２の発光確率モデルを算出の際、処理回路１５０は、算出機能１５０ｂにより、シンチレータ内での発光位置の不定性を考慮した計算を、式（１７）及び式（２１）～（２６）と同様の式を用いて行うこともできる。この場合、処理回路１５０は、算出機能１５０ｂにより、シンチレータ内での位置ｘごとに、発光確率密度ｐ_ｄｏｉ（ｘ）を、減弱係数μの値に基づいて算出する。続いて、処理回路１５０は、算出機能１５０ｂにより、分布関数モデルｐ_ｐｈ；ｍ（ｔ）と発光確率密度ｐ_ｄｏｉ（ｘ）とに基づいて、発光位置不確定性を考慮した発光の分布関数モデルｐ_ｐｈ；ｍ ^ｄｏｉ（ｔ）を算出する。続いて、処理回路１５０は、算出機能１５０ｂにより、発光位置不確定性を考慮した発光の分布関数モデルｐ_ｐｈ；ｍ ^ｄｏｉ（ｔ）に基づいて、発光位置不確定性を考慮した発光確率モデルｐ_ｔ；ｍ ^ｄｏｉ（t）を算出する。

このようにして、処理回路１５０は、算出機能１５０ｂにより、シンチレータ内での発光位置の不定性を考慮して、第２の発光確率モデルｐ_ｔ,ｍ(t)を算出する。

続いて、図２に戻り、ステップＳ５００において、処理回路１５０は、特定機能１５０ｃにより、第１の発光確率モデルｐ_ｔ；ｌ ^ｄｏｉ（t_１）と第２の発光確率モデルｐ_ｔ；ｍ ^ｄｏｉ（t_２）とに基づいて、第１の検出器ｌと第２の検出器ｍとによって規定されるＬＯＲの検出時間差に関する確率分布モデルｐ_{ＴＯＦ;l,m} ^doi(θ_l, θ_ｍ,t)を特定する。

図９において、曲線５０は、第１の検出器ｌにおける第１の発光確率モデルｐ_ｔ；ｌ ^ｄｏｉ（θ_l,t_l）を、遅延時間ｔ_ｌの関数としてプロットしたものを示している。また、曲線５１は、第２の検出器ｍにおける第２の発光確率モデルｐ_ｔ；ｍ ^ｄｏｉ（θ_ｍ,t_ｍ）を、遅延時間ｔ_ｍの関数としてプロットしたものを示している。処理回路１５０は、特定機能１５０ｃにより、検出器ｌにおける検出時刻ｔ_ｌにおける事象と、検出器ｍにおける検出時刻ｔ_ｍにおける積事象を、検出時間差の関数として表現した確率分布モデルｐ_{ＴＯＦ;l,m} ^doi(θ_l, θ_ｍ,t)を、特定する。曲線５２は、このような確率分布モデルｐ_{ＴＯＦ;l,m} ^doi(θ_l, θ_ｍ,t)の例を示している。

ここで、時間幅ｄｔの間における、検出器ｌにおける検出時刻ｔ_ｌにおける事象と、検出器ｍにおける検出時刻ｔ_ｍにおける積事象の数は、以下の式（２７）で与えられる。

また、検出器ｌにおける検出時刻ｔ_ｌと、検出器ｍにおける検出時刻ｔ_ｍとの検出時間の差は、以下の式（２８）で与えられる。

従って、ｐ_{ＴＯＦ;l,m} ^doi(θ_l, θ_ｍ,t)は、以下の式（２９）で与えられる。

以上のように、処理回路１５０は、算出機能１５０ｂにより、式（２９）に基づいて、第１の検出器ｌと第２の検出器ｍとによって規定されるＬＯＲの検出時間差に関する確率分布モデルｐ_{ＴＯＦ;l,m} ^doi(θ_l, θ_ｍ,t)を算出する。処理回路１５０は、この確率分布モデルを、ＰＥＴ画像再構成におけるＴＯＦカーネルとして利用することで、検出器ごとのＴＯＦ分解能の差をＰＥＴ画像再構成に反映させることができ、再構成画像の画質を向上させることができる。

例えば、図１０に示されているように、処理回路１５０は、算出機能１５０ｂにより、ＬＯＲごとに、ＬＯＲの検出時間差に関する確率分布モデルｐ_{ＴＯＦ;l,m} ^doi(θ_l, θ_ｍ,t)を算出する。例えば、処理回路１５０は、算出機能１５０ｂにより、検出器６０と検出器６１との間のＬＯＲ６３に関する確率分布モデルｐ_{ＴＯＦ;l,m'} ^doi(θ_l, θ_ｍ',t)を算出し、また、検出器６０と検出器６２との間のＬＯＲ６４に関する確率分布モデルｐ_{ＴＯＦ;l,m} ^doi(θ_l, θ_ｍ,t)を算出する。このような構成をとることで、検出器ごとに個別の確率分布モデルを生成することができ、例えば検出器ごとに性能差にばらつきがある場合などでも、再構成画像の画質を向上させることができる。

図２に戻り、ステップＳ６００において、処理回路１５０は、再構成機能１５０ｄにより、ステップＳ５００で特定された確率分布モデルｐ_{ＴＯＦ;l,m} ^doi(θ_l, θ_ｍ,t)をＴｏＦカーネルとして利用することで、ＬＯＲに対応するデータの再構成を行い、再構成画像を生成する。すなわち、処理回路１５０は、再構成機能１５０ｄにより、検出時間差に関する確率分布モデルをＴｏＦ（Ｔｉｍｅ Ｏｆ Ｆｌｉｇｈｔ）カーネルとして利用することで、ＬＯＲに対応するデータの再構成を行う。

図１１～図１２に、第１の実施形態に係る手法の検証結果が示されている。この検証では、１Ｄ－ＭＬ-ＥＭ（Ｍａｘｉｍｕｍ Ｌｉｋｅｌｉｈｏｏｄ Ｅｘｐｅｃｔａｔｉｏｎ Ｍａｘｉｍｉｚａｔｉｏn)法を用いて、第１の実施形態に係る手法により生成されたＴｏＦカーネルの性能の評価を行った。

位置ｊの放射体の同時計測が位置ｉで起こる場合を考え、空間座標を同時計数時刻で定義し、ＴｏＦカーネルをｃ_ｊｉとし、放射線物質の分布をλ（ｊ）とすると、コインシデンススペクトルｙ（ｉ）は、以下の式（３０）で与えられる。

ここで、ＭＬ－ＥＭ法を用いて、放射性物質の分布を求めると、ｋ＋１回目の更新における放射性物質の分布λ^ｋ＋１は、以下の式（３１）で与えられる。

ここで、ＴｏＦカーネルが妥当なものであれば、ＭＬ－ＥＭ法にて再構成される放射性物質の分布が、狭い領域にピークを持つと考え、種々のＴｏＦカーネルを用いて、放射性物質の分布の再構成を行った。

具体的には、検討するＴｏＦカーネルの種類として、４種類のＴｏＦカーネルを検討した。すなわち、チェレンコフ光のガウス分布とシンチレーション光のガウス分布の和を用いた複数のガウスシアンカーネル、シンチレーション光のみを仮定した単一ガウシアンカーネル、チェレンコフ光のみを仮定した単一ガウシアンカーネル、及び第１の実施形態に係る、発光モデルのモデルベースでチェレンコフ光とシンチレーション光の両方を考慮したモデルベースカーネルの４つのカーネルを検討し、図１１に示されているように、実測スペクトル７０に対してフィッティングを行った。曲線７３、曲線７２、曲線７１及び曲線７４は、それぞれ、複数のガウシアンカーネル、シンチレーション光のみを仮定した単一ガウシアンカーネル、チェレンコフ光のみを仮定した単一ガウシアンカーネル、モデルベースカーネルを用いてフィティングを行った結果を示している。

これら４種類のＴｏＦカーネルを用いて、ＭＬ－ＥＭ法により再構成を行い、放射性物質の分布の再構成を行った結果が、図１２に示されている。図１２において、曲線８３、曲線８２、曲線８１、曲線８４が、それぞれ、複数のガウシアンカーネル、シンチレーション光のみを仮定した単一ガウシアンカーネル、チェレンコフ光のみを仮定した単一ガウシアンカーネル、モデルベースカーネルを用いて再構成を行った場合の計算結果である。モデルベースカーネルを用いて再構成を行った場合、ＭＬ－ＥＭ法にて再構成される放射性物質のピークの幅が狭くなっており、第１の実施形態に係る手法を用いて生成したＴｏＦカーネルの妥当性が確認された。

以上のように、第１の実施形態に係る核医学診断装置においては、処理回路１５０は、シンチレータ内の発光モデルに基づいて、検出時間差に関する確率分布モデルを、検出器ごとに算出した。これにより、ＬＯＲごとに個別のＴｏＦカーネルを設定して画像再構成をすることができるので、例えば検出器ごとに性能にばらつきがある場合等でも、再構成後のＰＥＴ画像の画質を向上させることができる。

実施形態は、上述の例に限られない。上述の場合においては、検出器３がモノリシック検出器でない通常の検出器である場合について説明したが、実施形態はこれらに限られない。実施形態は、検出器３が、単結晶として構成され、１つのシンチレータ結晶が区分けされずに一体化したまま使用されるモノリシックシンチレータである場合であってもよい。当該シンチレータは、例えば、ＬＹＳＯ（Lutetium Yttrium Oxyorthosilicate）、ＬＳＯ（Lutetium Oxyorthosilicate）、ＬＧＳＯ（Lutetium Gadolinium Oxyorthosilicate）等の、シンチレータ結晶によって形成される。

検出器３がモノリシック検出器である場合、検出器３は一つのシンチレータ結晶が区分けされずに一体化したまま使用されるものであるが、検出器３が複数の仮想的な検出器からなると考えて、それぞれの仮想的な検出器の間のＴｏＦカーネルを生成することにより、同様の処理を行うことができる。例えば、検出器３が複数の検出素子アレイからなる場合、それぞれの検出素子アレイを１つの検出器と考えて、図２のフローチャートで示された処理と同様の処理を行うことができる。一例として、図２のステップＳ１００からステップＳ６００や、式（２０）、（２９）等における自然数ｌ及びｍを、検出器に付された番号の代わりに、検出素子アレイに付された番号であると読み替えて同様な処理を行うことで、処理回路１５０は、特定機能当該検出素子アレイの組み合わせに応じた検出時間差に関する確率分布モデルを特定することができる。

すなわち、検出器３がモノリシック検出器である場合、実施形態における核医学診断装置は、少なくとも第１の検出素子アレイｌと第２の検出素子アレイｍとが設けられたモノリシック検出器からなる検出器３と、処理回路１５０とを備える。処理回路１５０は、取得機能１５０ａにより、第１の検出素子アレイｌにて検出される第１の光子数情報と、第２の検出素子アレイｍにて検出される第２の光子数情報とを取得する。処理回路１５０は、算出機能１５０ｂにより、第１の光子数情報に基づいて第１の検出素子アレイｌに対応する第１の発光確率モデルｐ_ｔ；ｌ ^ｄｏｉ（θ_l,t_l）と、第２の光子数情報に基づいて第２の検出素子アレイｍに対応する第２の発光確率モデルｐ_ｔ；ｍ ^ｄｏｉ（θ_ｍ,t_ｍ）とを算出する。処理回路１５０は、特定機能１５０ｃにより、第１の発光確率モデルｐ_ｔ；ｌ ^ｄｏｉ（θ_l,t_l）と第２の発光確率モデルｐ_ｔ；ｍ ^ｄｏｉ（θ_ｍ,t_ｍ）とに基づいて、第１の検出素子アレイｌと第２の検出素子アレイｍとの組み合わせに応じた検出時間差に関する確率分布モデルｐ_{ＴＯＦ;l,m} ^doi(θ_l, θ_ｍ,t)を特定する。

また、実施形態に係るデータ処理方法は、少なくとも第１の検出素子アレイｌと第２の検出素子アレイｍとが設けられたモノリシック検出器の第１の検出素子アレイｌにて検出される第１の光子数情報と、第２の検出素子アレイｍにて検出される第２の光子数情報とを取得し、第１の光子数情報に基づいて第１の検出素子アレイｌに対応する第１の発光確率モデルｐ_ｔ；ｌ ^ｄｏｉ（θ_l,t_l）と、第２の光子数情報に基づいて第２の検出素子アレイｍに対応する第２の発光確率モデルｐ_ｔ；ｍ ^ｄｏｉ（θ_ｍ,t_ｍ）とを算出し、第１の発光確率モデルｐ_ｔ；ｌ ^ｄｏｉ（θ_l,t_l）と第２の発光確率モデルｐ_ｔ；ｍ ^ｄｏｉ（θ_ｍ,t_ｍ）とに基づいて、第１の検出素子アレイｌと第２の検出素子アレイｍとの組み合わせに応じた検出時間差に関する確率分布モデルｐ_{ＴＯＦ;l,m} ^doi(θ_l, θ_ｍ,t)を特定する。

また、実施形態に係るプログラムは、少なくとも第１の検出素子アレイｌと第２の検出素子アレイｍとが設けられたモノリシック検出器の第１の検出素子アレイｌにて検出される第１の光子数情報と、第２の検出素子アレイｍにて検出される第２の光子数情報とを取得し、第１の光子数情報に基づいて第１の検出素子アレイｌに対応する第１の発光確率モデルｐ_ｔ；ｌ ^ｄｏｉ（θ_l,t_l）と、第２の光子数情報に基づいて第２の検出素子アレイｍに対応する第２の発光確率モデルｐ_ｔ；ｍ ^ｄｏｉ（θ_ｍ,t_ｍ）とを算出し、第１の発光確率モデルｐ_ｔ；ｌ ^ｄｏｉ（θ_l,t_l）と第２の発光確率モデルｐ_ｔ；ｍ ^ｄｏｉ（θ_ｍ,t_ｍ）とに基づいて、第１の検出素子アレイｌと第２の検出素子アレイｍとの組み合わせに応じた検出時間差に関する確率分布モデルｐ_{ＴＯＦ;l,m} ^doi(θ_l, θ_ｍ,t)を特定する処理をコンピュータに実行させる。

（第２の実施形態）
はじめに、第２の実施形態の背景について説明する。ＰＥＴ画像の再構成を行う際において、ＴｏＦ（Ｔｉｍｅ Ｏｆ Ｆｌｉｇｈｔ）カーネルを用いて画像再構成が行われる。ＴｏＦカーネルとは、２つの検出器で検出された信号の検出時間差の関数として表された、検出イベントの確率密度関数である。ここで、画像再構成に用いられるＴｏＦカーネルの時間幅が、検出器のＴｏＦ性能に整合するものでない場合、再構成後の画像の画質が低下する場合がある。例えば、検出器のＴｏＦ性能と比較して、ＴｏＦカーネルの時間幅が小さい場合には、画像再構成後の放射性物質の位置の推定位置が、本来の分布よりも局所的に集積してしまう。一方、検出器のＴｏＦ性能と比較して、ＴｏＦカーネルの時間幅が大きい場合には、画像再構成後の放射線物質の位置の推定位置が、本来の分布より周辺部分に分布してしまう。従って、検出器のＴｏＦ性能と、画像再構成に用いるＴｏＦカーネルとが、できるだけ一致することが望ましい。

しかしながら、例えば、イベント数が少ない場合等、観測揺らぎが大きくなり、ＴｏＦカーネルの算出に誤差が生じる場合がある。そこで、第２の実施形態では、発光モデルに基づいて、同時計数イベント毎に確率分布モデルを生成したのち、それらの確率分布を加算することによりＴｏＦカーネルの算出を行う。このことで、イベントごとの観測値の平均からの揺らぎを考慮したＴｏＦカーネルを生成して画像再構成をすることができるので、例えば少ない同時計数イベントで観測揺らぎが大きく見えてしまう場合等でも、再構成後のＰＥＴ画像の画質を向上させることができる。

かかる背景に基づいて、第２の実施形態に係る核医学診断装置は、処理回路１５０を備える。処理回路１５０は、取得機能１５０ａにより、第１の検出器で検出される第１のイベントに関する第１の光子数情報と、第１の検出器とは異なる第２の検出器で検出される第２のイベントに関する第２の光子数情報とを取得する。処理回路１５０は、特定機能１５０ｃにより、第１の光子数情報と第２の光子数情報とに基づいて、第１のイベントと第２のイベントとに基づいて規定されるコインシデンスに関する確率分布モデルを特定する。

また、第２の実施形態に係るデータ処理方法は、第１の検出器で検出される第１のイベントに関する第１の光子数情報と、第１の検出器とは異なる第２の検出器で検出される第２のイベントに関する第２の光子数情報とを計測し、第１の光子数情報と第２の光子数情報とに基づいて、第１のイベントと第２のイベントとに基づいて規定されるコインシデンスに関する確率分布モデルを特定する。

また、第２の実施形態に係るプログラムは、第１の検出器で検出される第１のイベントに関する第１の光子数情報と、前記第１の検出器とは異なる第２の検出器で検出される第２のイベントに関する第２の光子数情報とを計測し、第１の光子数情報と第２の光子数情報とに基づいて、第１のイベントと第２のイベントとに基づいて規定されるコインシデンスに関する確率分布モデルを特定する処理をコンピューターに実行させる。

以下、図１３～図１７を用いて、第２の実施形態に係る核医学診断装置が行う処理について、詳しく説明する。

図１３に、第２の実施形態に係る核医学診断装置としてのＰＥＴ装置１００が行う処理の流れが示されている。以下のフローチャートでは、簡単のため、第１の検出器ｌと第２の検出器ｍに係る確率分布モデル（ＴｏＦカーネル）を生成する場合についてのみ記載している。しかしながら、典型的には、画像再構成においては、すべてのｌとｍの組について処理回路１５０は確率分布モデルを生成した上でステップＳ７０の画像再構成処理を行う。この場合、処理回路１５０は、すべてのｌとｍの組についてステップＳ１０からステップＳ６０の処理を繰り返した上で、ステップＳ７０の処理を行う。

ステップＳ１０において、処理回路１５０は、処理を行う検出器の組を示すｌとｍの値を選択する。すなわち、処理回路１５０は、第１のイベントを検出する検出器である第１の検出器ｌ（ｌ番目の検出器）と、第２のイベントを検出する検出器である第２の検出器ｍ（ｍ番目の検出器）とを特徴づけるｌとｍの値を選択する。なお、ここでいう第１のイベントと第２のイベントとで、一つの同時計数イベントが構成される。別の言い方をすれば、同時計数イベントの各々は、第１のイベントと第２のイベントの２つのイベントからなる。

続いて、ステップＳ２０において、処理回路１５０は、同時計数イベントの番号を示すｉの値に１をセットする。ここで、ｉ=1とは、処理回路１５０が、１番目の同時計数イベントの処理を行うことを意味する。

続いて、ステップＳ３０において、処理回路１５０は、特定機能１５０ｃにより、ｌ番目の検出器で検出されたイベントと、ｍ番目の検出器で検出されたイベントとの同時計数イベントのうちのｉ番目の同時計数イベントに対して、コインシデンスに関する確率分布モデルｐ_{ｃｔ；ｌ、ｍ} ^{ｉ、（ｄｏｉ）}（θ_ｌ ^ｉ,θ_ｍ ^ｉ、ｔ）を算出する。ここで、θ_ｌ ^ｉ、θ_ｍ ^ｉは、ｉ番目の同時計数イベントにおいて、それぞれｌ番目の検出器及びｍ番目の検出器について取得されるパラメータを示す。また、（ｄｏｉ）の記号は、シンチレータ内での発光位置の不定性を考慮してもよいし、しなくてもよいことを意味し、シンチレータ内での発光位置の不定性を考慮する実施形態においては、処理回路１５０は、特定機能１５０ｃにより、シンチレータ内での発光位置の不定性を考慮したコインシデンスに関する確率分布モデルｐ_{ｃｔ；ｌ、ｍ} ^{ｉ、ｄｏｉ}（θ_ｌ ^ｉ,θ_ｍ ^ｉ、ｔ）を算出し、シンチレータ内での発光位置の不定性を考慮しない実施形態においては、シンチレータ内での発光位置の不定性を考慮しないコインシデンスに関する確率分布モデルｐ_{ｃｔ；ｌ、ｍ} ^ｉ（θ_ｌ ^ｉ,θ_ｍ ^ｉ、ｔ）を算出することを意味する。

図１４に、図１３のステップＳ３０の処理について、より詳細な処理の流れが示されている。すなわち、図１４のステップＳ１００～ステップＳ５００は、図１３のステップＳ３０をより詳細に説明したものとなる。

はじめに、ステップＳ１００において、処理回路１５０は、取得機能１５０ａにより、ｉ番目の同時計数イベントについて、第１の検出器ｌで検出される第１のイベントに関する第１の光子数情報を取得する。

ここで、第１の光子数情報とは、例えば、当該イベントについてのシンチレーション光子数Ｎ^ｉ _Ｓ;lである。また、別の例として、第１の光子数情報とは、例えば、当該イベントについてのチェレンコフ光子数Ｎ^ｉ _Ｃ;lである。また、別の例として、第１の光子数情報とは、当該イベントについてのシンチレーション光子数Ｎ^ｉ _Ｓ;l及びシンチレーション光子数Ｎ^ｉ _Ｓ;lである。すなわち、処理回路１５０は、取得機能１５０ａにより、シンチレーション光子数Ｎ^ｉ _Ｓ;l、チェレンコフ光子数Ｎ^ｉ _Ｃ;l、またはシンチレーション光子数Ｎ^ｉ _Ｓ;l及びチェレンコフ光子数Ｎ^ｉ _Ｃ;lを、第１の光子数情報として取得する。ここで、処理回路１５０が取得機能１５０ａにより取得する第１の光子数情報は、イベントごとに取得されるものである。

ここで、ＰＥＴ装置１００は、第１の検出器ｌで取得されるエネルギースペクトルに基づいて、第１の光子数情報を取得してもよい。

ここで、第１の実施形態においてすでに述べたように、チェレンコフ光の光子密度ｎ_ｃ（ｔ）は、Ｎ_ｃをチェレンコフ光の光子数、σ_Ｃをチェレンコフ光スペクトルの半値幅、ｔ_０をチェレンコフ光のピーク時刻として、式（１）で与えられる。また、価電子の励起密度ｎ_ｅｘ（ｔ）は、Ｎ_ｓをシンチレーション光の光子数、σ_Ｓをシンチレーション光スペクトルの半値幅、ｔ_１をシンチレーション光のピーク時刻として、式（２）で与えられる。また、シンチレーション光の光子密度ｎ_ｓ（ｔ）は、緩和定数をτ_Sとして、価電子の励起密度ｎ_ｅｘ（ｔ）を用いて、式（３）で与えられる。従って、価電子の励起密度ｎ_ｅｘ（ｔ）の関数形を式（２）で仮定した場合、式（３）に式（２）を代入すると、式（４）が得られる。

すなわち、ＰＥＴ装置１００は、第１の検出器ｌ及び第２の検出器ｍの各々で取得されるエネルギースペクトルを取得し、当該エネルギースペクトルを基に曲線１７や曲線２２に示されるスペクトルを生成し、生成したスペクトルに基づいて、第１の光子数情報及び第２の光子数情報を取得してもよい。

なお、発光の分布関数ｐ_ｐｈ（ｔ）は、チェレンコフ光の光子密度ｎ_ｃ（ｔ）とシンチレーション光の光子密度ｎ_ｓ（ｔ）との和に比例するから、式（５）が成り立つ。

ここで、発光の分布関数ｐ_ｐｈ（ｔ）は、式（６）のように規格化されている。

なお、ステップＳ１００において、処理回路１５０は、取得機能１５０ａにより、光子数情報以外の情報を、第１の検出器ｌに係るパラメータセットとして、取得してもよい。一例として、処理回路１５０は、取得機能１５０ａにより、チェレンコフ光に対する光センサの時間分解能σ_Ｃを、第１の光子数情報に代えて、または第１の光子数情報に加えて、第１の検出器ｌにおいて取得するパラメータセットθ_ｌの一つとして、取得してもよい。また、別の例として、処理回路１５０は、取得機能１５０ａにより、シンチレーション光に対する光センサの時間分解能σ_Ｓを、第１の光子数情報に代えて、または第１の光子数情報に加えて、第１の検出器ｌにおいて取得するパラメータセットθ_ｌの一つとして、取得してもよい。また、別の例として、処理回路１５０は、取得機能１５０ａにより、過剰ノイズ発生確率を、第１の検出器ｌにおいて取得するパラメータセットθ_ｌの一つとして、取得してもよい。ここでいう過剰ノイズの発生確率とは、例えばオプティカルクロストークやアフターパルス等のＡＰＤ（Ａｖａｌａｎｃｈｅ ｐｈｏｔｏ ｄｉｏｄｅ）過剰ノイズの発生確率を言う。

ステップＳ２００において、処理回路１５０は、取得機能１５０ａにより、ｉ番目の同時計数イベントについて、第１の検出器ｌとは異なる第２の検出器ｍで検出される第２のイベントに関する第２の光子数情報を取得する。

ここで、第２の光子数情報とは、例えば、当該同時計数イベントについてのシンチレーション光子数Ｎ^ｉ _Ｓ;ｍである。また、別の例として、第２の光子数情報とは、例えば、当該同時計数イベントについてのチェレンコフ光子数Ｎ^ｉ _Ｃ;ｍである。また、別の例として、第２の光子数情報とは、当該同時計数イベントについてのシンチレーション光子数Ｎ^ｉ _Ｓ;ｍ及びシンチレーション光子数Ｎ^ｉ _Ｓ;ｍである。すなわち、処理回路１５０は、取得機能１５０ａにより、シンチレーション光子数Ｎ^ｉ _Ｓ;ｍ、チェレンコフ光子数Ｎ^ｉ _Ｃ;ｍ、またはシンチレーション光子数Ｎ^ｉ _Ｓ;ｍ及びチェレンコフ光子数Ｎ^ｉ _Ｃ;ｍを、第２の光子数情報として取得する。ここで、処理回路１５０が取得機能１５０ａにより取得する第２の光子数情報は、イベントごとに取得されるものである。

ＰＥＴ装置１００は、第２の検出器ｍで取得されるエネルギースペクトルに基づいて、第２の光子数情報を取得してもよい。

また、処理回路１５０は、取得機能１５０ａにより、チェレンコフ光に対する光センサの時間分解能σ_Ｃ、シンチレーション光に対する光センサの時間分解能σ_Ｓ、過剰ノイズ発生確率等を、第２の検出器ｍにおいて取得するパラメータセットθ_ｍとして、取得してもよい。

ステップＳ３００において、処理回路１５０は、算出機能１５０ｂにより、ステップＳ１００で取得された第１の光子数情報に基づいて、ｉ番目の同時計数イベントにおける、第１の検出器ｌに対する第１の発光確率モデルp_ｔ；ｌ ^{ｉ,(ｄｏｉ)}を算出する。ここで、処理回路１５０は、算出機能１５０ｂにより、第１の検出器ｌに対する第１の発光確率モデルp_ｔ；ｌ ^{ｉ,(ｄｏｉ)}を、同時計数イベントごとに定める。

図１５にステップＳ３００の処理の一例がより詳細に説明されている。図１５のフローチャートは、ステップＳ３００の処理の一例をより詳細に説明したフローチャートである。

ステップＳ３１０において、処理回路１５０は、算出機能１５０ｂにより、ステップＳ１００で取得された第１の検出器ｌにおける第１の光子数情報に基づいて、発光の分布関数モデルpⁱ _ｐｈ，ｌ (t)を生成する。発光の分布関数p_ｐｈ (t)は、時刻ｔに発光が起こる確率を示す分布関数である。発光の分布関数モデルpⁱ _ｐｈ,l (t)は、時刻に関して全時間で積分すると１になる。また、発光の分布関数モデルpⁱ _ｐｈ，ｌ(t)は、検出器ごと及び同時計数イベントごとに定められる。

処理回路１５０は、算出機能１５０ｂにより、第１の光子数情報に基づいて、例えばすでに第１の実施形態で述べたように、式（１）、（４）、（５）、（６）等を用いて発光の分布関数モデルpⁱ _ｐｈ,l (t)を生成する。第１の光子数情報が、検出器ごと及び同時計数イベントごとに取得されているので、処理回路１５０は、算出機能１５０ｂにより、発光の分布関数モデルp_ｐｈ,l ^ｉ (t)を、検出器ごと及び同時計数イベントごとに算出する。第１の実施形態においてすでに述べたように、図６に、発光の分布関数p_ｐｈ(t)の一例が示されている。図６において、曲線３０は、発光の分布関数p_ｐｈ(t)を示している。また、検出時刻を定義するためのフォトンの閾値（Ｎ’）を併記する。

続いて、ステップＳ３２０Ａにおいて、処理回路１５０は、算出機能１５０ｂにより、発光の分布関数モデルpⁱ _ｐｈ,l (θⁱ _ｌ,t)に基づいて、第１の発光確率モデルｐⁱ _ｔ,l(θⁱ _ｌ, t)を算出する。

ここで、第１の発光確率モデルｐⁱ _ｔ,l(≧Ｎ’，θⁱ _ｌ, t)は、検出事象に含まれる複数の光子の内の最初のＮ’個以上の光子を検出する確率密度に関するモデルである。換言すると、第１の発光確率モデルｐ_ｔ,l(≧Ｎ’，t)は、検出事象に含まれる複数の光子の内の所定の数の光子を検出する確率密度に関するモデルである。

第１の実施形態においてすでに述べたように、式（７）～式（１３）が成り立つ。

すなわち、処理回路１５０は、算出機能１５０ｂにより、式（１２）及び式（１３）を用いることで、発光の分布関数ｐ_ph(t)に基づいて、第１の発光確率モデルp_ｔ(t)を計算することができる。このように、処理回路１５０は、算出機能１５０ｂにより、式（１２）及び式（１３）に第１の検出器における発光の分布関数モデルp_ｐｈ (t)を代入し、第１の光子数情報を式（１３）の右辺のＮに代入して式（１３）の右辺を評価することで、第１の光子数情報及び同時計数イベントごとに定められた第１の検出器における発光の分布関数モデルp_ｐｈ (t)に基づいて、第１のイベントについての第１の発光確率モデルｐ_ｔ(t)を算出する。

なお、上述の発光の分布関数ｐ_ph(t)及び第１の第１の発光確率モデルp_t(t)は、各検出器ごと及び同時計数イベントごとに定めることができる。従って、検出器の添え字及び同時計数イベントに関する添え字を明示的に記載すると、式（１３）は、以下の式（３２）のように書くこともできる。

ここで、p^ｉ _ｔ；ｌ(t)はｉ番目の同時計数イベントの検出器ｌにおける第１の発光確率モデル、ｐ^ｉ _ｐｈ；ｌ(t)はｉ番目の同時計数イベントの検出器ｌにおける発光の分布関数、Ｎ^ｉ _ｓ；ｌはｉ番目の同時計数イベントの検出器ｌにおけるシンチレーション光の光子数、Ｎ^ｉ _ｃ；ｌはｉ番目の同時計数イベントの検出器ｌにおけるチェレンコフ光の光子数を示す。

なお、第１の発光確率モデルを算出の際、処理回路１５０は、算出機能１５０ｂにより、シンチレータ内での発光位置の不定性を考慮した計算を行うこともできる。

この計算については、第１の実施形態において、図８及び式（１５）～（１８）を用いて議論した。すでに述べたように、式（１５）～（１７）が成り立つ。

第１の実施形態と同様に、式（１３）を基に、シンチレータ内の発光位置ｘを考慮した定式化を考えると、シンチレータ内の発光位置ｘを考慮した取り扱いを考える場合、以下の式（３３）が成り立つ。

ここで、p_ｔ、ｌ ^ｉ,dｏｉ(t)は、同時計数イベントごとに定められた、シンチレータ内での発光位置の不定性を考慮した第１の発光確率モデルであり、ｐ_ｐｈ、ｌ ^{ｉ,ｄｏｉ}(t)は、同時計数イベントごとに定められた、シンチレータの位置の発光位置の不定性を考慮した発光の分布関数モデルである。

なお、式（３３）は、処理回路１５０が取得機能１５０ａによりステップＳ１００において取得するパラメータθ_ｌを関数の引数として明示的に表現すると、以下の式（３４）及び式（３５）のようにも書くことができる。

ここで、θ^ｉ _ｌは、同時計数イベントごとに定められた、検出器ｌにおいて取得されるパラメータセットを示しており、例えば、同時計数イベントごとに定められたシンチレーション光子数Ｎ^ｉ _ｓ；ｌ、チェレンコフ光子数Ｎ^ｉ _ｃ；ｌ、チェレンコフ光に対する光センサの時間分解能σ^ｉ _Ｃ、シンチレーション光に対する光センサの時間分解能σ^ｉ _Ｓ等のパラメータの組を表す。

式（３３）に戻り、シンチレータの位置の発光位置の不定性を考慮した場合のステップＳ３００の処理について説明する。図１６は、シンチレータ内での発光位置の不定性を考慮した場合の、発光確率モデルの算出の流れを示したフローチャートである。

図１５で示した場合と同様に、ステップＳ３１０において、処理回路１５０は、算出機能１５０ｂにより、第１の検出器ｌにおける第１の光子数情報等のパラメータθ_ｌに基づいて、同時計数イベントごとの発光の分布関数モデルｐ_ph,l(t)を生成する。

続いて、ステップＳ３１１において、処理回路１５０は、算出機能１５０ｂにより、シンチレータ内での位置ｘごとに、発光確率密度ｐ_ｄｏｉ（ｘ）を、減弱係数μの値に基づいて算出する。位置ｘにおける発光確率密度ｐ_ｄｏｉ（ｘ）は、式（２１）で表される。なお、発光確率密度ｐ_ｄｏｉ（ｘ）は、式（２２）に示されるように、規格化されている。

図１６に戻り、ステップＳ３１２において、処理回路１５０は、特定機能１５０ｃにより、ステップＳ３１０で生成された、同時計数イベントごとに定められた発光の分布関数である分布関数モデルｐ^ｉ _ｐｈ、ｌ（ｔ）と発光確率密度ｐ_ｄｏｉ（ｘ）とに基づいて、同時計数イベントごとに定められた発光位置不確定性を考慮した発光の分布関数モデルｐ_ｐｈ；ｌ ^i,doi(ｔ)を第１の発光確率モデルとして算出する。具体的には、処理回路１５０は、特定機能１５０ｃにより、発光位置不確定性を考慮した発光の分布関数モデルｐ^{ｉ，ｄｏｉ} _ｐｈ；ｌ（ｔ）を、ステップＳ３１０で生成された分布関数モデルｐⁱ _ｐｈ、ｌ（ｔ）と、ステップＳ３１１で算出された発光確率密度ｐ_ｄｏｉ（ｘ）とに基づいて、式（２３）に基づいて算出する。

続いて、ステップＳ３２０Ｂにおいて、処理回路１５０は、特定機能１５０ｃにより、同時計数イベントごとに定められた、発光位置不確定性を考慮した発光の分布関数モデルｐ^{i，ｄｏｉ} _ｐｈ；ｌ（ｔ）に基づいて、同時計数イベントごとに定められた、発光位置不確定性を考慮した発光確率モデルｐ^{ｉ,ｄｏｉ} _ｔ；ｌ（θ_l;t）を算出する。具体的には、式（２４）及び式（２５）に示されるように、処理回路１５０は、特定機能１５０ｃにより、式（１３）と同様に、同時計数イベントごとに定められた発光位置不確定性を考慮した発光の分布関数モデルｐ^i,doi _ｐｈ（ｔ）に基づいて、同時計数イベントごとに定められた発光位置不確定性を考慮した発光確率モデルｐ^i,dｏｉ _t（t）を算出する。

なお、式（２３）において、式（１７）を基に被積分関数のｘをｔ’に変数変換すると、式（２６）が得られるので、処理回路１５０は、特定機能１５０ｃにより、式（２４）の右辺の評価において、式（２６）を用いて算出を行ってもよい。

このようにして、処理回路１５０は、特定機能１５０ｃにより、シンチレータ内での発光位置の不定性を考慮して、第１の発光確率モデルおよび第２の発光確率モデルを算出する。

図１４に戻り、ステップＳ４００において、処理回路１５０は、特定機能１５０ｃにより、ステップＳ３００と同様の処理を行うことにより、ステップＳ２００で取得された第２の光子数情報に基づいて、第２の検出器ｍに対する、同時計数イベントごとに定められた第２の発光確率モデルp^{ｉ、（ｄｏｉ）} _ｔ；ｍを算出する。

具体的には、処理回路１５０は、特定機能１５０ｃにより、ステップＳ２００で取得された第２の検出器ｍにおける同時計数イベントごとに定められた第２の光子数情報に基づいて、例えば式（１）、（４）、（５）、（６）等に相当する式を用いて、同時計数イベントごとに定められた発光の分布関数モデルp^ｉ _ｐｈ；ｍ (t)を生成する。

続いて、処理回路１５０は、特定機能１５０ｃにより、第２の光子数情報及び第２の検出器ｍにおける発光の分布関数モデルp^ｉ _ｐｈ；ｍ(t)に基づいて、同時計数イベントごとに定められた第２の発光確率モデルｐ^ｉ _ｔ；ｍ(t)を、例えば式（１４）等に相当する式を用いて算出する。

また、第２の発光確率モデルを算出の際、処理回路１５０は、算出機能１５０ｂにより、シンチレータ内での発光位置の不定性を考慮した計算を、式（１７）及び式（２１）～（２６）と同様の式を用いて行うこともできる。この場合、処理回路１５０は、算出機能１５０ｂにより、シンチレータ内での位置ｘごとに、発光確率密度ｐ_ｄｏｉ（ｘ）を、減弱係数μの値に基づいて算出する。続いて、処理回路１５０は、算出機能１５０ｂにより、同時計数イベントごとに定められた分布関数モデルｐ^ｉ _ｐｈ；ｍ（ｔ）と発光確率密度ｐ_ｄｏｉ（ｘ）とに基づいて、同時計数イベントごとに定められた発光位置不確定性を考慮した発光の分布関数モデルｐ^{ｉ、ｄｏｉ} _ｐｈ；ｍ（ｔ）を算出する。続いて、処理回路１５０は、特定機能１５０ｃにより、同時計数イベントごとに定められた発光位置不確定性を考慮した発光の分布関数モデルｐ^{ｉ、ｄｏｉ} _ｐｈ；ｍ（ｔ）に基づいて、同時計数イベントごとに定められた発光位置不確定性を考慮した発光確率モデルｐ^{ｉ、ｄｏｉ} _ｔ；ｍ（t）を算出する。

このようにして、処理回路１５０は、特定機能１５０ｃにより、シンチレータ内での発光位置の不定性を考慮して、同時計数イベントごとに定められた第２の発光確率モデルｐ_ｔ,ｍ(t)を算出する。

続いて、図１３に戻り、ステップＳ５００において、処理回路１５０は、特定機能１５０ｃにより、ステップS３００で算出された第１の発光確率モデルｐ_ｔ；ｌ ^ｄｏｉ（t_１）と算出された第２の発光確率モデルｐ_ｔ；ｍ ^ｄｏｉ（t_２）とに基づいて、同時計数イベントごとに特定された、第１のイベントと第２のイベントとに基づいて規定されるコインシデンスに関する確率分布モデルｐ^{ｉ、（ｄｏｉ）} _ｃｔ;l,m(θ_l, θ_ｍ,t)を特定する。

図１７において、曲線２５３は、第１の検出器ｌのシンチレータ２５１における同時計数イベントごとに定められた第１の発光確率モデルｐ^{ｉ，ｄｏｉ} _ｔ；ｌ（θ_l,t_l）を、遅延時間t_lの関数としてプロットしたものを示している。また、曲線２５４は、第２の検出器ｍのシンチレータ２５２における第２の発光確率モデルｐ^{ｉ，ｄｏｉ} _ｔ；ｍ（θ_ｍ,t_ｍ）を、遅延時間ｔ_ｍの関数としてプロットしたものを示している。処理回路１５０は、特定機能１５０ｃにより、検出器ｌにおける検出時刻ｔ_ｌにおける事象と、検出器ｍにおける検出時刻ｔ_ｍにおける積事象を、検出時間差の関数として表現した、同時計数イベントごとに定められた確率分布モデルｐ^{ｉ，ｄｏｉ} _ｃｔ;l,m(θ_l, θ_ｍ,t)を、特定する。曲線２５６は、このような確率分布モデルｐ^{ｉ，ｄｏｉ} _ｃｔ;l,m(θ_l, θ_ｍ,t)の例を示している。

ここで、時間幅ｄｔの間における、検出器ｌにおける検出時刻ｔ_ｌにおける事象と、検出器ｍにおける検出時刻ｔ_ｍにおける積事象の数は、以下の式（３６）で与えられる。

ここで、θ^ｉ _ｌは、ｉ番目の同時計数イベントに係る、検出器ｌにおいて取得されるパラメータセットを示している。また、パラメータθを、より具体的な形で表現すると、以下の式（３７）が成り立つ。

ここで、Ｎ^ｉ _ｓ；ｌ、Ｎ^ｉ _ｓ；ｍは、ｉ番目の同時計数イベントに係る、検出器ｌ及びｍにおけるシンチレーション光子数であり、Ｎ^ｉ _ｃ；ｌ、Ｎ^ｉ _ｃ；ｍは、ｉ番目の同時計数イベントに係る、検出器ｌ及びｍにおけるチェレンコフ光子数を示す。

また、検出器ｌにおける検出時刻ｔ_ｌと、検出器ｍにおける検出時刻ｔ_ｍとの検出時間の差は、すでに述べた式（２８）で与えられる。

従って、ｐ^{ｉ，ｄｏｉ} _ｃｔ;l,m(θ_l, θ_ｍ,t)は、以下の式（３８）で与えられる。

また、パラメータθを、より具体的な形で表現すると、以下の式（３９）が成り立つ。

以上のように、図１３のステップＳ３０の処理により、処理回路１５０は、特定機能１５０ｃにより、第１の検出器ｌと第２の検出器ｍとによって規定されるＬＯＲの検出時間差に関する確率分布モデルｐ^{ｉ、ｄｏｉ} _ｃｔ;l,m(θ^ｉ _l, θ^ｉ _ｍ,t)を、ｉ番目の同時計数イベントについて算出する。処理回路１５０は、この処理を、全ての同時計数イベントについて算出し、それらを加算することで、最終的な、第１の検出器ｌと第２の検出器ｍとによって規定されるＬＯＲの検出時間差に関する確率分布モデルｐ^ｄｏｉ _ｃｔ;l,m(θ^ｉ _l, θ^ｉ _ｍ,t)を算出し、それを基に、データの再構成を行うことができる。

ステップＳ４０において、処理回路１５０は、特定機能１５０ｃにより、すべての同時計数イベントについて、ステップＳ３０の処理が終わったか否かを判定する。処理回路１５０が、特定機能１５０ｃにより、すべての同時計数イベントについて、ステップＳ３０の処理が終わっていないと判定した場合（ステップＳ４０ Ｎｏ）、処理はステップＳ５０に進み、ｉ＋１番目のデータについて、処理が行われる。一方、処理回路１５０が、特定機能１５０ｃにより、すべての同時計数イベントについて、ステップＳ３０の処理が終わったと判定した場合（ステップＳ４０ Ｙｅｓ）、処理はステップＳ６０に進む。

ステップＳ６０において、処理回路１５０は、特定機能１５０ｃにより、ステップＳ３０で算出された、各同時計数イベントごとの確率分布モデルｐ^{ｉ、ｄｏｉ} _ｃｔ;l,m(θ^ｉ _l, θ^ｉ _ｍ,t)を加算して、第１の検出器ｌと第２の検出器ｍにおけるコインシデンスに関する確率分布モデルｐ^ｔｏｔ _{ｃｔ；ｌ、ｍ}を生成する。

具体的には、第１の検出器ｌと第２の検出器ｍにおけるコインシデンスに関する確率分布モデルｐ^ｔｏｔ _{ｃｔ；ｌ、ｍ}は、以下の式（４０）で表される。

ここで、第１の検出器ｌと第２の検出器ｍにおけるコインシデンスに関する確率分布モデルｐ^ｔｏｔ _{ｃｔ；ｌ、ｍ}は、後述のステップＳ７０において、ＴｏＦカーネルとして用いられることで、画像再構成に用いられる。すなわち、処理回路１５０は、特定機能１５０ｃにより、各イベントごとに特定された確率分布モデルｐ^{ｉ、ｄｏｉ} _ｃｔ;l,m(θ^ｉ _l, θ^ｉ _ｍ,t)を、各イベントについて加算することで、第１の検出器ｌと第２の検出器ｍにおけるコインシデンスに関する確率分布モデルｐ^ｔｏｔ _{ｃｔ；ｌ、ｍ}を、第１の検出器ｌと第２の検出器ｍとの間のＴｏＦカーネルとして生成する。このように、処理回路１５０は、特定機能１５０ｃにより、第１の光子数情報と第２の光子数情報とに基づいて、第１のイベントと第２のイベントとに基づいて規定されるコインシデンスに関する確率分布モデルｐ^ｔｏｔ _{ｃｔ；ｌ、ｍ}を特定する。

ここで、第２の実施形態では、処理回路１５０は、特定機能１５０ｃにより、ステップＳ３０において同時計数イベントごとに確率分布モデルｐ^{ｉ、ｄｏｉ} _ｃｔ;l,m(θ^ｉ _l, θ^ｉ _ｍ,t)を生成した上で、これらをステップＳ６０において加算して、最終的な確率分布モデルｐ^ｔｏｔ _{ｃｔ；ｌ、ｍ}を生成している。このことにより、処理回路１５０は、特定機能１５０ｃにより、事象ごとの観測値の平均からの揺らぎを考慮したＴｏＦカーネルを生成することができる。この結果、処理回路１５０は、特定機能１５０ｃにより、例えば少ない事象で観測揺らぎが大きく見えてしまう場合等でも、妥当なＴｏＦカーネルを生成することができる。

続いて、ステップＳ７０において、処理回路１５０は、再構成機能１５０ｄにより、ステップＳ６０で特定された、第１の光子数情報と第２の光子数情報とに基づいて規定されたコインシデンスに関する確率分布モデルｐ^ｔｏｔ _{ｃｔ；ｌ、ｍ}をＴＯＦカーネルとして利用することで、第１の検出器ｌと第２の検出器ｍとの間のＬＯＲに対応するデータの再構成を行う。

図１８～図２０に、第２の実施形態に係る手法の検証結果が示されている。この検証では、リストモードデータから、本手法を用いて確率分布モデルｐ^ｔｏｔ _{ｃｔ；ｌ、ｍ}を生成し、実測値と比較した。

図１８に、ＢＧＯ結晶を用いて測定されたリストモードの散布図が示されており、図中の丸一つ一つ、例えばイベント２６０及びイベント２６１が一つのイベントに対応する。処理回路１５０は、算出機能１５０ｂにより、イベント２６０を基に、イベント２６０に対応する、イベントごとの確率分布モデルｐ_ｃｔ（ｔ）を算出した。曲線２６３は、イベント２６０に対応するｐ_ｃｔ（ｔ）である。また、処理回路１５０は、算出機能１５０ｂにより、イベント２６１を基に、イベント２６１に対応する、イベントごとの確率分布モデルｐ_ｃｔ（ｔ）を算出した。曲線２６２は、イベント２６１に対応するｐ_ｃｔ（ｔ）である。処理回路１５０は、算出機能１５０ｂにより、同様な処理をすべてのイベントについて行い、イベントごとの確率分布モデルｐ_ｃｔ（ｔ）を加算することにより、確率分布モデルｐ^ｔｏｔ（ｔ）を生成した。

図１９は、ＬＹＳＯにおける実測値と確率分布モデルｐ^ｔｏｔ _ｃｔ（ｔ）との比較を示している。曲線２６５が実測値であり、曲線２６４が、計算された確率分布モデルｐ^ｔｏｔ _ｃｔ（ｔ）である。また、図２０は、ＢＧＯにおける実測値と確率分布モデルｐ^ｔｏｔ _ｃｔ（ｔ）との比較を示している。曲線２６６が実測値であり、曲線２６７が、計算された確率分布モデルｐ^ｔｏｔ _ｃｔ（ｔ）である。どちらの場合でも、実測値と確率分布モデルｐ^ｔｏｔ（ｔ）とが概ね一致しており、同時計数イベントごとに発光モデルをもとに確率分布モデルｐ_ｃｔ（ｔ）を生成し、それらを加算して確率分布モデルｐ^ｔｏｔ _ｃｔ（ｔ）を生成する第２の実施形態の手法の妥当性が確認された。

実施形態は、上述の例に限られない。上述の場合においては、検出器３がモノリシック検出器でない通常の検出器である場合について説明したが、実施形態はこれらに限られない。実施形態は、検出器３が、単結晶として構成され、１つのシンチレータ結晶が区分けされずに一体化したまま使用されるモノリシックシンチレータである場合であってもよい。当該シンチレータは、例えば、ＬＹＳＯ（Lutetium Yttrium Oxyorthosilicate）、ＬＳＯ（Lutetium Oxyorthosilicate）、ＬＧＳＯ（Lutetium Gadolinium Oxyorthosilicate）等の、シンチレータ結晶によって形成される。

検出器３がモノリシック検出器である場合、検出器３は一つのシンチレータ結晶が区分けされずに一体化したまま使用されるものであるが、検出器３が複数の仮想的な検出器からなると考えて、それぞれの仮想的な検出器の間のＴｏＦカーネルを生成することにより、同様の処理を行うことができる。例えば、検出器３が複数の検出素子アレイからなる場合、それぞれの検出素子アレイを１つの検出器と考えて、図１３のフローチャートで示された処理と同様の処理を行うことができる。一例として、図１３のステップＳ１０からステップＳ７０や図１４のステップＳ１００からステップＳ３００、式（３５）、（３８）等における自然数ｌ及びｍを、検出器に付された番号の代わりに、検出素子アレイに付された番号であると読み替えて同様な処理を行うことで、処理回路１５０は、特定機能当該検出素子アレイの組み合わせに応じた検出時間差に関する確率分布モデルを特定することができる。

すなわち、検出器３がモノリシック検出器である場合、第２の実施形態における核医学診断装置は、処理回路１５０を備える。処理回路１５０は、取得機能１５０ａにより、少なくとも第１の検出素子アレイｌと第２の検出素子アレイｍとが設けられたモノリシック検出器と、第１の検出素子アレイｌにて検出される第１のイベントに対応する第１の光子数情報と、第２の検出素子アレイｍにて検出され第１のイベントに対応する第２の光子数情報とを取得する。処理回路１５０は、特定機能１５０ｃにより、第１の光子数情報と第２の光子数情報とに基づいて、第１の検出素子アレイｌと第２の検出素子アレイｍとの組み合わせに応じた検出時間差に関する確率分布モデルｐ^ｔｏｔ _{ｃｔ；ｌ、ｍ}を特定する。

また、第２の実施形態に係るデータ処理方法は、少なくとも第１の検出素子アレイｌと第２の検出素子アレイｍとが設けられたモノリシック検出器の第１の検出素子アレイｌにて検出される第１のイベントに対応する第１の光子数情報と、第２の検出素子アレイｍにて検出され第１のイベントに対応する第２の光子数情報とを計測し、第１の光子数情報と第２の光子数情報とに基づいて、第１の検出素子アレイｌと第２の検出素子アレイｍとの組み合わせに応じた検出時間差に関する確率分布モデルｐ^ｔｏｔ _{ｃｔ；ｌ、ｍ}を特定する。

また、第２の実施形態に係るプログラムは、少なくとも第１の検出素子アレイｌと第2の検出素子アレイｍとが設けられたモノリシック検出器の第１の検出素子アレイｌにて検出される第１のイベントに対応する第１の光子数情報と、第２の検出素子アレイｍにて検出され第１のイベントに対応する第２の光子数情報とを計測し、第１の光子数情報と第２の光子数情報とに基づいて、第１の検出素子アレイｌと第２の検出素子アレイｍとの組み合わせに応じた検出時間差に関する確率分布モデルｐ^ｔｏｔ _{ｃｔ；ｌ、ｍ}を特定する処理をコンピューターに実行させる。

以上のように、第２の実施形態に係る核医学診断装置においては、処理回路１５０は、シンチレータ内の発光モデルに基づいて、検出時間差に関する確率分布モデルを、同時計数イベントごとに算出した上でそれらを加算して検出時間差に関する確率分布モデルを生成した。これにより、事象ごとの観測値の平均からの揺らぎを考慮したＴｏＦカーネルを生成して画像再構成をすることができるので、例えば少ない事象で観測揺らぎが大きく見えてしまう場合等でも、再構成後のＰＥＴ画像の画質を向上させることができる。

以上説明した少なくとも１つの実施形態によれば、画質を向上させることができる。

以上の実施形態に関し、発明の一側面および選択的な特徴として以下の付記を開示する。

(付記１)
本発明の一つの側面において提供される核医学診断装置は、取得部と、算出部と、特定部とを備える。取得部は、第１の検出器で検出される第１の光子数情報と、前記第１の検出器とは異なる第２の検出器で検出される第２の光子数情報とを取得する。算出部は、前記第１の光子数情報に基づいて前記第１の検出器に対応する第１の発光確率モデルと、前記第２の光子数情報に基づいて前記第２の検出器に対応する第２の発光確率モデルとを算出する。特定部は、前記第１の発光確率モデルと前記第２の発光確率モデルとに基づいて、前記第１の検出器と前記第２の検出器とによって規定されるＬＯＲ（Ｌｉｎｅ Ｏｆ Ｒｅｓｐｏｎｓｅ）の検出時間差に関する確率分布モデルを特定する。

（付記２）
前記検出時間差に関する確率分布モデルをＴｏＦ（Ｔｉｍｅ Ｏｆ Ｆｌｉｇｈｔ）カーネルとして利用することで、前記ＬＯＲに対応するデータの再構成を行う再構成処理部を更に備えてもよい。

（付記３）
前記第１の発光確率モデル及び前記第２の発光確率モデルは、検出事象に含まれる複数の光子の内の最初の光子を検出する確率密度に関するモデルであってもよい。

（付記４）
前記第１の発光確率モデル及び前記第２の発光確率モデルは、検出事象に含まれる複数の光子の内の所定の数の光子を検出する確率密度に関するモデルであってもよい。

（付記５）
前記第１の光子数情報及び前記第２の光子数情報は、所定のキャリブレーションスキャンを行うことで特定されてもよい。

（付記６）
前記第１の光子数情報及び前記第２の光子数情報は、前記キャリブレーションスキャンを行うことで前記第１の検出器及び前記第２の検出器の各々で取得されるエネルギースペクトルに基づいて特定されてもよい。

（付記７）
前記第１の発光確率モデル及び前記第２の発光確率モデルは、前記第１の検出器及び前記第２の検出器の各々における発光の分布関数モデルに基づいてもよい。

（付記８）
前記第１の発光確率モデル及び前記第２の発光確率モデルは、検出事象に対応する光子数の期待値に基づいてもよい。

（付記９）
前記第１の光子数情報及び前記第２の光子数情報は、各々シンチレーション光子数を含んでもよい。

（付記１０）
前記第１の光子数情報及び前記第２の光子数情報は、各々チェレンコフ光子数を更に含んでもよい。

（付記１１）
前記算出部は、シンチレータ内での発光位置の不定性を考慮して、前記第１の発光確率モデルおよび前記第２の発光確率モデルを算出してもよい。

（付記１２）
前記算出部は、前記シンチレータ内での位置ごとに発光確率密度を減弱係数に基づいて算出し、発光の分布関数及び前記発光確率密度に基づいて、前記第１の発光確率モデルおよび前記第２の発光確率モデルを算出してもよい。

（付記１３）
本発明の一つの側面において提供されるデータ処理方法は、
第１の検出器で検出される第１の光子数情報と、前記第１の検出器とは異なる第２の検出器で検出される第２の光子数情報とを取得し、
前記第１の光子数情報に基づいて前記第１の検出器に対応する第１の発光確率モデルと、前記第２の光子数情報に基づいて前記第２の検出器に対応する第２の発光確率モデルとを算出し、
前記第１の発光確率モデルと前記第２の発光確率モデルとに基づいて、前記第１の検出器と前記第２の検出器とによって規定されるＬＯＲ（Ｌｉｎｅ Ｏｆ Ｒｅｓｐｏｎｓｅ）の検出時間差に関する確率分布モデルを特定する。

（付記１４）
本発明の一つの側面において提供されるプログラムは、第１の検出器で検出される第１の光子数情報と、前記第１の検出器とは異なる第２の検出器で検出される第２の光子数情報とを取得し、前記第１の光子数情報に基づいて前記第１の検出器に対応する第１の発光確率モデルと、前記第２の光子数情報に基づいて前記第２の検出器に対応する第２の発光確率モデルとを算出し、前記第１の発光確率モデルと前記第２の発光確率モデルとに基づいて、前記第１の検出器と前記第２の検出器とによって規定されるＬＯＲの検出時間差に関する確率分布モデルを特定する処理をコンピュータに実行させる。

（付記１５）
本発明の一つの側面において提供される核医学診断装置は、少なくとも第１の検出素子アレイと第２の検出素子アレイとが設けられたモノリシック検出器と、前記第１の検出素子アレイにて検出される第１の光子数情報と、前記第２の検出素子アレイにて検出される第２の光子数情報とを取得する取得部と、前記第１の光子数情報に基づいて前記第１の検出素子アレイに対応する第１の発光確率モデルと、前記第２の光子数情報に基づいて前記第２の検出素子アレイに対応する第２の発光確率モデルとを算出する算出部と、前記第１の発光確率モデルと前記第２の発光確率モデルとに基づいて、前記第１の検出素子アレイと前記第２の検出素子アレイとの組み合わせに応じた検出時間差に関する確率分布モデルを特定する特定部と、を備える。

（付記１６）
本発明の一つの側面において提供されるデータ処理方法は、少なくとも第１の検出素子アレイと第２の検出素子アレイとが設けられたモノリシック検出器の前記第１の検出素子アレイにて検出される第１の光子数情報と、前記第２の検出素子アレイにて検出される第２の光子数情報とを取得し、前記第１の光子数情報に基づいて前記第１の検出素子アレイに対応する第１の発光確率モデルと、前記第２の光子数情報に基づいて前記第２の検出素子アレイに対応する第２の発光確率モデルとを算出し、前記第１の発光確率モデルと前記第２の発光確率モデルとに基づいて、前記第１の検出素子アレイと前記第２の検出素子アレイとの組み合わせに応じた検出時間差に関する確率分布モデルを特定する。

（付記１７）
本発明の一つの側面において提供されるプログラムは、少なくとも第１の検出素子アレイと第２の検出素子アレイとが設けられたモノリシック検出器の前記第１の検出素子アレイにて検出される第１の光子数情報と、前記第２の検出素子アレイにて検出される第２の光子数情報とを取得し、前記第１の光子数情報に基づいて前記第１の検出素子アレイに対応する第１の発光確率モデルと、前記第２の光子数情報に基づいて前記第２の検出素子アレイに対応する第２の発光確率モデルとを算出し、前記第１の発光確率モデルと前記第２の発光確率モデルとに基づいて、前記第１の検出素子アレイと前記第２の検出素子アレイとの組み合わせに応じた検出時間差に関する確率分布モデルを特定する処理をコンピュータに実行させる。

(付記１８)
本発明の一つの側面において提供される核医学診断装置は、取得部と、特定部とを備える。取得部は、第１の検出器で検出される第１のイベントに関する第１の光子数情報と、前記第１の検出器とは異なる第２の検出器で検出される第２のイベントに関する第２の光子数情報とを取得する。特定部は、前記第１の光子数情報と前記第２の光子数情報とに基づいて、前記第１のイベントと前記第２のイベントとに基づいて規定されるコインシデンスに関する確率分布モデルを特定する。

（付記１９）
前記コインシデンスに関する確率分布モデルをＴｏＦ（Ｔｉｍｅ Ｏｆ Ｆｌｉｇｈｔ）カーネルとして利用することで、前記第１の検出器と前記第２の検出器との間のＬＯＲ（Ｌｉｎｅ Ｏｆ Ｒｅｓｐｏｎｓｅ）に対応するデータの再構成を行う再構成処理部を更に備えてもよい。

（付記２０）
前記取得部は、前記第１の光子数情報及び前記第２の光子数情報として、シンチレーション光子数を計測してもよい。

（付記２１）
前記取得部は、前記第１の光子数情報及び前記第２の光子数情報として、チェレンコフ光子数を更に計測してもよい。

（付記２２）
前記特定部は、各イベントごとに特定された前記確率分布モデルを加算して、前記第１の検出器と前記第２の検出器との間のＴｏＦカーネルを生成してもよい。

（付記２３）
前記特定部は、前記第１の光子数情報に基づいて、前記第１のイベントについての第１の発光確率モデルを生成し、前記第２の光子数情報に基づいて、前記第２のイベントについての第２の発光確率モデルを生成し、前記第１の発光確率モデルと前記第２の発光確率モデルとに基づいて、前記確率分布モデルを特定してもよい。

（付記２４）
前記特定部は、イベントごとに定められた発光分布関数モデルに基づいて、前記第１の発光確率モデルおよび前記第２の発光確率モデルを生成してもよい。

（付記２５）
前記取得部により計測される前記第１の光子数情報及び前記第２の光子数情報は、イベントごとに取得されるものであってもよい。

（付記２６）
前記特定部は、シンチレーター内での発光位置の不定性を考慮して、前記第１の発光確率モデルおよび前記第２の発光確率モデルを算出してもよい。

（付記２７）
前記特定部は、前記シンチレーター内での位置ごとに発光確率密度を減弱係数に基づいて算出し、発光の分布関数及び前記発光確率密度に基づいて、前記第１の発光確率モデルおよび前記第２の発光確率モデルを算出してもよい。

（付記２８）
前記第１の発光確率モデルは、検出事象に含まれる複数の光子の内の所定の数の光子を検出する確率密度に関するモデルであってもよい。

（付記２９）
本発明の一つの側面において提供されるデータ処理方法は、第１の検出器で検出される第１のイベントに関する第１の光子数情報と、前記第１の検出器とは異なる第２の検出器で検出される第２のイベントに関する第２の光子数情報とを計測し、前記第１の光子数情報と前記第２の光子数情報とに基づいて、前記第１のイベントと前記第２のイベントとに基づいて規定されるコインシデンスに関する確率分布モデルを特定する。

（付記３０）
本発明の一つの側面において提供されるプログラムは、第１の検出器で検出される第１のイベントに関する第１の光子数情報と、前記第１の検出器とは異なる第２の検出器で検出される第２のイベントに関する第２の光子数情報とを計測し、前記第１の光子数情報と前記第２の光子数情報とに基づいて、前記第１のイベントと前記第２のイベントとに基づいて規定されるコインシデンスに関する確率分布モデルを特定する処理をコンピューターに実行させる。

（付記３１）
本発明の一つの側面において提供される核医学診断装置は、取得部と、特定部とを備える。取得部は、少なくとも第１の検出素子アレイと第２の検出素子アレイとが設けられたモノリシック検出器と、前記第１の検出素子アレイにて検出される第１のイベントに対応する第１の光子数情報と、第２の検出素子アレイにて検出され前記第１のイベントに対応する第２の光子数情報とを取得する。特定部は、前記第１の光子数情報と前記第２の光子数情報とに基づいて、前記第１の検出素子アレイと前記第２の検出素子アレイとの組み合わせに応じた検出時間差に関する確率分布モデルを特定する。

（付記３２）
本発明の一つの側面において提供されるデータ処理方法は、少なくとも第１の検出素子アレイと第２の検出素子アレイとが設けられたモノリシック検出器の前記第１の検出素子アレイにて検出される第１のイベントに対応する第１の光子数情報と、第２の検出素子アレイにて検出され前記第１のイベントに対応する第２の光子数情報とを計測し、前記第１の光子数情報と前記第２の光子数情報とに基づいて、前記第１の検出素子アレイと前記第２の検出素子アレイとの組み合わせに応じた検出時間差に関する確率分布モデルを特定する。

（付記３３）
本発明の一つの側面において提供されるプログラムは、少なくとも第１の検出素子アレイと第２の検出素子アレイとが設けられたモノリシック検出器の前記第１の検出素子アレイにて検出される第１のイベントに対応する第１の光子数情報と、第２の検出素子アレイにて検出され前記第１のイベントに対応する第２の光子数情報とを計測し、前記第１の光子数情報と前記第２の光子数情報とに基づいて、前記第１の検出素子アレイと前記第２の検出素子アレイとの組み合わせに応じた検出時間差に関する確率分布モデルを特定する処理をコンピューターに実行させる。

いくつかの実施形態を説明したが、これらの実施形態は、例として提示したものであり、発明の範囲を限定することは意図していない。これら実施形態は、その他の様々な形態で実施されることが可能であり、発明の要旨を逸脱しない範囲で、種々の省略、置き換え、変更、実施形態同士の組み合わせを行うことができる。これら実施形態やその変形は、発明の範囲や要旨に含まれると同様に、特許請求の範囲に記載された発明とその均等の範囲に含まれるものである。

１５０ 処理回路
１５０ａ 取得機能
１５０ｂ 算出機能
１５０ｃ 特定機能
１５０ｄ 再構成機能
１５０ｅ 制御機能
１５０ｆ 受付機能
１５０ｇ 画像生成機能
１５０ｈ 表示制御機能

Claims (23)

1. 第１の検出器で検出される第１の光子数情報と、前記第１の検出器とは異なる第２の検出器で検出される第２の光子数情報とを取得する取得部と、
   前記第１の光子数情報に基づいて前記第１の検出器に対応する第１の発光確率モデルと、前記第２の光子数情報に基づいて前記第２の検出器に対応する第２の発光確率モデルとを算出する算出部と、
   前記第１の発光確率モデルと前記第２の発光確率モデルとに基づいて、前記第１の検出器と前記第２の検出器とによって規定されるＬＯＲ（Ｌｉｎｅ Ｏｆ Ｒｅｓｐｏｎｓｅ）の検出時間差に関する確率分布モデルを特定する特定部と、を備える、核医学診断装置。
2. 前記検出時間差に関する確率分布モデルをＴｏＦ（Ｔｉｍｅ Ｏｆ Ｆｌｉｇｈｔ）カーネルとして利用することで、前記ＬＯＲに対応するデータの再構成を行う再構成処理部を更に備える、請求項１に記載の核医学診断装置。
3. 前記第１の発光確率モデル及び前記第２の発光確率モデルは、検出事象に含まれる複数の光子の内の最初の光子を検出する確率密度に関するモデルである、請求項１に記載の核医学診断装置。
4. 前記第１の発光確率モデル及び前記第２の発光確率モデルは、検出事象に含まれる複数の光子の内の所定の数の光子を検出する確率密度に関するモデルである、請求項１に記載の核医学診断装置。
5. 前記第１の光子数情報及び前記第２の光子数情報は、所定のキャリブレーションスキャンを行うことで特定される、請求項１に記載の核医学診断装置。
6. 前記第１の光子数情報及び前記第２の光子数情報は、前記キャリブレーションスキャンを行うことで前記第１の検出器及び前記第２の検出器の各々で取得されるエネルギースペクトルに基づいて特定される、請求項５に記載の核医学診断装置。
7. 前記第１の発光確率モデル及び前記第２の発光確率モデルは、前記第１の検出器及び前記第２の検出器の各々における発光の分布関数モデルに基づく、請求項１に記載の核医学診断装置。
8. 前記第１の発光確率モデル及び前記第２の発光確率モデルは、検出事象に対応する光子数の期待値に基づく、請求項１に記載の核医学診断装置。
9. 前記第１の光子数情報及び前記第２の光子数情報は、各々シンチレーション光子数を含む、請求項１に記載の核医学診断装置。
10. 前記第１の光子数情報及び前記第２の光子数情報は、各々チェレンコフ光子数を更に含む、請求項８に記載の核医学診断装置。
11. 前記算出部は、シンチレータ内での発光位置の不定性を考慮して、前記第１の発光確率モデルおよび前記第２の発光確率モデルを算出する、請求項１に記載の核医学診断装置。
12. 前記算出部は、前記シンチレータ内での位置ごとに発光確率密度を減弱係数に基づいて算出し、発光の分布関数及び前記発光確率密度に基づいて、前記第１の発光確率モデルおよび前記第２の発光確率モデルを算出する、請求項１１に記載の核医学診断装置。
13. 前記取得部は、前記第１の検出器で検出される第１のイベントに関する前記第１の光子数情報と、前記第２の検出器で検出される第２のイベントに関する前記第２の光子数情報とを取得し、
    前記特定部は、前記第１の光子数情報と前記第２の光子数情報とに基づいて、前記第１のイベントと前記第２のイベントとに基づいて規定されるコインシデンスに関する前記確率分布モデルを特定する、請求項１に記載の核医学診断装置。
14. 前記コインシデンスに関する確率分布モデルをＴｏＦ（Ｔｉｍｅ Ｏｆ Ｆｌｉｇｈｔ）カーネルとして利用することで、前記第１の検出器と前記第２の検出器との間のＬＯＲ（Ｌｉｎｅ Ｏｆ Ｒｅｓｐｏｎｓｅ）に対応するデータの再構成を行う再構成処理部を更に備える、請求項１３に記載の核医学診断装置。
15. 前記特定部は、イベントごとに特定された前記確率分布モデルを加算して、前記第１の検出器と前記第２の検出器との間のＴｏＦカーネルを生成する、請求項１３に記載の核医学診断装置。
16. 前記算出部は、前記第１の光子数情報に基づいて、前記第１のイベントについての第１の発光確率モデルを生成し、前記第２の光子数情報に基づいて、前記第２のイベントについての第２の発光確率モデルを生成し、
    前記特定部は、前記第１の発光確率モデルと前記第２の発光確率モデルとに基づいて、前記確率分布モデルを特定する、請求項１３に記載の核医学診断装置。
17. 前記算出部は、イベントごとに定められた発光分布関数モデルに基づいて、前記第１の発光確率モデルおよび前記第２の発光確率モデルを生成する、請求項１６に記載の核医学診断装置。
18. 前記取得部により計測される前記第１の光子数情報及び前記第２の光子数情報は、イベントごとに取得されるものである、請求項１３に記載の核医学診断装置。
19. 第１の検出器で検出される第１の光子数情報と、前記第１の検出器とは異なる第２の検出器で検出される第２の光子数情報とを取得し、
    前記第１の光子数情報に基づいて前記第１の検出器に対応する第１の発光確率モデルと、前記第２の光子数情報に基づいて前記第２の検出器に対応する第２の発光確率モデルとを算出し、
    前記第１の発光確率モデルと前記第２の発光確率モデルとに基づいて、前記第１の検出器と前記第２の検出器とによって規定されるＬＯＲ（Ｌｉｎｅ Ｏｆ Ｒｅｓｐｏｎｓｅ）の検出時間差に関する確率分布モデルを特定する、データ処理方法。
20. 第１の検出器で検出される第１の光子数情報と、前記第１の検出器とは異なる第２の検出器で検出される第２の光子数情報とを取得し、
    前記第１の光子数情報に基づいて前記第１の検出器に対応する第１の発光確率モデルと、前記第２の光子数情報に基づいて前記第２の検出器に対応する第２の発光確率モデルとを算出し、
    前記第１の発光確率モデルと前記第２の発光確率モデルとに基づいて、前記第１の検出器と前記第２の検出器とによって規定されるＬＯＲの検出時間差に関する確率分布モデルを特定する処理をコンピュータに実行させる、プログラム。
21. 少なくとも第１の検出素子アレイと第２の検出素子アレイとが設けられたモノリシック検出器と、
    前記第１の検出素子アレイにて検出される第１の光子数情報と、前記第２の検出素子アレイにて検出される第２の光子数情報とを取得する取得部と、
    前記第１の光子数情報に基づいて前記第１の検出素子アレイに対応する第１の発光確率モデルと、前記第２の光子数情報に基づいて前記第２の検出素子アレイに対応する第２の発光確率モデルとを算出する算出部と、
    前記第１の発光確率モデルと前記第２の発光確率モデルとに基づいて、前記第１の検出素子アレイと前記第２の検出素子アレイとの組み合わせに応じた検出時間差に関する確率分布モデルを特定する特定部と、を備える、核医学診断装置。
22. 少なくとも第１の検出素子アレイと第２の検出素子アレイとが設けられたモノリシック検出器の前記第１の検出素子アレイにて検出される第１の光子数情報と、前記第２の検出素子アレイにて検出される第２の光子数情報とを取得し、
    前記第１の光子数情報に基づいて前記第１の検出素子アレイに対応する第１の発光確率モデルと、前記第２の光子数情報に基づいて前記第２の検出素子アレイに対応する第２の発光確率モデルとを算出し、
    前記第１の発光確率モデルと前記第２の発光確率モデルとに基づいて、前記第１の検出素子アレイと前記第２の検出素子アレイとの組み合わせに応じた検出時間差に関する確率分布モデルを特定する、データ処理方法。
23. 少なくとも第１の検出素子アレイと第２の検出素子アレイとが設けられたモノリシック検出器の前記第１の検出素子アレイにて検出される第１の光子数情報と、前記第２の検出素子アレイにて検出される第２の光子数情報とを取得し、
    前記第１の光子数情報に基づいて前記第１の検出素子アレイに対応する第１の発光確率モデルと、前記第２の光子数情報に基づいて前記第２の検出素子アレイに対応する第２の発光確率モデルとを算出し、
    前記第１の発光確率モデルと前記第２の発光確率モデルとに基づいて、前記第１の検出素子アレイと前記第２の検出素子アレイとの組み合わせに応じた検出時間差に関する確率分布モデルを特定する処理をコンピュータに実行させる、プログラム。

Images