JP7352382B2 - 画像処理装置、画像処理方法及びプログラム - Google Patents

Description

本発明は、画像処理装置、画像処理方法及びプログラムに関する。

近年、放射線のエネルギー情報を利用した画像化技術であるスペクトラル・イメージング技術では複数の異なる放射線エネルギー画像から複数の物質の厚さ画像を求めたり、あるいは面密度と実効原子番号画像の画像を求めたりすることができる。こうした厚さ画像、面密度画像、実効原子番号画像(「物質にかかわる平面分布」の一例)は物質透過後の放射線が厚さなどの指数関数として減衰すると仮定して、透過後の放射線の画素値から逆算することで求めることができる。

特許文献1では軟部組織の画像を平滑化し、その画像を蓄積画像から減算することで、骨部画像の画質を改善する技術が開示されている。

特開平３－２８５４７５号公報

しかしながら、物質にかかわる平面分布の推定演算においては、システムノイズや放射線に由来する量子ノイズが計算過程で伝播し、出力された画像はノイズを含んだものとなり得る。特許文献１の技術においても、単純な四則演算を用いて平面分布の推定演算をしており、放射線の線量が低くなる場合などにおいては、ノイズを除去しきれず、また、平滑化により画像が不鮮明になる場合が生じ得る。

本発明は、上記の課題に鑑み、確率的に分布するノイズを低減し、より優れた画質の画像を取得することが可能な画像処理技術を提供することを目的とする。

本発明の目的を達成するために、本発明の一態様による画像処理装置は、異なる放射線エネルギーで撮影した複数の放射線画像を入力画像とし、前記入力画像に含まれる物質について平面分布を示す出力画像を取得する平面分布取得手段を備え、前記平面分布取得手段は、尤度と、注目画素の画素値または周囲の画素の画素値を入力とする事前確率から、入力が前記入力画像であるときに、出力が前記平面分布の値となる確率が最大となるように前記出力画像を決定する。

本発明によれば、ノイズを低減し、より優れた画質の画像を取得することが可能になる。

添付図面は明細書に含まれ、その一部を構成し、本発明の実施形態を示し、その記述と共に本発明の原理を説明するために用いられる。
実施形態に係る放射線撮影装置の構成例を示す図。 実施形態に係る撮影処理の流れを説明する図。 物質の平面分布を取得する処理の流れを説明する図。 画素値の差分と事前確率の分布を模式的に示す図。

以下、添付図面を参照して本発明の実施形態を詳しく説明する。尚、以下の実施形態は特許請求の範囲に係る本発明を限定するものでなく、また本実施形態で説明されている特徴の組み合わせの全てが本発明の解決手段に必須のものとは限らない。以下の実施形態及び特許請求の範囲において、放射線は、Ｘ線の他、α線、β線、γ線、及び各種粒子線なども含む。

図１は、本発明の実施形態に係る放射線撮影装置１００の構成例を示す図である。放射線撮影装置１００は、放射線管球１０１、ＦＰＤ１０３（Flat Panel Detector（放射線検知部））、コンピュータ１０９（画像処理装置）、操作盤１０４及び同期装置１１８、照射指示部１０８を有する。

放射線管球１０１は、同期装置１１８から送信された放射線照射の指示信号に基づいて、放射線を発生させ、被写体１０２に放射線を照射する。放射線管球１０１には、低エネルギーの放射線を除去するためのフィルタを設けることがある。放射線管球１０１から放射線が被写体１０２に照射されると、放射線は被写体１０２を透過してＦＰＤ１０３に入射する。

ＦＰＤ１０３は、放射線を可視光に変換する蛍光体と、可視光を電荷、電圧に変換する撮像素子（素子は複数個あり、1つの素子と1つの画素が対応する）と、変換された電圧をさらにデジタル値の画像へ変換する画像形成部と、画像に補正（オフセット補正、ゲイン補正、欠陥補正）を行う画像補正部を有する。撮像素子は、入射光に応じた信号を出力する画素がアレイ状（二次元の領域）に配置されている。各画素の画像形成部は蛍光体により可視光に変換された放射線を電気信号に変換し、画像信号として出力する。尚、放射線を電圧に直接変換する方式のＦＰＤ、例えば、テルル化カドミウム(ＣｄＴｅ)や非晶質セレン(ａ-Ｓｅ)を用いたＦＰＤにおいても本実施形態を適用することが可能である。ＦＰＤ１０３は画像信号に基づく電荷の蓄積を行って放射線画像を取得し、放射線画像をコンピュータ１０９に出力する。

操作盤１０４は、放射線管球１０１から照射する放射線の照射条件を指定することが可能であり、操作盤１０４は、条件指定部１０５、撮影部位指定部１０６、条件表示部１０７を有する。条件指定部１０５は、ボタン等の操作部を介した操作者からの操作入力に基づいて、管電圧、管電流、照射時間を指定可能である。操作盤１０４は画質調整を行うための画質調整部(ダイヤルなど)を保持しており、操作者はダイヤルを操作することで医師などの画像確認者が放射線画像を確認しやすい画質に設定することができる。

撮影部位指定部１０６には、撮影部位と放射線の照射条件が予め関連付けられていて、操作者は、被写体の撮影部位を指定すると、撮影部位指定部１０６は、指定されたさ撮影部位に関連付けられている放射線の照射条件を呼び出し、設定することができるようになっている。条件表示部１０７には、指定されている放射線の照射条件が表示されるように構成されている。

照射指示部１０８はフットペダルで、操作者はペダルを踏むことで、照射指示信号が同期装置１１８に入力され、同期装置１１８の制御に基づいて放射線管球１０１から放射線を照射することが可能になる。

ＦＰＤ１０３から出力された画像（放射線画像）は、コンピュータ１０９へ転送される。コンピュータ１０９は、ＦＰＤ１０３から取得した放射線画像を処理する画像処理装置として機能することが可能であり、機能構成として、ＦＰＤ１０３から画像（放射線画像）を取得する画像取得部１１０、事前情報取得部１１１、平面分布取得部１１２、事前情報保存部１１３を有している。これらの機能構成は、例えば、図示されていないコンピュータ１０９の中央演算処理装置（ＣＰＵ）、主記憶装置、ハードディスクなどの記憶装置、高速計算のためのグラフィック・プロセッシング・ユニット、ＬＡＮ(Local Area Network)アダプタなどのハードウェアにより構成することが可能である。

また、各機能はソフトウェアとして実装することが可能である。例えば、中央演算処理装置と、記憶装置から読み込んだプログラムを用いて、各部の機能を構成することが可能である。同様の機能を果たすのであれば、それらは集積回路（同等の電気回路）などで構成してもよい。尚、ＦＰＤ１０３に搭載されている画像補正部を、コンピュータ１０９に設けてもよい。また、操作盤１０４の代わりに、コンピュータ１０９に予め放射線の照射条件を保存しておくことも可能である。この場合、操作盤１０４で放射線の照射条件を指定する際に、操作盤１０４が保存されている照射条件をコンピュータ１０９から読み出して、設定すればよい。

画像取得部１１０は、ＦＰＤ１０３から画像（放射線画像）を取得する。また、事前情報取得部１１１は、事前情報として、例えば、放射線に関する確率分布の情報（例えばポアソン分布、正規分布の関数やポアソン分布や正規分布のパラメタである平均値、標準偏差などの情報）、管球情報(管電圧、管電流、照射時間、放射線フィルタの有無を含む情報)、ＦＰＤ１０３に関する情報(蓄積時間など)、撮影部位情報(例えば、被写体の胸部、頭部などの部位)、被写体情報(入力装置１１４で入力される被写体の識別番号、氏名、性別、年齢、身長、体重など)、造影剤情報(造影剤の種類、濃度など)、Ｃアーム１１６の情報(材質、位置、角度など)、架台１１７の情報(架台の位置など)が含まれる。

平面分布取得部１１２は、異なる放射線エネルギーで撮影した複数の放射線画像を入力画像とし、前記入力画像に含まれる物質について平面分布を示す出力画像を取得する。平面分布取得部１１２による計算プロセスは後に説明する。

コンピュータ１０９に接続する入力装置１１４は、コンピュータ１０９やＦＰＤ１０３に対する指示を入力するユーザーインターフェイスとして機能する。入力装置１１４は、例えば、コンピュータ１０９のキーボードなどにより構成することが可能である。操作者は、入力装置１１４を使用して、撮影部位情報や被写体情報を入力することが可能である。撮影部位情報や被写体情報は、入力装置１１４による入力の他に、例えば、コンピュータ１０９に実装されるＬＡＮアダプタを通じ、ネットワークを介して外部の装置から入力してもよい。

また、コンピュータ１０９に接続する表示装置１１５は、コンピュータ１０９がＦＰＤ１０３から受信した放射線画像（デジタル画像）やコンピュータ１０９から出力された画像（画像処理された画像）を表示することが可能である。表示装置１１５として、特にカラー・ディスプレイが用いられた場合は、物質にかかわる平面分布を色として表現することが可能になり、画像診断性能の向上を図ることが可能になる。

また、操作盤１０４の機能の一部を入力装置１１４、表示装置１１５が備えてもよい。例えば、操作盤１０４の条件指定部１０５あるいは撮影部位指定部１０６を通じて行った放射線の照射条件の入力操作を入力装置１１４にて行ってもよいし、入力された管電圧などの放射線の照射条件の表示を、条件表示部１０７の代わりに表示装置１１５で行ってもよい。

同期装置１１８は、放射線管球１０１、ＦＰＤ１０３、操作盤１０４、照射指示部１０８及びコンピュータ１０９に接続され、ＦＰＤ１０３の状態と、照射指示部１０８の押下状態と、コンピュータ１０９の処理状態から放射線の照射可否を決定する。同期装置１１８は照射制御部として機能し、決定したタイミングに基づいて放射線の照射及び撮影を制御する。

本発明の実施形態に係る撮影処理の流れを図２のフローチャートを参照して説明する。放射線撮影装置１００の操作者は、ステップＳ２０１において、放射線の照射条件を入力する。放射線の照射条件入力は、例えば、操作盤１０４の条件指定部１０５あるいは撮影部位指定部１０６を通じて行うことが可能である。また、入力装置１１４を通じて、被写体の識別番号などの情報も併せて入力する。

ステップＳ２０２において、操作者は、照射指示部１０８を通じて放射線の照射を指示する。その後、ステップＳ２０３において、同期装置１１８は、放射線管球１０１、ＦＰＤ１０３、操作盤１０４、照射指示部１０８、コンピュータ１０９の状態から照射可否を判断する。照射不可の場合（Ｓ２０３－Ｎｏ）、ステップＳ２０４において、同期装置１１８は表示装置１１５に警告を表示させる。一方、ステップＳ２０３において、同期装置１１８が照射可と判定した場合（Ｓ２０３－Ｙｅｓ）、処理をステップＳ２０５に進める。

ステップＳ２０５において、事前情報取得部１１１は、事前情報として、放射線に関する確率分布の情報、管球情報、ＦＰＤ１０３に関する情報、撮影部位情報、被写体情報等を取得する。この時、識別物質数などの情報に応じて、画像取得枚数(例えば、エネルギーレベルの異なる２つのエネルギー（高エネルギー、低エネルギー）による物質分離を行うなら2枚)を決定する。そして、事前情報取得部１１１は、決定した画像取得枚数やＦＰＤ１０３における情報(蓄積時間など)等を同期装置１１８に通知する。

ステップＳ２０６において、同期装置１１８は、事前情報取得部１１１から取得した情報に基づいて、放射線照射の指示信号を生成し、放射線管球１０１とＦＰＤ１０３に放射線照射の指示信号を送信する。これにより、放射線管球１０１による放射線の照射とＦＰＤ１０３による画像取得が制御される。ＦＰＤ１０３は、異なる放射線エネルギーレベルで撮影された複数の放射線画像（例えば、低エネルギー放射線画像、高エネルギー放射線画像）を撮影する。

ＦＰＤ１０３により必要な画像取得枚数が取得されたら、ステップＳ２０７において、撮影を終了する。ステップＳ２０７において、撮影が終了すると、画像取得部１１０は、ＦＰＤ１０３から画像（放射線画像）を取得する。すなわち、画像取得部１１０は、異なる放射線エネルギーレベルで撮影された複数の放射線画像（低エネルギー放射線画像、高エネルギー放射線画像）を取得する。

ステップＳ２０８において、平面分布取得部１１２は、ステップＳ２０５で求められた放射線に関する確率分布の情報等の事前情報と、ステップＳ２０６にて取得された画像とに基づいて、物質にかかわる平面分布を算出する。平面分布の詳細な算出処理については後述する。

ステップＳ２０９において、コンピュータ１０９のＣＰＵは平面分布の算出結果をコンピュータ１０９内の記憶装置し、また、ＣＰＵは、平面分布の算出結果を表示装置１１５に表示させる。ＣＰＵは、例えば、分離した物質の厚さ、実効原子番号、面密度などの平面分布の算出結果と可視光の色とを対応させた画像（物質識別画像）を表示装置１１５にカラーで表示するような表示制御を行うことが可能である。

（平面分布算出処理）
次に、ステップＳ２０８の物質にかかわる平面分布の算出処理方法を詳細に説明する。ここでは、物質にかかわる平面分布算出の一例として、２色(高エネルギー放射線画像と低エネルギー放射線画像の２種類)の画像から、物質ごとの厚さの分布を示す画像（物質分離画像）を出力するプロセスについて説明するが、画像が３色以上の場合においても同様に平面分布の算出処理を適用することができる。

ＦＰＤ１０３内に一時保存された放射線分布情報は、サンプリングホールド実施後、読み出し可能となり、画像取得部１１０は、ＦＰＤ１０３から異なるタイミングで、放射線分布情報（Ｘ_Ｌ）と放射線分布情報（Ｘ_Ｌ＋Ｘ_Ｈ）の読み出しを実施する。画像取得部１１０は、放射線分布情報（Ｘ_Ｌ＋Ｘ_Ｈ）から放射線分布情報（Ｘ_Ｌ）を差し引くことで、放射線分布情報（Ｘ_Ｈ）を得ることができる。ここで、低エネルギーの放射線分布情報（Ｘ_Ｌ）が低エネルギーの放射線による画像（低エネルギー放射線画像）になり、高エネルギーの放射線分布情報（Ｘ_Ｈ）が高エネルギーの放射線による画像（高エネルギー放射線画像）になる。

平面分布取得部１１２は、異なる放射線エネルギーで撮影した複数の放射線画像に含まれる物質についての物質特性画像を取得する。平面分布取得部１１２は、１回の放射線照射で取得された複数の放射線画像を用いて、物質分離画像、物質識別画像といった物質特性画像を取得することが可能である。

ここで物質分離画像とは、複数の放射線画像（低エネルギー放射線画像、高エネルギー放射線画像）に含まれる複数の物質を分離した画像であり、被写体を特定の２以上の物質で表した場合に、その物質の厚さで構成された２以上の画像に分離した画像をいう。

また、物質識別画像とは、複数の放射線画像（低エネルギー放射線画像、高エネルギー放射線画像）に含まれる物質の実効原子番号と面密度の分布を示す画像であり、被写体を特定の１物質で表した場合に、その物質の実効原子番号と面密度の分布に分解した画像をいう。

また、物質にかかわる平面分布は、分離した物質の厚さの分布を示す画像（厚さ画像）、特定の物質に関する実効原子番号の分布を示す画像（実効原子番号画像）と面密度の分布を示す画像（面密度）の２種類に分ける場合や、実効原子番号画像と面密度画像と厚さ画像の３種類以上の画像に分ける場合においても同様に平面分布の算出処理を適用することができる。

入力画像として、異なる放射線エネルギーの画像として、２つの放射線画像（低エネルギー放射線画像、高エネルギー放射線画像）をＩ=(Ｉ(Ｅ₁、Ｉ(Ｅ₂））とする。ここで、Ｅ₁、Ｅ₂は放射線のエネルギーを表す。また、出力画像として、物質の厚さ画像をｆ（ｄ）=(ｄ₁、ｄ₂）とする。ノイズのない場合の入力画像（Ｉ）と出力画像（ｆ（ｄ））との関係は、以下の数１式により表すことができる。

実際に得られる画像において、放射線の統計的性質に由来する量子ノイズや、検出器の電気的な雑音に由来するシステムノイズが存在する場合には、入力画像（Ｉ）と出力画像（ｆ（ｄ））との関係において、数１式の関係は必ずしも成立しない場合が生じ得る。一般にノイズをｎとすれば、入力画像（Ｉ）、出力画像（ｆ（ｄ））とノイズの関係は、以下の数２式のように、ノイズｎは、入力画像（Ｉ）と出力画像（ｆ（ｄ））との差分として表される。

ノイズｎの従う確率分布は、物質の厚さがｄの時、画像がIとなる確率と考えることができる。よって、ノイズｎの従う確率分布はｐ（Ｉ｜ｄ）と表すことができる。

求めようとしている物質にかかわる平面分布は、入力画像がI（=(Ｉ(Ｅ₁、Ｉ(Ｅ₂））)の時、入力画像に含まれ得る各物質の厚さｄの平面分布である。各物質の厚さｄについて、求める厚さは入力画像がIであるとき、確率が最大となる物質の厚さｄであると考えることができる。つまり以下の数３式を満たす物質の厚さｄとなる。平面分布取得部１１２は、入力が入力画像（Ｉ）であるときに、出力が平面分布の値（ｄ）となる確率（事後確率）が最大となるように出力画像を決定する。

また、求める物質の厚さは確率分布の対数が最大となる物質の厚さｄと考えてもよく、その場合、式の形は以下の数４式、あるいは数５式で表すことができる。

数５式をベイズの定理を用いて変形すると、以下の数６式のように変形することができる。確率（事後確率）は、尤度と事前確率に基づいて取得することができる。すなわち、平面分布取得部１１２は、尤度と、注目画素の画素値または周囲の画素の画素値を入力とする事前確率とから、入力が入力画像（Ｉ）であるときに、出力が平面分布の値（ｄ）となる確率（事後確率）が最大となるように出力画像を決定する。ここで、尤度は、入力される値が平面分布の値（物質の厚さｄ）のときに、入力に対応する画像が入力画像となる確率を示すものである。

ここで、数６式において、ｌｏｇｐ（Ｉ｜ｄ）は、尤度と呼ばれる項であり、ｌｏｇｐ_prior（ｄ）は事前確率と呼ばれる項である。事前確率の項の形式を変更することで、多様な事前情報を反映した確率を計算することが可能になる。この事前確率は、数６式の解をより頑健にする効果を持ち、正則化項と呼ばれることもある。

数６式を求める方法はいくつか存在するが、極値が最小値の場合、例えば、最急降下法で求めることが可能である。この他、ニュートン法やその他一般的に知られている方法を用いることが可能である。最急降下法を適用する場合、以下の数７式の形でｄ_ｎ＋１の逐次近似解を求め、求めた逐次近似解から数６式の解を求めることができる。平面分布取得部１１２は、出力画像を反復演算による逐次近似解により取得することができる。

数７式において、ｄ_ｎはn回目の反復で求められた物質の厚さである。αは反復の速度を表す係数である。ｇｒａｄはベクトルの微分（勾配）を表す記号である。つまり、入力画像がIであるとき、確率が最大となる物質の厚さｄ_ｎ＋１は、数７式のように反復を繰り返すことで逐次的に求めることができる。

図３は、図２のステップＳ２０８の詳細な処理（物質の平面分布を取得する処理）の流れを説明するフローチャートである。ステップＳ３０１において、画像取得部１１０は、ＦＰＤ１０３から画像（放射線画像）を取得する。図２のステップＳ２０７において撮影が終了すると、画像取得部１１０は、異なる放射線エネルギーレベルで撮影された複数の放射線画像（低エネルギー放射線画像、高エネルギー放射線画像）を取得する。低エネルギー放射線画像、高エネルギー放射線画像の２つの放射線画像が、ステップＳ３０１における入力画像Ｉとなる。

ステップＳ３０２において、事前情報取得部１１１は、ステップＳ２０５で求められた放射線に関する確率分布の情報等の事前情報を取得する。

次にステップＳ３０３において、平面分布取得部１１２は、逐次的演算を行う前の初期値（ｄ_１）の設定を行う。設定する初期値は、例えば、系にノイズがないと仮定したときの平面分布（物質の厚さ）とすることができる。ノイズのない場合の入力画像（Ｉ）と物質の厚さｄに関する出力画像（ｆ（ｄ））との関係は、数１式のように表すことができるので、平面分布（物質の厚さ）の初期値ｄ_１は、以下の数８式から求めることができる。数８式に基づいて、平面分布取得部１１２は、ノイズがないと仮定したときの平面分布の値を初期値として反復演算を行う。

別の初期値の設定方法として、数８式で得られた画像に更に画像処理を加えた画像を物質の厚さの初期値として入力することも可能である。あるいは、エネルギーサブトラクション法により、異なる放射線エネルギーにより撮影された複数の放射線画像に含まれる複数の物質に分離した出力画像（例えば、軟部組織の厚さ画像や骨部の厚さ画像）について平滑化する手法で物質にかかわる平面分布を計算し、その出力を初期値として入力してもよい。すなわち、平面分布取得部１１２は、出力画像に対する平滑化処理により取得した平面分布の値を初期値として反復演算を行うことも可能である。

ステップＳ３０４以降の処理が反復演算（逐次近似）の処理となり、平面分布取得部１１２は、予め設定された繰り返し回数ｎに到達するまで（Ｓ３０４－Ｙｅｓ）、ステップＳ３０５において、数７式における逐次近似解の更新を繰り返す。逐次近似解の更新処理において、数７式の微分項(ｇｒａｄ)の処理は数値的に求めてもよいし、解析的に求めてもよい。微分項(ｇｒａｄ)の処理方法の詳細は後述する。

そして、ステップＳ３０６において、平面分布取得部１１２は、反復の速度を表す係数α_ｎを更新する。一般に繰り返しの最初の数回における解は最終的な解と比較して差分が大きいはずであるから、平面分布取得部１１２は、反復の速度を表す係数α_ｎを大きい値（第１係数値）に設定する。反復演算を行うにつれて誤差が小さくなるはずだから、平面分布取得部１１２は、反復の速度を表す係数α_ｎを小さい値、すなわち、第１係数値に比べて小さい第２係数値に設定する。反復の速度を表す係数α_ｎは、例えば、ｋ×α_ｎ→α_ｎ＋１という形式となる。ここで、「→」の表記は、ｋ×α_ｎをα_ｎ＋１に代入するという意味であり、係数ｋは0より大きく1より小さい実数の定数である（０＜ｋ（実数の定数）＜１）。第１係数値をα_ｎとすると、第２係数値（α_ｎ＋１）は、ｋ×第１係数値で求めることができる。平面分布取得部１１２は、反復演算（逐次近似）の繰り返し回数に応じて、反復の速度を表す係数α_ｎが小さくなるように係数α_ｎの設定を更新する。

ステップＳ３０４の判定処理で、設定された繰り返し回数ｎに到達した場合（Ｓ３０４－Ｎｏ）、平面分布取得部１１２は、本処理を終了する。

（具体的な数式を用いた処理の例）
（物質の厚さと画像（出力画像）の関係ｆ（ｄ）について）
次に、具体的な数式を用いた処理の例について説明する。まず、数１式、数２式などに現れる物質の厚さと画像（出力画像）の関係ｆ（ｄ）については、以下の数９式、あるいは数１０式に示した透過に関するランベルト・ベールの法則の式を用いることができる。数９式、数１０式において、Ｉ'／Ｉ_０は減弱率、μは線減弱係数、ｄは物質の厚さであり、Ｅは放射線のエネルギー、Φ（Ｅ）は放射線のスペクトルを示す。Ｉは入力画像であり、異なる２つの放射線エネルギーで撮影した画像（低エネルギー放射線画像、高エネルギー放射線画像）である。独立した２つの連立方程式を解くことで、未知の変数である物質の厚さｄ_１、ｄ_２の値を求めることができる。

（数６式の尤度の項について）
入力される値が平面分布の値（物質の厚さｄ）のときに、入力に対応する画像が入力画像となる確率を示す尤度の項について説明する。平面分布取得部１１２は、尤度の確率分布として、正規分布またはポアソン分布のいずれかを用いる。

物質の厚さｄの時、画像がIとなる確率を示す尤度（ｌｏｇｐ（Ｉ｜ｄ）：数６式）は、確率分布を正規分布で近似した以下の数１１式や、確率分布をポアソン分布で近似した以下の数１２式を用いることができる。ここで、数１１式において、σは分散を示している。また、数１２式において、所定の時間中に平均でλ_１回発生する事象がｋ_１回（ｋは正の整数）発生する確率、λ₂回発生する事象がｋ₂回（k は正の整数）発生する確率を示し、Ｉ_１＝ｋ_１となる確率と、Ｉ_２＝ｋ_２となる確率は独立である。

（数６式の事前確率の項について）
平面分布取得部１１２は、異なる放射線エネルギーで撮影した複数の放射線画像を入力画像とし、入力画像に含まれる物質について平面分布を示す出力画像を取得し、出力画像における注目画素の画素値と周囲の画素の画素値との差分により平面分布の連続性を評価する。

出力画像における連続性に関する事前確率は、平面分布の連続性の評価に基づいて取得される。数６式で示した事前確率の項（ｌｏｇｐ_prior（ｄ））については、例えば、注目画素と周囲の画素との比較（差分）を行う項を導入することができる。一例として、例えば、以下の数１３式を事前確率の項として用いることが可能である。ここで、平面分布取得部１１２は、事前確率を、出力画像における注目画素の画素値と注目画素に対する周囲の画素の画素値との差分に基づいて取得する。

数１３式において、「｜ ｜」の表記は絶対値を示す。数１３式における、ｋ、ｌは、注目画素に対して上下、左右、斜めの隣接画素の項を表現するための添え字であり、一例として、ｋ＝±１、ｌ＝±１の組み合わせが挙げられる。比較（差分）の範囲は、隣接画素に限られず、ｋ＝±１、ｌ＝±１の組み合わせの隣接画素を含む、より広い範囲にわたった周囲の画素との比較（差分）の項を、数１３式の事前確率の項に追加することも可能である。

数１３式で示したような事前確率の項(正則化項)を導入することで、平面分布取得部１１２は、出力画像（物質分離画像）における注目画素の画素値と周囲の画素の画素値との比較（差分）により、出力画像における平面分布の空間的な連続性を評価することができる。そして、平面分布取得部１１２は、平面分布の空間的な連続性の評価結果に基づいて、出力画像に対するノイズ低減処理を行う。これにより、確率的に分布するノイズが低減され、より優れた画質の画像を取得することが可能になる。

数１３式で示したような事前確率の項(正則化項)における正則係数μ_i、η_iは事前確率の項(正則化項)の影響の度合い（強さ）を表す係数である。正則係数μ_i、η_iが大きい値をとるほど事前確率で導入した仮定が強く影響することになる。図４は、画素値の差分と事前確率の分布を模式的に示す図である。注目画素４０１及び周囲の画素４０２のように画素値の差分が小さい場合、事前確率は大きくなる傾向を示す。一方、注目画素４０１及び周囲の画素４０３のように画素値の差分が大きい場合、事前確率は小さくなる傾向を示す。

すなわち、数１３式の事前確率の項(正則化項)は、連続性の評価値（注目画素の画素値と周囲の画素の画素値との差分）が小さければ小さいほど事前確率が大きくなることを意味する。

また、正則係数μ_i、η_iに大きい値を設定しても事前確率の項(正則化項)の出力を大きくすることができる。これにより、事前確率で導入した仮定（画素値の連続性により出力画像は滑らかであるという事前の仮定）が強く影響する。正則係数μ_i、η_iが大きい値であればあるほど、出力画像の取得演算処理（逐次近似演算）により注目画素と周囲の画素との画素値の差分が小さくなり、より滑らかな出力画像（物質分離画像）が得られるようになる。滑らかな出力画像が得られる傾向にあるから、事前確率の項(正則化項)を用いることで、物質にかかわる平面分布（物質の厚さの分布を示す物質分離画像）のノイズを低減することができる。

平面分布取得部１１２は、事前情報取得部１１１が取得した事前情報に含まれる撮影条件に基づいて、事前確率の項(正則化項)の影響の度合い（強さ）を表す正則係数μ_i、η_iを変更することができる。例えば、放射線の照射条件(線量)が所定の閾値に比べて小さい場合、平面分布取得部１１２は、正則係数μ_i、η_iを所定の基準値に比べて大きい値に変更して事前確率の項(正則化項)を大きくすることで、ノイズ除去（ノイズ低減）をより優先した設定にすることができる。

平面分布取得部１１２は、放射線の照射条件(線量)に関する情報を、事前情報取得部１１１が取得した事前情報に基づいて取得することができる。例えば、事前情報に含まれる管球情報(管電圧、管電流、照射時間)は線量に関わる情報であり、平面分布取得部１１２は、放射線の照射条件(線量)に基づいて、事前確率の項(正則化項)の影響の度合いを表す正則係数μ_i、η_iを変更することができる。尚、放射線撮影装置１００の入力装置１１４に画質調整を行うための調整装置が付属している場合は、平面分布取得部１１２は、調整装置が示す値に基づいて事前確率の項(正則化項)の正則係数μ_i、η_iの設定を調整することができる。

また、事前情報に含まれる被写体情報(身長、体重など)は被写体の体形に関わる情報であり、線量に関わる情報である。平面分布取得部１１２は、放射線の照射条件(線量)に関する情報として被写体情報（被写体の体形に関する情報）に基づいて、事前確率の項(正則化項)の影響の度合いを表す正則係数μ_i、η_iを変更することができる。例えば、被写体の体重が大きいほど透過する線量は低くなるので、平面分布取得部１１２は、正則係数μ_i、η_iを所定の基準値に比べて大きい値に変更して事前確率の項(正則化項)を大きくすることで、ノイズ除去（ノイズ低減）をより優先した設定にすることができる。

また、Ｃアーム１１６や、架台１１７の位置設定の情報も線量に関わる情報である。例えば、Ｃアーム１１６の位置や角度、すなわち、放射線発生装置（放射線管球１０１）と放射線検知部（ＦＰＤ１０３）との相対的な位置または角度によって、放射線が被写体１０２を通過する距離が長くなる場合には、距離が短い場合に比べて透過する線量は相対的に低くなることが予想される。また、放射線が被写体を通過する距離が短くなる場合は、距離が長い場合に比べて透過する線量は相対的に高くなることが予想される。

平面分布取得部１１２は、放射線の照射条件(線量)に関する情報として放射線発生装置（放射線管球１０１）と放射線検知部（ＦＰＤ１０３）との相対的な位置または角度に基づいて、事前確率の項(正則化項)の影響の度合いを表す正則係数μ_i、η_iを変更することができる。

また、架台１１７は放射線を遮る要素の一つであるため、架台１１７の有無も線量に関わる情報であり、平面分布取得部１１２は、放射線の照射条件(線量)に関する情報として架台１１７の情報に基づいて、事前確率の項(正則化項)の影響の度合いを表す正則係数μ_i、η_iを変更することができる。

発明者は、物質分離画像のうち、特に脂肪などの軟部組織に関する物質の厚さの分布を示す画像（厚さ画像）に対して、事前確率の項(正則化項)を導入することが、ノイズ低減等の画質向上に効果的であることを見出した。そのため、着目する部位が脂肪などの軟部組織である場合には、軟部組織に係る厚さ画像に対する事前確率の項(正則化項)の正則係数を、骨などの硬部組織に係る厚さ画像に対する事前確率の項(正則化項)の正則係数に比べて大きく設定して事前確率の項(正則化項)を大きくすることで、ノイズが除去（ノイズ低減）され、より優れた画質の画像（物質分離画像）を提供することが可能になる。

この場合、平面分布取得部１１２は、軟部組織の厚さ画像における、注目画素の画素値と周囲の画素の画素値との差分により軟部組織の厚さ画像における連続性の評価を行う。そして、平面分布取得部１１２は、平面分布の空間的な連続性の評価結果に基づいて、出力画像に対するノイズ低減処理を行う。これにより、確率的に分布するノイズが低減され、より優れた画質の画像を取得することが可能になる。

尚、着目する部位が、骨などの硬部組織である場合には、骨などの硬部組織に係る厚さ画像に対する事前確率の項(正則化項)の正則係数を、脂肪などの軟部組織に係る厚さ画像に対する事前確率の項(正則化項)の正則係数に比べて大きく設定してもよい。

数１３式を解く場合に最急降下法においては、事前確率の項(正則化項)の厚さに関する微分を計算することになるが、数１３式を微分した場合、数１３式のうち、以下の数１４式、数１５式に示す２つの項は不連続な解となる。

数１４式及び数１５式に示す２つの項は、具体的には以下の数１６式でそれぞれ表すことができる。

数１４式及び数１５式における不連続性を回避するためには、絶対値を含む項の代わりに微分が連続になる関数、例えば、シグモイド関数を含む以下の数１７式の関数の積分を、数１４式で示す項の代わりに用いることで不連続性を回避することができる。

以上、事前確率の項(正則化項)は注目画素と周囲の画素との差分の項（例えば、数１３式）を持つ例を説明したが、事前確率の項(正則化項)は必ずしも注目画素と周囲の画素の差分であるとは限らない。平面分布取得部１１２は、事前確率を、出力画像（物質の厚さ画像）における注目画素の画素値に基づいて取得することも可能である。例えば、物質の厚さ画像(ｄ₁やｄ₂)はゼロ以上の有限の値となるから、物質の厚さを入力として、以下の数１８式に示すような正則化項を導入することができる。数１８式において、κは正の定数で、正則係数の一例である。

物質の厚さがゼロより小さくなると、数１８式に示す事前確率の項(正則化項)であるｌｏｇ_prior(ｄ_1i,j、ｄ_2i,j)もゼロより小さくなり、ゆえに尤度と事前確率の項(正則化項)とを合わせた全体の事後確率は小さくなる。物質の厚さ画像(ｄ₁，ｄ₂)はゼロより大きい画素値を持つ画像となる可能性が向上する。この事前確率も、正負の値にばらつくはずだった厚さ画像を正の値だけにする効果を持つから、周囲の画素との差分の項を持つ事前確率の項(正則化項)と同様に平面分布の連続性を評価する項の一例であり、またノイズを低減する結果をもたらす。すなわち、平面分布取得部１１２は、出力画像（物質の厚さ画像）における注目画素の画素値により平面分布の連続性を評価し、評価結果に基づいて出力画像におけるノイズを低減する。これにより、確率的に分布するノイズが低減され、より優れた画質の画像を取得することが可能になる。

（処理の高速化）
数１式、数８式から数１３式を、数７式に代入することで、数７式の逐次近似解（更新解）を求めることができるが、式の変形により解析解をもとめることが可能となり、結果演算をより高速化することができる。

数１１式は、物質の厚さｄの時、画像がIとなる確率を示す尤度の確率分布を正規分布で近似した式であるが、数１１式の対数をとると、数１１式は以下の数１９式のように表すことができる。

尤度に関しては数１９式のように展開することができ、事前確率の項(正則化項)に関しては数１３式のように展開できるので、数７式中の－logp(I｜ｄ)－logｐ_prior（ｄ）という数式部分は、以下の数２０式で表すことができる。

数２０式を数７式に代入すると、ベクトルの微分（勾配）を表すｇｒａｄの項は、以下の数２１式のように表すことができる。

数２１式の解を解析的に求めることにより、演算精度を維持しつつ演算処理をより高速化することができる。

数１２式は尤度の確率分布をポアソン分布で近似した式であるが、数１２式の対数をとると、数１２式は数２２式のように表すことができる。

尤度に関しては数２２式のように展開することができ、事前確率の項(正則化項)に関しては数１３式のように展開できるので、数７式中の－logp(I｜ｄ)－logｐ_prior（ｄ）という数式部分は、以下の数２３式で表すことができる。

数２３式を数７式に代入すると、ベクトルの微分（勾配）を表すｇｒａｄの項は、以下の数２４式のように表すことができる。

数２４式の解を解析的に求めることにより、演算精度を維持しつつ演算処理をより高速化することができる。

（単色近似）
尤度の確率分布を正規分布で近似した数１１式に現れるｆ_i（ｄ）、あるいは尤度の確率分布をポアソン分布で近似した数１２式に現れるλ_i（ｄ）（i＝１、２）は平均画素値となる。これはノイズがないときの平均画素値であり、数９式に基づいて決定することができる。数９式あるいは数１０式を、数２１式、数２４式に代入することで最急降下法による逐次近似解を数７式から求めることができる。数９式及び数１０式における放射線のスペクトル（Φ（Ｅ））の幅を考慮することで、放射線が透過したときのエネルギー変化の影響(ビームハードニング)を抑制することが可能になる。

一方、以下の数２５式に示すように放射線スペクトルを、デルタ関数（ｅ＝１、２）を用いて単色近似することで、より簡潔に解を表現することができる。計算精度が必要な場合は数９式や数１０式を用いて計算を行い、計算速度が必要である場合、数９式及び数１０式における放射線のスペクトル（Φ（Ｅ））に単色近似を導入にして計算を行えばよい。平面分布取得部１１２は、異なる放射線エネルギーに対応した放射線スペクトルを単色の放射線スペクトルとして近似して出力画像を取得する。

単色近似を導入したとき、数９式、あるいは、数１０式は、以下の数２６式となり、物質の厚さｄの座標も考慮すると、数２６式は数２７式のように表すことができる。

数２７式を微分すると数２８式、数２９式のように表すことができる。

数２７式、数２８式、数２９式を正規分布に関する数２１式に代入すると、数２１式の各成分は、以下の数３０式及び数３１式のように表すことができる。

同様に数２７式、数２８式、数２９式をポアソン分布に関する数２４式に代入することで、単色近似を導入した場合におけるポアソン分布に関する数２４式の各成分を求めることができる。尚、ここでは、放射線スペクトルの単色近似について説明したが、正規分布やポアソン分布の再生性、つまり独立な確率変数Ｘ_１とＸ_２がそれぞれ正規分布やポアソン分布(パラメタλ_１とλ_２)に従うとき、確率変数の和Ｘ_１＋Ｘ_２も正規分布やポアソン分布に従うことを考慮すると、放射線スペクトルの幅を考慮したときの画素値Ｉ＝Σ_ｉi（Ｅ_i）についても同様の議論を展開することができる。この場合、平面分布取得部１１２は、異なる放射線エネルギーに対応した放射線スペクトルに基づいて取得された各画像の画素値を加算して出力画像を取得し、取得した出力画像に対して本実施形態の処理を適用することができる。

（その他の実施形態）
本発明は、上述の実施形態の１以上の機能を実現するプログラムを、ネットワーク又は記憶媒体を介してシステム又は装置に供給し、そのシステム又は装置のコンピュータにおける１つ以上のプロセッサーがプログラムを読出し実行する処理でも実現可能である。また、１以上の機能を実現する回路（例えば、ＡＳＩＣ）によっても実現可能である。

本発明は上記実施形態に制限されるものではなく、本発明の精神及び範囲から離脱することなく、様々な変更及び変形が可能である。

１００：放射線撮影装置、１０１：放射線管球（放射線発生装置）、１０３：ＦＰＤ（放射線検知部）：１０４：操作盤、１０８：照射指示部、１０９：コンピュータ（画像処理装置）、１１０：画像取得部、１１１：事前情報取得部、１１２：平面分布取得部、１１３：事前情報保存部、１１４：入力装置、１１５：表示装置、１１６：Ｃアーム、１１７：架台、１１８：同期装置

Claims (16)

1. 異なる放射線エネルギーで撮影した複数の放射線画像を入力画像とし、前記入力画像に含まれる物質について平面分布を示す出力画像を取得する平面分布取得手段を備え、
   前記平面分布取得手段は、尤度と、注目画素の画素値または周囲の画素の画素値を入力とする事前確率とから、入力が前記入力画像であるときに、出力が前記平面分布の値となる確率が最大となるように前記出力画像を決定する画像処理装置。
2. 前記平面分布取得手段は、前記事前確率を、前記出力画像における前記注目画素の画素値と前記周囲の画素の画素値との差分に基づいて取得する、又は、前記出力画像における前記注目画素の画素値に基づいて取得する請求項１に記載の画像処理装置。
3. 前記平面分布取得手段は、前記差分により前記平面分布の連続性を評価し、前記評価に基づいて前記出力画像におけるノイズを低減する請求項２に記載の画像処理装置。
4. 前記平面分布取得手段は、前記出力画像を反復演算による逐次近似解により取得する請求項１乃至３のいずれか１項に記載の画像処理装置。
5. 前記平面分布取得手段は、前記反復演算の繰り返し回数に応じて、反復の速度を表す係数が小さくなるように前記係数の設定を更新する請求項４に記載の画像処理装置。
6. 前記平面分布取得手段は、ノイズがないと仮定したときの平面分布の値を初期値として前記反復演算を行う、又は、前記出力画像に対する平滑化処理により取得した平面分布の値を初期値として前記反復演算を行う請求項４または５に記載の画像処理装置。
7. 前記平面分布取得手段は、前記異なる放射線エネルギーに対応した放射線スペクトルに基づいて取得された各画像の画素値を加算して前記出力画像を取得する、又は、前記異なる放射線エネルギーに対応した放射線スペクトルを単色の放射線スペクトルとして近似して前記出力画像を取得する請求項１乃至６のいずれか１項に記載の画像処理装置。
8. 前記平面分布取得手段は、前記尤度の確率分布として、正規分布またはポアソン分布のいずれかを用いる請求項１乃至７のいずれか１項に記載の画像処理装置。
9. 前記平面分布取得手段は、前記出力画像として、前記複数の放射線画像に含まれる複数の物質に分離した画像を取得する請求項１乃至８のいずれか１項に記載の画像処理装置。
10. 前記平面分布取得手段は、前記出力画像として、被写体を構成する軟部組織の厚さ画像と、前記被写体を構成する硬部組織の厚さ画像を取得する請求項９に記載の画像処理装置。
11. 前記平面分布取得手段は、前記軟部組織の厚さ画像における、注目画素の画素値と周囲の画素の画素値との差分により前記軟部組織の厚さ画像における連続性の評価を行う請求項１０に記載の画像処理装置。
12. 前記平面分布取得手段は、前記事前確率の影響の度合いを表す正則係数を、撮影条件と、放射線の照射条件と、被写体の体形に関する情報との少なくとも一つに基づいて変更する請求項１乃至１１のいずれか１項に記載の画像処理装置。
13. 異なる放射線エネルギーで撮影した被写体の複数の放射線画像における注目画素の画素値または周囲の画素の画素値を入力とする事前確率と尤度とを用いて、前記被写体に含まれる物質について平面分布を取得する平面分布取得手段を備える画像処理装置。
14. 異なる放射線エネルギーで撮影した複数の放射線画像を入力画像とし、前記入力画像に含まれる物質について平面分布を示す出力画像を取得する工程と、
    尤度と、注目画素の画素値または周囲の画素の画素値を入力とする事前確率とから、入力が前記入力画像であるときに、出力が前記平面分布の値となる確率が最大となるように前記出力画像を決定する工程と、を有する画像処理方法。
15. 異なる放射線エネルギーで撮影した被写体の複数の放射線画像における注目画素の画素値または周囲の画素の画素値を入力とする事前確率と尤度とを用いて、前記被写体に含まれる物質について平面分布を取得する工程を有する画像処理方法。
16. コンピュータを、請求項１乃至１３のいずれか１項に記載の画像処理装置の手段として機能させるプログラム。

Images