特許 6238565 平滑化方法及び平滑化装置

(19)【発行国】日本国特許庁(JP)

(12)【公報種別】特許公報(B2)

(11)【特許番号】6238565

(24)【登録日】2017年11月10日

(45)【発行日】2017年11月29日

(54)【発明の名称】平滑化方法及び平滑化装置

(51)【国際特許分類】

   G01T 1/161 20060101AFI20171120BHJP

【ＦＩ】

   G01T1/161 A

   G01T1/161 C

【請求項の数】16

【全頁数】16

(21)【出願番号】特願2013-105112(P2013-105112)

(22)【出願日】2013年5月17日

(65)【公開番号】特開2013-242313(P2013-242313A)

(43)【公開日】2013年12月5日

【審査請求日】2016年2月12日

(31)【優先権主張番号】13/475,584

(32)【優先日】2012年5月18日

(33)【優先権主張国】US

(73)【特許権者】

【識別番号】594164542

【氏名又は名称】東芝メディカルシステムズ株式会社

(74)【代理人】

【識別番号】110001771

【氏名又は名称】特許業務法人虎ノ門知的財産事務所

(72)【発明者】

【氏名】シャオフォン ニュウ

(72)【発明者】

【氏名】ユン ドン

(72)【発明者】

【氏名】ウェンリー・ワン

【審査官】 原 俊文

(56)【参考文献】

【文献】 特開昭６１−１３２８８７（ＪＰ，Ａ）

【文献】 特表２０１１−５０７６０３（ＪＰ，Ａ）

【文献】 特表２０１２−５２１２４７（ＪＰ，Ａ）

【文献】 特表２０１２−５０６０９２（ＪＰ，Ａ）

【文献】 向井祐一朗，大池洋史，和田俊和，“ビジョンアルゴリズムとノルムの選定に関する一考察”，電子情報通信学会技術研究報告，２０１１年 ５月３０日，第１１１巻，第７７号，ｐ．７３−３８

【文献】 河村圭，石井大祐，渡辺裕，“全変動最小化の高速計算手法”，電子情報通信学会論文誌，２０１０年 ３月 １日，第Ｊ９３−Ｄ巻，第３号，ｐ．３２６−３３５

【文献】 Hongbin Guo, Rosemary A. Renaut, Kewei Chen, Eric Reiman，"FDG-PET Parametric Imaging by Total Variation Minimization"，Computerized Medical Imaging and Graphics，２００９年 ６月，Vol. 33, No. 4，p. 295-303

【文献】 Xinhui Duan, Li Zhang, Yuxiang Xing, Zhiqiang Chen, Jianping Cheng，"Few-View Projection Reconstruction With an Iterative Reconstruction-Reprojection Algorithm and TV Constraint"，IEEE Transactions on Nuclear Science，２００９年 ６月，Vol. 56, No. 3，p. 1377-1382

(58)【調査した分野】（Int.Cl.，ＤＢ名）

Ｇ０１Ｔ １／１６１−１／１６６

(57)【特許請求の範囲】

【請求項1】

ポジトロン放射断層撮影（Positron Emission Tomography：ＰＥＴ）のＰＥＴスキャナで得られるデータから、遅延同時計数により偶発イベントデータを取得する取得ステップと、
Ｌ２ノルムを最小化する基準を用いた反復処理により、前記偶発イベントデータを平滑化処理する平滑化ステップと、
を含み、
前記取得ステップは、再構成断面上を示す座標を（ｓ、φ）とする、前記反復処理の最初の反復時の時間ｔ＝０に対応する初期偶発イベントデータであるｕ（ｓ、φ、ｔ＝０）を、ｕ_０（ｓ、φ）として取得し、
前記平滑化ステップは、
ｓ及びφに対する２次中心差分ｕ_ｓｓ及びｕ_φφを算出する２次中心差分算出ステップと、
一定のパラメータであるλにより定義される勾配ｕ_ｔを、ｕ_ｔ＝２（ｕ_ｓｓ＋ｕ_φφ）−λ（ｕ−ｕ_０）を用いて算出する勾配算出ステップと、
勾配の降下方向をパラメータ化するための理論的な時間である第１の時間ｔ_１と、第２の時間ｔ_２との時間差であり、正数であるΔｔをΔｔ＝ｔ_２−ｔ_１とし、最初の反復時のｔ_１を０とした場合、ｕ（ｓ、φ、ｔ_２）＝ｕ（ｓ、φ、ｔ_１）＋Δｔｕ_ｔを用いて、偶発イベントデータを更新する更新ステップと、
前記２次中心差分算出ステップと、前記勾配算出ステップと、前記更新ステップとを、前回の反復からのｕ（ｓ、φ、ｔ）の変化が、所定の閾値より小さくなるまで、反復する反復ステップと、
を含む、平滑化方法。

【請求項2】

ポジトロン放射断層撮影（Positron Emission Tomography：ＰＥＴ）のＰＥＴスキャナで得られるデータから、遅延同時計数により偶発イベントデータを取得する取得ステップと、
Ｌ２ノルムを最小化する基準を用いた反復処理により、前記偶発イベントデータを平滑化処理する平滑化ステップと、
を含み、
前記取得ステップは、再構成空間におけるアキシャル面上の位置を示す座標を（ｓ、φ）とし、再構成空間における軸方向の位置を示す座標を（ｚ、θ）とする、前記反復処理の最初の反復時の時間ｔ＝０に対応する初期偶発イベントデータであるｕ（ｓ、φ、ｚ、θ、ｔ＝０）を、ｕ_０（ｓ、φ、ｚ、θ）として取得し、
前記平滑化ステップは、
ｓ、φ、ｚ及びθに対する２次中心差分ｕ_ｓｓ、ｕ_φφ、ｕ_ｚｚ及びｕ_θθを算出する２次中心差分算出ステップと、
一定のパラメータであるλにより定義される勾配ｕ_ｔを、ｕ_ｔ＝２（ｕ_ｓｓ＋ｕ_φφ＋ｕ_ｚｚ＋ｕ_θθ）−λ（ｕ−ｕ_０）を用いて算出する勾配算出ステップと、
勾配の降下方向をパラメータ化するための理論的な時間である第１の時間ｔ_１と、第２の時間ｔ_２との時間差であり、正数であるΔｔをΔｔ＝ｔ_２−ｔ_１とし、最初の反復時のｔ_１を０とした場合、ｕ（ｓ、φ、ｚ、θ、ｔ_２）＝ｕ（ｓ、φ、ｚ、θ、ｔ_１）＋Δｔｕ_ｔを用いて、偶発イベントデータを更新する更新ステップと、
前記２次中心差分算出ステップと、前記勾配算出ステップと、前記更新ステップとを、前回の反復からのｕ（ｓ、φ、ｚ、θ、ｔ）の変化が、所定の閾値より小さくなるまで、反復する反復ステップと、
を含む、平滑化方法。

【請求項3】

λは、０．０５と２．０との間の値である、請求項１又は２に記載の平滑化方法。

【請求項4】

λは、０と０．０１との間の値である、請求項１又は２に記載の平滑化方法。

【請求項5】

Δｔは、０．０５と０．１との間の値である、請求項１又は２に記載の平滑化方法。

【請求項6】

ｕ（ｓ、φ、ｔ）は、以下の式（１）により得られる値が、前記所定の閾値より小さい時、収束する、請求項１に記載の平滑化方法。

【数1】

【請求項7】

ｕ（ｓ、φ、ｚ、θ、ｔ）は、以下の式（２）により得られる値が、前記所定の閾値より小さい時、収束する、請求項２に記載の平滑化方法。

【数2】

【請求項8】

前記初期偶発イベントデータが、前記ＰＥＴスキャナで測定可能な全ての応答線が含まれる領域に限定されるように、マスクを適用するマスクステップ、
を更に含む、請求項２に記載の平滑化方法。

【請求項9】

ポジトロン放射断層撮影（Positron Emission Tomography：ＰＥＴ）のＰＥＴスキャナで得られるデータから、遅延同時計数により偶発イベントデータを取得する取得処理と、
Ｌ２ノルムを最小化する基準を用いた反復処理により、前記偶発イベントデータを平滑化する平滑化処理と、を実行する処理回路、
を備え、
前記処理回路は、
前記取得処理において、再構成断面上を示す座標を（ｓ、φ）とする、前記反復処理の最初の反復時の時間ｔ＝０に対応する初期偶発イベントデータであるｕ（ｓ、φ、ｔ＝０）を、ｕ_０（ｓ、φ）として取得し、
前記平滑化処理において、
ｓ及びφに対する２次中心差分ｕ_ｓｓ及びｕ_φφを算出する２次中心差分算出処理と、
一定のパラメータであるλにより定義される勾配ｕ_ｔを、ｕ_ｔ＝２（ｕ_ｓｓ＋ｕ_φφ）−λ（ｕ−ｕ_０）を用いて算出する勾配算出処理と、
勾配の降下方向をパラメータ化するための理論的な時間である第１の時間ｔ_１と、第２の時間ｔ_２との時間差であり、正数であるΔｔをΔｔ＝ｔ_２−ｔ_１とし、最初の反復時のｔ_１を０とした場合、ｕ（ｓ、φ、ｔ_２）＝ｕ（ｓ、φ、ｔ_１）＋Δｔｕ_ｔを用いて、偶発イベントデータを更新する更新処理と、
前記２次中心差分算出処理と、前記勾配算出処理と、前記更新処理とを、前回の反復からのｕ（ｓ、φ、ｔ）の変化が、所定の閾値より小さくなるまで、反復する反復処理と
を実行する、平滑化装置。

【請求項10】

ポジトロン放射断層撮影（Positron Emission Tomography：ＰＥＴ）のＰＥＴスキャナで得られるデータから、遅延同時計数により偶発イベントデータを取得する取得処理と、
Ｌ２ノルムを最小化する基準を用いた反復処理により、前記偶発イベントデータを平滑化する平滑化処理と、を実行する処理回路、
を備え、
前記処理回路は、
前記取得処理において、再構成空間におけるアキシャル面上の位置を示す座標を（ｓ、φ）とし、再構成空間における軸方向の位置を示す座標を（ｚ、θ）とする、前記反復処理の最初の反復時の時間ｔ＝０に対応する初期偶発イベントデータであるｕ（ｓ、φ、ｚ、θ、ｔ＝０）を、ｕ_０（ｓ、φ、ｚ、θ）として取得し、
前記平滑化処理において、
ｓ、φ、ｚ及びθに対する２次中心差分ｕ_ｓｓ、ｕ_φφ、ｕ_ｚｚ及びｕ_θθを算出する２次中心差分算出処理と、
一定のパラメータであるλにより定義される勾配ｕ_ｔを、ｕ_ｔ＝２（ｕ_ｓｓ＋ｕ_φφ＋ｕ_ｚｚ＋ｕ_θθ）−λ（ｕ−ｕ_０）を用いて算出する勾配算出処理と、
勾配の降下方向をパラメータ化するための理論的な時間である第１の時間ｔ_１と、第２の時間ｔ_２との時間差であり、正数であるΔｔをΔｔ＝ｔ_２−ｔ_１とし、最初の反復時のｔ_１を０とした場合、ｕ（ｓ、φ、ｚ、θ、ｔ_２）＝ｕ（ｓ、φ、ｚ、θ、ｔ_１）＋Δｔｕ_ｔを用いて、偶発イベントデータを更新する更新処理と、
前記２次中心差分算出処理と、前記勾配算出処理と、前記更新処理とを、前回の反復からのｕ（ｓ、φ、ｚ、θ、ｔ）の変化が、所定の閾値より小さくなるまで、反復する反復処理と
を実行する、平滑化装置。

【請求項11】

λは、０．０５と２．０との間の値である、請求項９又は１０に記載の平滑化装置。

【請求項12】

λは、０と０．０１との間の値である、請求項９又は１０に記載の平滑化装置。

【請求項13】

Δｔは、０．０５と０．１との間の値である、請求項９又は１０に記載の平滑化装置。

【請求項14】

ｕ（ｓ、φ、ｔ）は、以下の式（３）により得られる値が、前記所定の閾値より小さい時、収束する、請求項９に記載の平滑化装置。

【数3】

【請求項15】

ｕ（ｓ、φ、ｚ、θ、ｔ）は、以下の式（４）により得られる値が、前記所定の閾値より小さい時、収束する、請求項１０に記載の平滑化装置。

【数4】

【請求項16】

前記処理回路は、前記初期偶発イベントデータが、前記ＰＥＴスキャナで測定可能な全ての応答線が含まれる領域に限定されるように、マスクを適用するマスク処理、
を更に実行する、請求項１０に記載の平滑化装置。

【発明の詳細な説明】

【技術分野】

【0001】

本発明の実施形態は、平滑化方法及び平滑化装置に関する。

【背景技術】

【0002】

ポジトロン放射断層撮影（Positron Emission Tomography：ＰＥＴ）の再構成において、偶発推定（random estimation）のために用いられる最も正確で一般的な方法は、遅延同時計数ウィンドウ推定（delayed coincidence window estimation）である。この推定では、１つの対消滅（annihilation）で発生する１対の光子間（一対のガンマ線間）の相互関係は、完全に排除され、その結果、偶発イベント（random event）のみを残すことができる。すなわち、ＰＥＴにおいて、同時計数用のタイムウィンドウを用いて生成される同時計数情報には、偶発的に同時計数された情報である偶発同時計数情報が含まれる。この偶発同時計数情報の補正では、一般的に、遅延同時計数を基礎とする方法が用いられている。遅延同時計数とは、一方のガンマ線検出イベントの情報（計数情報）に含まれる時間をずらして、タイムウィンドウ内に属する２つのガンマ線検出イベントを、擬似的に同時計数する方法である。遅延同時計数により同時計数された遅延同時計数情報には、対消滅ガンマ線の２つの計数情報を組み合わせた真の同時計数情報が含まれない。このため、遅延同時計数情報を用いて、同時計数情報に含まれる偶発同時計数情報を見積もることができる。ＰＥＴスキャナから得られるデータ（個々のイベントの計数情報）から遅延同時計数により取得されるデータは、偶発イベントデータとなる。しかし、限定された同時計数タイムウィンドウと、短い収集時間のために、統計学的には、測定された遅延データは、真の偶発データの分布の１つの場面（不連続で離散的な分布）しか具現化していない。従って、適切に処理されない場合、遅発同時計数ウィンドウデータ（遅延同時計数情報により得られるサイノグラムの画像や再構成画像のデータ）を減算した後のプロンプトデータ（同時計数用のタイムウィンドウを用いて生成される同時計数情報により得られるサイノグラムの画像や再構成画像のデータ）の中には、差異（すなわち、ノイズ）が増加する。すなわち、遅延同時計数情報は、離散的で有限数のデータであるため、統計ノイズが大きい。従来の遅延同時計数では、偶発同時計数の連続分布を見積もるために、平滑化処理が必要となる。

【0003】

高周波ノイズを削除するために、フーリエ解析のように、均一にサンプリングされるデータ空間で展開される画像平滑化処理技術は多く知られている。しかし、直接取得される生の遅延サイノグラムは、均一にサンプリングされる空間にない。生の遅延サイノグラムを、均一にサンプリングされる補間されたサイノグラムに変換できるならば、様々な標準的な平滑化処理技術を使用できるであろう。均一に補間されたサイノグラム空間での偶発サイノグラムの平滑化処理の後、不規則にサンプリングされている生のサイノグラム空間へ逆補間することができるであろう。

【0004】

ノイズは、ノイズを有する偶発サイノグラムのスペクトルの高周波成分を表わす。偶発サイノグラムの平滑化処理は、フーリエ領域の高周波をフィルタ処理することによって、適用可能になる。すなわち、偶発サイノグラムの平滑化処理は、フーリエ変換した周波数空間で、高周波除去することで実行できる。しかし、平滑化処理は、また、偶発サイノグラム分布に存在する信号に大きな変化も引き起こす。例えば、平滑化した結果は、高周波除去に対して鋭敏で、変化が大きい。更に、フーリエ方式は、信号のグローバルな表現であり、補間された偶発サイノグラムが、測定可能な全ての応答線（Line-Of-Response：ＬＯＲ、同時計数線とも呼ばれる）を含むマスクの内側に限定される偶発サイノグラムの平滑化処理の場合のように、不規則なマスクを用いた画像に対しては失敗する。マスクは、画像における測定可能部位を表す。従って、マスクの外側の任意の部位は、測定不可であり、ピクセル強度はゼロである。すなわち、フーリエ変換は、グローバル変換であるため、再構成用の有効視野（Field-Of-View：ＦＯＶ）等の領域に対して、マスクを適用している場合、利用できない。

【先行技術文献】

【特許文献】

【0005】

【特許文献1】特開２００９−４２０２９号公報

【発明の概要】

【発明が解決しようとする課題】

【0006】

本発明が解決しようとする課題は、偶発同時計数の連続分布を推定することができる平滑化方法及び平滑化装置を提供することである。

【課題を解決するための手段】

【0007】

実施形態の平滑化方法は、取得ステップと、平滑化ステップとを含む。取得ステップは、ポジトロン放射断層撮影（Positron Emission Tomography：ＰＥＴ）のＰＥＴスキャナで得られるデータから、遅延同時計数により偶発イベントデータを取得する。平滑化ステップは、Ｌ２ノルムを最小化する基準を用いた反復処理により、前記偶発イベントデータを平滑化処理する。取得ステップは、再構成断面上を示す座標を（ｓ、φ）とする、前記反復処理の最初の反復時の時間ｔ＝０に対応する初期偶発イベントデータであるｕ（ｓ、φ、ｔ＝０）を、ｕ_０（ｓ、φ）として取得し、平滑化ステップは、ｓ及びφに対する２次中心差分ｕ_ｓｓ及びｕ_φφを算出する２次中心差分算出ステップと、一定のパラメータであるλにより定義される勾配ｕ_ｔを、ｕ_ｔ＝２（ｕ_ｓｓ＋ｕ_φφ）−λ（ｕ−ｕ_０）を用いて算出する勾配算出ステップと、勾配の降下方向をパラメータ化するための理論的な時間である第１の時間ｔ_１と、第２の時間ｔ_２との時間差であり、正数であるΔｔをΔｔ＝ｔ_２−ｔ_１とし、最初の反復時のｔ_１を０とした場合、ｕ（ｓ、φ、ｔ_２）＝ｕ（ｓ、φ、ｔ_１）＋Δｔｕ_ｔを用いて、偶発イベントデータを更新する更新ステップと、前記２次中心差分算出ステップと、前記勾配算出ステップと、前記更新ステップとを、前回の反復からのｕ（ｓ、φ、ｔ）の変化が、所定の閾値より小さくなるまで、反復する反復ステップと、を含む。

【図面の簡単な説明】

【0008】

【図1】図１は、偶発差異の削減前後の補間された偶発サイノグラムを示す図である。

【図2】図２は、図１に示す偶発差異の削減前後の偶発サイノグラムのφ＝６４におけるプロファイルを示す図である。

【図3】図３は、本開示の一実施形態に係る方法のフローチャートである。

【図4】図４は、本開示の別の一実施形態に係る方法のフローチャートである。

【図5】図５は、本開示の実施形態を実行するコンピュータシステムを示す図である。

【発明を実施するための形態】

【0009】

以下、添付図面を参照して、平滑化方法及び平滑化装置の実施形態を詳細に説明する。本明細書に記載する実施形態とこれに付随する多くの利点は、以下の詳細な説明を参照し添付図面と関連付けて考えれば、より完全に理解できる。ただし、添付図面は、一例として提供され、本開示を限定するものではない。

【0010】

一実施形態では、ポジトロン放射断層撮影（Positron Emission Tomography：ＰＥＴ）のＰＥＴスキャナから得られる偶発イベントデータを平滑化処理するための方法（平滑化方法）及び装置（平滑化装置）が提供される。一実施形態に係る方法は、ＰＥＴスキャナで得られるデータから、遅延同時計数により偶発イベントデータを取得する取得ステップと、Ｌ２ノルムを最小化する基準を用いた反復処理により、前記偶発イベントデータを平滑化処理する平滑化ステップと、を含む。２次元の平滑化処理を行なう場合の実施形態に係る方法における取得ステップは、ｔ＝０に対応する初期偶発イベントデータであるｕ（ｓ、φ、ｔ＝０）を、ｕ_０（ｓ、φ）として取得する。そして、２次元の平滑化処理を行なう場合の実施形態に係る方法における平滑ステップは、ｓ及びφに対する２次中心差分ｕ_ｓｓ及びｕ_φφを算出する２次中心差分算出ステップと、一定のパラメータであるλにより定義される勾配ｕ_ｔを、ｕ_ｔ＝２（ｕ_ｓｓ＋ｕ_φφ）−λ（ｕ−ｕ_０）を用いて算出する勾配算出ステップと、正数であるΔｔをΔｔ＝ｔ_２−ｔ_１とし、最初の反復時のｔ_１を０とした場合、ｕ（ｓ、φ、ｔ_２）＝ｕ（ｓ、φ、ｔ_１）＋Δｔｕ_ｔを用いて、偶発イベントデータを更新する更新ステップと、「２次中心差分算出ステップと、勾配算出ステップと、更新ステップ」とを、前回の反復からのｕ（ｓ、φ、ｔ）の変化が、所定の閾値より小さくなるまで、反復する反復ステップと、を含む。

【0011】

更に、上記の方法の実施形態では、λは、０．０５と２．０との間の値である。或いは、上記の方法の実施形態では、λは、０と０．０１との間の値である。また、更に、上記の方法では、Δｔは、０．０５と０．１との間の値である。また、上記の方法では、ｕ（ｓ、φ、ｔ）は、以下の式（１）により得られる値が、前記所定の閾値より小さい時、収束する。

【0012】

【数1】

【0013】

また、３次元の平滑化処理を行なう場合の別の実施形態に係る方法における取得ステップは、ｔ＝０に対応する初期偶発イベントデータであるｕ（ｓ、φ、ｚ、θ、ｔ＝０）を、ｕ_０（ｓ、φ、ｚ、θ）として取得する。そして、３次元の平滑化処理を行なう場合の別の実施形態に係る方法における平滑ステップは、ｓ、φ、ｚ及びθに対する２次中心差分ｕ_ｓｓ、ｕ_φφ、ｕ_ｚｚ及びｕ_θθを算出する２次中心差分算出ステップと、一定のパラメータであるλにより定義される勾配ｕ_ｔを、ｕ_ｔ＝２（ｕ_ｓｓ＋ｕ_φφ＋ｕ_ｚｚ＋ｕ_θθ）−λ（ｕ−ｕ_０）を用いて算出する勾配算出ステップと、正数であるΔｔをΔｔ＝ｔ_２−ｔ_１とし、最初の反復時のｔ_１を０とした場合、ｕ（ｓ、φ、ｚ、θ、ｔ_２）＝ｕ（ｓ、φ、ｚ、θ、ｔ_１）＋Δｔｕ_ｔを用いて、偶発イベントデータを更新する更新ステップと、「２次中心差分算出ステップと、勾配算出ステップと、更新ステップ」とを、前回の反復からのｕ（ｓ、φ、ｚ、θ、ｔ）の変化が、所定の閾値より小さくなるまで、反復する反復ステップと、を含む。

【0014】

更に、上記の方法の実施形態では、λは、０．０５と２．０との間の値である。或いは、上記の方法の実施形態では、λは、０と０．０１との間の値である。また、更に、上記の方法では、Δｔは、０．０５と０．１との間の値である。また、上記の方法では、ｕ（ｓ、φ、ｚ、θ、ｔ）は、以下の式（２）により得られる値が、前記所定の閾値より小さい時、収束する。

【0015】

【数2】

【0016】

更に、３次元の平滑化処理を行なう場合の別の実施形態に係る方法では、初期偶発イベントデータが、ＰＥＴスキャナで測定可能な全ての応答線が含まれる領域に限定されるように、マスクを適用するマスクステップ、を更に含む。

【0017】

本開示の一実施形態では、Ｌ２ノルムのサイノグラム画像の勾配（gradient）を反復的に最小化するための偶発サイノグラム平滑化処理方法が記載される。サイノグラム画像の勾配は、サイノグラム空間のピクセルごとの局所的な勾配の総和として計算される。

【0018】

２次元（２Ｄ）の平滑化処理では、ｕ_０（ｓ、φ）は、ノイズを有する画像の各ピクセルの強度値を示す。ここで、「ｓ及びφ」は、ノイズを有する画像の全てのピクセルを示すΩ（積分範囲）に属する。すなわち、「ｓ及びφ」は、再構成断面（アキシャル面）の位置を示す座標再構成断面（アキシャル面）上を示す座標となる。また、２次元（２Ｄ）の平滑化処理では、ｕ（ｓ、φ）は、真のノイズの無い画像を示す。イメージングモデルは、「ｕ_０（ｓ，φ）＝ｕ（ｓ，φ）＋ｎ（ｓ，φ）」として記述される。

【0019】

ただし、ｕ（ｓ、φ）は、加法性雑音（例えば、ガウス分布）であり、Ｅ［ｕ_０（ｓ、φ）］＝ｕ（ｓ、φ）である。

【0020】

目的は、画像勾配のＬ２ノルムを最小化する基準の下で、ｕ_０（ｓ、φ）からｕ（ｓ、φ）を推定することである。下記の導出を明確にするために、「ｕ及びｕ_０」を、「ｕ（ｓ、φ）及びｕ_０（ｓ、φ）」を示す記号とすると、最適化問題は、以下の式（３）となる。

【0021】

【数3】

【0022】

従って、最小化される目的関数は、以下の式（４）となる。

【0023】

【数4】

【0024】

オイラーラグランジュ方程式（Euler-Lagrange equation）を使用することにより、以下の式（５）が得られる。

【0025】

【数5】

【0026】

ここで、変量「ｔ」は、勾配の降下方向をパラメータ化するための理論的な時間として導入される。「ｕ_ｔ」は、ｔの方向に沿った勾配を表わす。「ｕ_０」は、「ｔ＝０」における進化画像（evolving image）を表わし、オリジナルのノイズを有する画像「ｕ（ｓ、φ、０）＝ｕ_０（ｓ、φ）」である。上記の「ｕ_ｔ＝２（ｕ_ｓｓ＋ｕ_φφ）−λ（ｕ−ｕ_０）」では、「λ」は、平滑化画像が、オリジナルのノイズを有する画像と大きく異なる場合、異なる量の制約を加えるのに用いられる制約パラメータである。強いノイズの場合（例えば、ノイズが信号又は画像を圧倒する場合）、λは「０」、又は、λは非常に小さな値（例えば、略０．０１より小さく、且つ、０より大きい値）に設定される。弱いノイズの場合、λは、約０．５と２との間に設定される。すなわち、λの値が大きいほど、偶発同時計数の離散的な分布に追従した平滑化が行なわれ、λの値が小さいほど、偶発同時計数の分布を優先的に連続分布にする平滑化が行なわれる。

【0027】

次に、画像は、以下の式（６）で更新される。

【0028】

【数6】

【0029】

式（６）では、Δｔ＝ｔ_２−ｔ_１は、安定性の限界を表わし、略０．０５と０．１との間の小さな値に設定される。「ｔ」が増加するにつれて、ノイズ除去画像が得られる。

【0030】

図３は、本開示の一実施形態に係る方法のフローチャートである。図３は、ＰＥＴスキャナから得られるイベントデータの偶発差異を削減する方法のフローチャートである。

【0031】

ステップ１００では、ｔ＝０で、初期偶発イベントデータである最初の偶発イベントデータ「ｕ（ｓ、φ、ｔ＝０）＝ｕ_０（ｓ、φ）」を取得する。初期偶発イベントデータは、遅延同時計数により得られた遅延同時計数情報から生成された偶発サイノグラムの画像、又は、偶発サイノグラムから再構成された再構成画像である。

【0032】

ステップ１１０では、ｓ、φに対して２次中心差分ｕ_ｓｓ、ｕ_φφを計算する。具体的には、２次中心差分（２回微分）は、以下の式（７）により算出される。なお、式（７）では、２つのパラメータを「ｘ及びｙ」で示している。

【0033】

【数7】

【0034】

次に、ステップ１２０では、勾配ｕ_ｔを、「ｕ_ｔ＝２（ｕ_ｓｓ＋ｕ_φφ）−λ（ｕ−ｕ_０）」を用いて計算する。ただし、λは、一定のパラメータである。

【0035】

ステップ１３０では、「ｕ（ｓ、φ、ｔ_２）＝ｕ（ｓ、φ、ｔ_１）＋Δｔｕ_ｔ」を用いて、偶発イベントデータを更新する。ただし、Δｔ＝ｔ_２−ｔ_１であり、最初の反復時は「ｔ_１＝０」、Δｔは、正である。

【0036】

ステップ１４０では、ｕ（ｓ、φ、ｔ）は、収束したか否かが判定される。すなわち、ステップＳ１４０では、適切なノルムを使ってｕ（ｓ、φ）の変化が、閾値より小さいか否かが、決定される。ｕ（ｓ、φ、ｔ）が収束した場合、処理は終了する。しかし、ｕ（ｓ、φ、ｔ）が収束しなかった場合、ステップ１１０〜１３０が、反復される。ｕ（ｓ、φ、ｔ）の収束は、上記の数（１）により算出される値が、所定の閾値より小さいか否かを調べることによって、決定することができる。

【0037】

図１は、偶発差異の削減前後の補間された偶発サイノグラムを示す図であり、図２は、図１に示す偶発差異の削減前後の偶発サイノグラムのφ＝６４におけるプロファイルを示す図である。図１及び２は、テストデータで上記の図３に示す方法を実行した結果の一例を示す図である。図１は、１枚の（ｓ、φ）スライスを用いて、偶発差異の削減の前後に補間された偶発サイノグラムを示す図である。図１の上図は、平滑化処理前の偶発サイノグラムを示し、図１の下図は、平滑化処理後の偶発サイノグラムを示す。図１の上図に示す平滑化処理前の偶発サイノグラムは、複数の点が描出されており、偶発同時計数が離散的な分布となっているのに対し、図３に示す平滑化処理により、図１の下図に示す平滑化処理後の偶発サイノグラムは、偶発同時計数が連続的で滑らかな分布として描出されている。図１は、図３に示す平滑化処理が、勾配の曲率にペナルティを課した局所的な平滑化処理であることを示している。また、図２は、偶発差異の削減前後における図１のφ＝６４のプロファイルを示している。図２からも、離散的な偶発同時計数の分布が、平滑化処理による推定結果から、連続的で滑らかな偶発同時計数の分布となることが示されている。すなわち、図３に示す方法によって、偶発同時計数の連続分布を推定することができる。

【0038】

上記の方法を、３次元再構成を行なう場合、すなわち、パラメータが４つの４次元（４Ｄ）平滑化処理に適用した場合、４Ｄデータ「ｕ（ｓ、φ、ｚ、θ）」に対する最適化問題は、以下の式（８）となる。なお、「ｓ及びφ」は、再構成空間におけるアキシャル面上の位置を示す座標であり、「ｚ及びθ」は、再構成空間における軸方向（アキシャル方向）の位置を示す座標である。

【0039】

【数8】

【0040】

ただし、式（８）のΩは、マスクにより定義された不規則な形状の部位を表す。マスクにより定義される部位は、矩形の立方体に限定されるものではない。

【0041】

画像は、式「ｕ（ｓ、φ、ｚ、θ、ｔ_２）＝ｕ（ｓ、φ、ｚ、θ、ｔ_１）＋（ｔ_２−ｔ_１）ｕ_ｔ」により更新される。ただし、「ｕ_ｔ＝２（ｕ_ｓｓ＋ｕ_φφ＋ｕ_ｚｚ＋ｕ_θθ）−λ（ｕ−ｕ_０）」である。

【0042】

図４は、本開示の別の一実施形態に係る方法のフローチャートである。図４は、ＰＥＴスキャナから得られるイベントデータの偶発差異を削減する方法のフローチャートである。

【0043】

ステップ２００では、ｔ＝０で、初期偶発イベントデータである最初の偶発イベントデータ「ｕ（ｓ、φ、ｚ、θ、ｔ＝０））＝ｕ_０（ｓ、φ、ｚ、θ）」を取得する。初期偶発イベントデータは、ＰＥＴスキャナが収集した３次元空間での遅延同時計数により得られた遅延同時計数情報から生成された偶発サイノグラムの画像、又は、偶発サイノグラムから再構成された再構成画像である。

【0044】

ステップ２１０では、ｓ、φ、ｚ、θに対して２次中心差分ｕ_ｓｓ、ｕ_φφ、ｕ_ｚｚ及びｕ_θθを計算する。具体的には、２次中心差分は、以下の式（９）により算出される。

【0045】

【数9】

【0046】

次に、ステップ２２０では、勾配ｕ_ｔを、「ｕ_ｔ＝２（ｕ_ｓｓ＋ｕ_φφ＋ｕ_ｚｚ＋ｕ_θθ）−λ（ｕ−ｕ_０）」を用いて計算する。ただし、λは、一定のパラメータである。

【0047】

ステップ２３０では、偶発イベントデータが、「ｕ（ｓ、φ、ｚ、θ、ｔ_２）＝ｕ（ｓ、φ、ｚ、θ、ｔ_１）＋Δｔｕ_ｔ」を用いて、偶発イベントデータを更新する。ただし、Δｔ＝ｔ_２−ｔ_１であり、最初の反復時は「ｔ_１＝０」、Δｔは、正である。

【0048】

ステップ２４０では、ｕ（ｓ、φ、ｚ、θ、ｔ）は、収束したか否かが判定される。すなわち、ステップＳ２４０では、適切なノルムを使ってｕ（ｓ、φ、ｚ、θ、ｔ）の変化が、閾値より小さいか否かが、決定される。ｕ（ｓ、φ、ｚ、θ、ｔ）が収束した場合、処理は終了する。しかし、ｕ（ｓ、φ、ｚ、θ、ｔ）が収束しなかった場合、ステップ２１０〜２３０が、反復される。ｕ（ｓ、φ、ｚ、θ、ｔ）の収束は、上記の数（２）により算出される値が、所定の閾値より小さいか否かを調べることによって、決定できる。方法は、取得した最初の偶発イベントデータが、ＰＥＴスキャナの測定可能な応答線全てを含む部位に限定されるように、マスクを適用するマスクステップを更に含むことができる。図４に示す方法によっても、偶発同時計数の連続分布を推定することができる。

【0049】

上述した２次偏導関数（second-order partial derivaties）であるｕ_ｓｓ、ｕ_φφ、ｕ_ｚｚ及びｕ_θθは、マスクによって限定されるボクセルで計算される。ボクセルが、マスクの境界線である場合、偏導関数は、境界線を横切る方向でゼロに設定される。

【0050】

上記の実施形態の利点は、マスクを適用する方法を、柔軟に適用可能であることである。すなわち、上記の実施形態では、再構成用のＦＯＶに任意のマスクを適用したデータに対しても、安定した平滑化処理を行なうことができる。図３及び図４に示す処理の後は、プロンプトデータから平滑化処理後の偶発イベントデータを差し引くことで、ノイズの無いＰＥＴ画像が再構成される。なお、図３及び図４に示す本実施形態の平滑化処理は、所謂、汎関数の変分法により停留値を求める処理と類似した処理である。

【0051】

上述した平滑化方法は、コンピュータシステム又はプログラム可能な論理を使って実行可能である。図５は、本開示の実施形態を実行できるコンピュータシステム１２０１を示す図である。コンピュータシステム１２０１は、例えば、上述の処理を実行する処理回路、モジュール、又は、ユニットを含むことができる。

【0052】

コンピュータシステム１２０１は、バス１２０２に連結されたディスクコントローラ１２０６を含む。ディスクコントローラ１２０６は、磁気のハードディスク１２０７、及び、取り外し可能なメディアドライブ１２０８（例えば、フロッピー（登録商標）ディスクドライブ、読取専用コンパクトディスクドライブ、読込み／書込みコンパクトディスクドライブ、コンパクトディスクジュークボックス、テープドライブ、及び、取り外し可能な光磁気ドライブ）のように、情報及び指令を格納するための１つ又は複数の格納デバイスを制御する。格納デバイスは、適切なデバイスインターフェースを使って、コンピュータシステム１２０１に加えることができる。適切なデバイスインターフェースは、例えば、ＳＣＳＩ（Small Computer System Interface）、ＩＤＥ（Integrated Device Electronics）、ＥＩＤＥ（Enhanced IDE）、ＤＭＡ（Direct Memory Access）、又は、ｕｌｔｒａ―ＤＭＡである。

【0053】

また、コンピュータシステム１２０１は、特定目的論理デバイス（例えば、ＡＳＩＣ（Application Specific Integrated Circuit））、又は、特定目的論理デバイス（例えば、ＳＰＬＤ（Simple Programmable Logic Device）、ＣＰＬＤ（Complex Programmable Logic Device）、ＦＰＧＡ（Field-Programmable Gate Array））を含むこともできる。

【0054】

また、コンピュータシステム１２０１は、バス１２０２に連結された表示部コントローラ１２０９を含むことができる。表示部コントローラ１２０９は、情報をコンピュータユーザーに表示するためのディスプレイ、又は、液晶ディスプレイ（Liquid Crystal Display：ＬＣＤ）等の表示部１２１０を制御する。コンピュータシステムは、コンピュータユーザーと対話したり、処理装置１２０３へ情報を提供したりするために、キーボード１２１１及びポインティングデバイス１２１２等の入力デバイスを含む。ポインティングデバイス１２１２は、例えば、マウス、トラックボール、タッチスクリーンセンサー用のフィンガー、又は、処理装置１２０３に指令情報及び命令選択を伝えたり、表示部１２１０上のカーソル運動を制御したりするためのポインティングスティックであっても良い。

【0055】

コンピュータシステム１２０１は、メインメモリ１２０４等のメモリに含まれる１つ又は複数の指令の１つ又は複数のシーケンスを実行する処理装置１２０３に応答して、本開示の処理ステップの一部又は全てを実行する。そのような指令は、ハードディスク１２０７、又は、取り外し可能なメディアドライブ１２０８等の別のコンピュータ可読媒体からメインメモリ１２０４に読み込むことができる。メインメモリ１２０４に含まれる指令のシーケンスを実行するために、多重処理配列（multi-processing arrangement）内の１つ又は複数のプロセッサを採用しても良い。別の実施形態では、配線による回路をソフトウェア指令の代わりに、又は、それとの組み合わせに使用してもよい。このように、実施形態は、ハードウェア回路とソフトウェアとの如何なる特有の組み合わせにも制限されない。

【0056】

上記のように、コンピュータシステム１２０１は、本開示の教示によってプログラムされた指令を格納し、且つ、本明細書に記載のデータ構造、表、記録、又は他のデータを含む少なくとも１つのコンピュータ可読媒体又はメモリを含む。コンピュータ可読媒体は、コンパクトディスク、ハードディスク、フロッピー（登録商標）ディスク、テープ、光磁気ディスク、「ＥＰＲＯＭ（Electrically Programmable Read Only Memory）、ＥＥＰＲＯＭ（Electrically Erasable Programmable Read Only Memory）、フラッシュＥＰＲＯＭ」等のＰＲＯＭ（Programmable Read Only Memory）、ＤＲＡＭ（Dynamic Random Access Memory）、ＳＲＡＭ（Static RAM）、ＳＤＲＡＭ（Synchronous DRAM）、又は、任意の他の磁気媒体、コンパクトディスク（例えば、ＣＤ−ＲＯＭ）、又は、任意の他の光学媒体、パンチカード、紙テープ、又は、穴のパターンをもつ他の物理的媒体である。

【0057】

本開示は、コンピュータ可読媒体のうちのいずれか１つ又はそれらの組み合わせに格納された、コンピュータシステム１２０１を制御するためのソフトウェアを含む。該ソフトウェアは、本開示を実行するためのデバイス、又は、デバイス群を駆動し、コンピュータシステム１２０１に人間のユーザーと対話できるようにするためのものである。そのようなソフトウェアには、デバイスドライバ、オペレーティングシステム、及びアプリケーションソフトが含まれるが、これらに限定されるものではない。そのようなコンピュータ可読媒体は、本開示を実施する際に実行される処理工程の全て又は一部（処理工程が区分される場合）を実行するための本開示のコンピュータプログラム製品を更に含む。

【0058】

本実施形態のコンピュータ符号化デバイスは、如何なる機械言語に翻訳処理可能な、又は、実行可能な符号化機構であってもよく、スクリプト、機械言語に翻訳可能なプログラム、ＤＬＬ（Dynamic Link Libraries）、Ｊａｖａ（登録商標）クラス、及び、完全な実行可能プログラムを含むが、これらに限定されるものではない。更に、本実施形態の処理工程の部分は、より良い性能、信頼性、及び費用の全て、又は、いずれかのために、区分けされてもよい。

【0059】

本明細書で使用される「コンピュータ可読媒体」は、処理装置１２０３に実行のための指令を提供する際に関係する任意の非一過性媒体を指す。コンピュータ可読媒体は、不揮発性媒体又は揮発性媒体を含む多くの形状をとることができるが、これらに限定されるものではない。不揮発性媒体は、例えば、ハードディスク１２０７、又は、取り外し可能なメディアドライブ１２０８等の光学ディスク、磁気ディスク、及び光磁気ディスクを含む。揮発性媒体は、メインメモリ１２０４等の動的メモリを含む。一方、伝送媒体は、バス１２０２を作るワイヤーを含む、同軸ケーブル、銅線、及び光ファイバーを含む。また、伝送媒体は、電波や赤外線のデータ通信中に生成される波等の音波又は光波の形をとっても良い。

【0060】

実施のために、処理装置１２０３への１つ又は複数の指令の１つ又は複数のシーケンスを実行する際に、様々な形式のコンピュータ可読媒体を含むことができる。例えば、指令は、最初、リモートコンピュータの磁気ディスクで処理されるものとする。リモートコンピュータは、本開示の全て又は一部を実行するための指令を、遠隔で動的メモリに読み込み、モデムを使って電話線を通して指令を送信できる。コンピュータシステム１２０１へのモデムローカルは、電話線でデータを受信し、データをバス１２０２に乗せることができる。バス１２０２は、データをメインメモリ１２０４へ伝達する。処理装置１２０３は、メインメモリ１２０４から指令を取り出して実行する。メインメモリ１２０４が受信した指令は、処理装置１２０３による実行の前か後かに、ハードディスク１２０７又は取り外し可能なメディアドライブ１２０８に任意に格納されることになる。

【0061】

また、コンピュータシステム１２０１は、バス１２０２に連結された通信インターフェース１２１３も含む。通信インターフェース１２１３は、ネットワークリンク１２１４に連結された双方向データ通信を備える。ネットワークリンク１２１４は、例えば、ローカルエリアネットワーク（Local Area Network：ＬＡＮ）１２１５、又は、インターネット等の別の通信ネットワーク１２１６に接続される。例えば、通信インターフェース１２１３は、任意のパケット交換ＬＡＮに取り付けられるネットワークインターフェイスカードであってもよい。別の例として、通信インターフェース１２１３は、ＩＳＤＮ（Integrated Services Digital Network）カードであっても良い。無線リンクが、実装されてもよい。そのような実装例では、通信インターフェース１２１３は、様々な種類の情報を表わすデジタルデータストリームを運ぶ電気、電磁気、又は光の信号を送受信する。

【0062】

ネットワークリンク１２１４は、一般的に、１つ又は複数のネットワークを通して別のデータデバイスへのデータ通信を提供する。例えば、ネットワークリンク１２１４は、ローカルエリアネットワーク１２１５を介して、又は、通信ネットワーク１２１６を介して、通信サービスを提供するサービスプロバイダによって運営される設備を通して別のコンピュータへの接続を提供できる。ローカルエリアネットワーク１２１５及び通信ネットワーク１２１６は、例えば、デジタルデータストリームを運ぶ電気、電磁気、又は光の信号と、関連する物理層（例えば、ＣＡＴ５ケーブル、同軸ケーブル、光ファイバー等）とを使用する。デジタルデータをコンピュータシステム１２０１に搬入搬出する、様々なネットワークを通る信号、ネットワークリンク１２１４上の信号、及び、通信インターフェース１２１３を通る信号は、ベースバンド信号、又は、搬送波ベース信号で実装されても良い。ベースバンド信号は、デジタルデータビットストリームを記述する非変調電気パルスとしてデジタルデータを伝達する。ここで、「ビット」という語は、シンボルを意味すると広く解釈されており、各シンボルは、少なくとも１つ又は複数の情報ビットを運ぶものである。また、デジタルデータは、導電媒体一面に伝播される、又は、伝播媒体の中を通る電磁波のように伝達される、振幅、位相、及び周波数の全て又はいずれかの偏移変調信号をともなうような搬送波を変調するために使用されてもよい。従って、デジタルデータは、有線通信路を通して非変調のベースバンドデータとして伝送、及び搬送波を変調することによって、ベースバンドとは異なる所定の周波数バンド内で伝送の両方又は一方が行われても良い。コンピュータシステム１２０１は、プログラムコードを含むデータを（単数又は複数の）ネットワークであるＬＡＮ１２１５及び通信ネットワーク１２１６、ネットワークリンク１２１４、及び通信インターフェース１２１３を経由して送受信できる。さらに、ネットワークリンク１２１４は、ＬＡＮ１２１５を経由して、ＰＤＡ（Personal Digital Assistant）ラップトップコンピュータ、又は、携帯電話等のモバイルデバイス１２１７への接続を提供できる。

【0063】

以上、説明したとおり本実施形態によれば、偶発同時計数の連続分布を推定することができる。

【0064】

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

【符号の説明】

【0065】

１２０１ コンピュータシステム
１２０２ バス
１２０３ 処理装置
１２０４ メインメモリ
１２０５ 読み出し専用メモリ（Read Only Memory：ＲＯＭ）
１２０６ ディスクコントローラ
１２０７ ハードディスク
１２０８ 取り外し可能なメディアドライブ
１２０９ 表示部コントローラ
１２１０ 表示部
１２１１ キーボード
１２１２ ポインティングデバイス
１２１３ 通信インターフェース
１２１４ ネットワークリンク
１２１５ ローカルエリアネットワーク（Local Area Network：ＬＡＮ）
１２１６ 通信ネットワーク
１２１７ モバイルデバイス

【図1】

【図2】

【図3】

【図4】

【図5】