(57)【特許請求の範囲】
【請求項１】
同一の被写体を異なるタイミングでそれぞれ撮像して得られた第１画像および第２画像を取得する画像取得部と、
複数の画像変形量を用いて前記第１画像をそれぞれ変形させて、該変形した第１画像のそれぞれを所定の分割条件に従って、複数の第１分割画像に分割するとともに、前記第２画像を前記第１分割画像に対応する複数の第２分割画像に分割し、前記第１分割画像および該第１分割画像に対応する前記第２分割画像の一対ごとに、該一対の分割画像のそれぞれの画素値の分布の類似性を表す分割画像類似度を算出し、複数の該分割画像類似度に基づいて前記変形した第１画像と前記第２画像の画像全体の類似性を評価し、評価の高い前記画像変形量を前記第１画像中の前記被写体の位置を前記第２画像中の前記被写体の位置に類似させる前記第１画像の画像変形量として推定する画像変形量推定部と、
前記推定された変形量に基づいて、前記第１画像を変形させた画像を生成する画像生成部とを備えることを特徴とする画像処理装置。
【請求項２】
前記所定の分割条件に従った分割が、所定の第１の空間パラメータに基づいて前記変形した第１画像を複数の第１分割画像に分割し、前記第１の空間パラメータに対応する第２の空間パラメータに基づいて前記第２画像を前記第１分割画像に対応する複数の第２分割画像に分割するものであり、
前記第１の空間パラメータが、前記第１画像内の所定の図形からの距離を表すものであり、前記第２の空間パラメータが、前記所定の図形に対応する前記第２画像内の図形からの距離を表すものであることを特徴とする請求項１記載の画像処理装置。
【請求項３】
前記所定の図形が点であり、
前記分割条件に従った分割が、前記第１の空間パラメータに応じて、該第１の空間パラメータの範囲ごとに前記第１画像を分割するとともに、前記第１の空間パラメータの範囲に対応する前記第２の空間パラメータの範囲ごとに前記第２画像を分割するものであることを特徴とする請求項２記載の画像処理装置。
【請求項４】
前記所定の分割条件に従った分割が、所定の第１の空間パラメータに基づいて前記変形した第１画像を複数の第１分割画像に分割し、前記第１の空間パラメータに対応する第２の空間パラメータに基づいて前記第２画像を前記第１分割画像に対応する複数の第２分割画像に分割するものであり、
前記第１の空間パラメータが、前記第１画像内の所定の図形に対する角度を表すものであり、前記第２の空間パラメータが、前記所定の図形に対応する前記第２画像内の図形に対する角度を表すものであることを特徴とする請求項１記載の画像処理装置。
【請求項５】
前記所定の図形が、前記第１画像内の所定の軸であり、
前記第１の空間パラメータが、前記第１画像内の前記所定の軸からの角度を表すものであり、
前記第２の空間パラメータが、前記第２画像内の前記所定の軸に対応する軸からの角度を表すものであることを特徴とする請求項４記載の画像処理装置。
【請求項６】
前記画像変形量推定部が、前記複数の分割画像類似度の総和が最大となる前記画像変形量を、前記第１画像中の前記被写体の位置を前記第２画像中の前記被写体の位置に類似させる前記第１画像の画像変形量として推定するものであることを特徴とする請求項１から５のいずれか１項記載の画像処理装置。
【請求項７】
前記分割画像類似度が、前記第１分割画像の画素値および前記第２分割画像の画素値を離散確率変数とした相互相関情報量を表すものであることを特徴とする請求項１から６のいずれか１項記載の画像処理装置。
【請求項８】
前記所定の図形が、ユーザが入力した位置により特定される点であることを特徴とする請求項３記載の画像処理装置。
【請求項９】
前記所定の図形が、自動認識技術により抽出された解剖学構造物の特徴的な位置により特定される点であることを特徴とする請求項３記載の画像処理装置。
【請求項１０】
画像取得部と、画像変形量推定部と、画像生成部を備えた画像処理装置に実行させる画像処理方法において、
同一の被写体を異なるタイミングでそれぞれ撮像して得られた第１画像および第２画像を取得する工程と、
複数の画像変形量を用いて前記第１画像をそれぞれ変形させて、該変形した第１画像のそれぞれを所定の分割条件に従って、複数の第１分割画像に分割するとともに、前記第２画像を前記第１分割画像に対応する複数の第２分割画像に分割し、前記第１分割画像および該第１分割画像に対応する前記第２分割画像の一対ごとに、該一対の分割画像のそれぞれの画素値の分布の類似性を表す分割画像類似度を算出し、複数の該分割画像類似度に基づいて前記変形した第１画像と前記第２画像の画像全体の類似性を評価し、評価の高い前記画像変形量を前記第１画像中の前記被写体の位置を前記第２画像中の前記被写体の位置に類似させる前記第１画像の画像変形量として推定する工程と、
前記推定された変形量に基づいて、前記第１画像を変形させた画像を生成する工程を備えたことを特徴とする画像処理方法。
【請求項１１】
コンピュータを、
同一の被写体を異なるタイミングでそれぞれ撮像して得られた第１画像および第２画像を取得する画像取得部と、
複数の画像変形量を用いて前記第１画像をそれぞれ変形させて、該変形した第１画像のそれぞれを所定の分割条件に従って、複数の第１分割画像に分割するとともに、前記第２画像を前記第１分割画像に対応する複数の第２分割画像に分割し、前記第１分割画像および該第１分割画像に対応する前記第２分割画像の一対ごとに、該一対の分割画像のそれぞれの画素値の分布の類似性を表す分割画像類似度を算出し、複数の該分割画像類似度に基づいて前記変形した第１画像と前記第２画像の画像全体の類似性を評価し、評価の高い前記画像変形量を前記第１画像中の前記被写体の位置を前記第２画像中の前記被写体の位置に類似させる前記第１画像の画像変形量として推定する画像変形量推定部と、
前記推定された変形量に基づいて、前記第１画像を変形させた画像を生成する画像生成部として機能させることを特徴とする画像処理プログラム。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、同一の被写体を撮影して取得した２つ画像において、両画像中の被写体の空間的な位置が一致するように、一方の画像の画像空間を変形させた画像を生成する画像処理装置および画像処理方法、並びに、画像処理プログラムに関する。
【背景技術】
【０００２】
同一の被写体を同一または異なる撮像装置（モダリティ）を用いて異なるタイミングで撮影した２つの３次元画像を用いた画像診断において、両画像を重ねたとき、両画像中の被写体の空間的な位置が一致する変換関数を推定し、推定した変換関数を用いて一方の画像を変形させて２種類の画像を位置合わせする非剛体レジストレーション技術が注目されている。この非剛体レジストレーション技術においては、画像空間を所定間隔で区切る制御点が設定され、この制御点の位置を変位させて、変形した一方の画像と他方の画像のそれぞれの画素値の類似性を評価する評価関数が最大となる制御点の画像変形量を決定し、このときの制御点の画像変形量に基づいて変換関数が推定される。
【０００３】
非特許文献１は、種類の異なるモダリティであるＰＥＴ装置（Positron emission tomography）およびＣＴ装置(Computed tomography)により同一の被写体をそれぞれ撮影して得られた画像間に非剛体レジストレーション処理を適用するものであり、評価関数に、２種類の画像の画素値の類似性の尺度として相互相関情報量（Mutual information）を用いている。
【先行技術文献】
【非特許文献】
【０００４】
【非特許文献１】David Mattes, David R. Haynor, Hubert Vesselle, Thomas K. Lewellen, William Eubank, “Nonrigid multimodality image registration”, Proceedings of the SPIE, volume 4322, p.1609-1620, 2001
【発明の概要】
【発明が解決しようとする課題】
【０００５】
しかしながら、非特許文献１に記載された方法では、２つの画像の画素値の分布の相関性のみに基づいて２つの画像の類似性を判断しているため、画像の空間的な特徴を判別できず、類似性の判定を誤る場合がある。
【０００６】
例えば、非特許文献１に記載された方法では、同じ画素値範囲に属する複数の被写体を撮影した２つの画像間において被写体の数や空間的位置が異なっていても、２つの画像間で同じ画素値範囲に属する画素の総数が同じ場合には、類似していると判断してしまう。具体的には、患者の胸部を撮影した一つの画像において、膵臓と肝臓を表す各画素の画素値は、所定の同じ画素値範囲に属するものとなる。もし、同一の患者の胸部を撮影した他方の画像に、この所定の範囲の画素値範囲を示す１つの陰影があり、この陰影の体積が、先述の一つの画像の肝臓と膵臓の体積の総和と等しい場合、非特許文献１の方法では両画像が類似していると判断してしまう可能性がある。
【０００７】
本発明は、かかる問題点に鑑みて、同一の画素値範囲に属する複数の被写体を含む画像間において、被写体の空間的な特徴を類似度の評価に反映させることにより精度良く類似度の評価を行って、より正確に被写体が一致するように一方の画像を変形させた画像を生成することができる画像処理装置、画像処理方法、及び画像処理プログラムを提供することを目的とする。
【課題を解決するための手段】
【０００８】
上記目的を達成するために、本発明に係る画像処理装置は、同一の被写体を異なるタイミングでそれぞれ撮像して得られた第１画像および第２画像を取得する画像取得部と、前記第１画像を変形させて、該変形した第１画像と前記第２画像との類似度を、前記変形した第１画像の画素値と前記第２画像の対応する画素値の分布の相関性を評価する評価関数により評価して前記変形した第１画像と前記第２画像が類似する前記第１画像の画像変形量を推定する画像変形量推定部と、前記推定された変形量に基づいて、第１画像を変形させた画像を生成する画像生成部とを備え、前記評価関数が、所定の分割条件に従って、前記変形した第１画像を複数の第１分割画像に分割するとともに前記第２画像を前記第１分割画像に対応する複数の第２分割画像に分割し、前記第１分割画像およびこれに対応する前記第２分割画像の一対ごとに、該複数の一対の画像のそれぞれの画素値の分布の類似性を表す分割画像類似度に基づいて、前記第１画像を変形した画像と前記第２画像の類似性を評価するものであることを特徴とする。
【０００９】
本発明に係る画像処理方法は、上記画像処理装置に実行させる画像処理方法において、画像取得部が同一の被写体を異なるタイミングでそれぞれ撮像して得られた第１画像および第２画像を取得する工程と、画像変形量推定部が前記第１画像を変形させて、該変形した第１画像と前記第２画像との類似度を、前記変形した第１画像の画素値と前記第２画像の対応する画素値の分布の相関性を評価する評価関数により評価して前記変形した第１画像と前記第２画像が類似する前記第１画像の画像変形量を推定する工程と、画像生成部が前記推定された変形量に基づいて、第１画像を変形させた画像を生成する工程を備えた画像処理方法であって、前記評価関数が、所定の分割条件に従って、前記変形した第１画像を複数の第１分割画像に分割するとともに前記第２画像を前記第１分割画像に対応する複数の第２分割画像に分割し、前記第１分割画像およびこれに対応する前記第２分割画像の一対ごとに、該複数の一対の画像のそれぞれの画素値の分布の類似性を表す分割画像類似度に基づいて、前記第１画像を変形した画像と前記第２画像の類似性を評価するものであることを特徴とする。
【００１０】
本発明に係る画像処理プログラムは、コンピュータを、同一の被写体を異なるタイミングでそれぞれ撮像して得られた第１画像および第２画像を取得する画像取得部と、前記第１画像を変形させて、該変形した第１画像と前記第２画像との類似度を、前記変形した第１画像の画素値と前記第２画像の対応する画素値の分布の相関性を評価する評価関数により評価して前記変形した第１画像と前記第２画像が類似する前記第１画像の画像変形量を推定する画像変形量推定部と、前記推定された変形量に基づいて、第１画像を変形させた画像を生成する画像生成部として機能させる画像処理プログラムであって、前記評価関数が、所定の分割条件に従って、前記変形した第１画像を複数の第１分割画像に分割するとともに前記第２画像を前記第１分割画像に対応する複数の第２分割画像に分割し、前記第１分割画像およびこれに対応する前記第２分割画像の一対ごとに、該複数の一対の画像のそれぞれの画素値の分布の類似性を表す分割画像類似度に基づいて、前記第１画像を変形した画像と前記第２画像の類似性を評価するものであることを特徴とする。
【００１１】
上記本発明に係る画像処理装置および方法並びにプログラムにおいて、第１画像および第２画像は、同じ被写体を撮影して得られたものであれば、同じモダリティにより撮影して得られたものであってもよく、異なるモダリティにより撮影して得られたものであってもよい。また、第１および第２のモダリティとして適用可能なモダリティとして、ＰＥＴ、ＣＴ、ＭＲＩ、ＳＰＥＣＴ、超音波画像などがあげられる。
【００１２】
また、本発明に係る画像処理装置において、前記評価関数は、所定の分割条件に従って、前記変形した第１画像を複数の第１分割画像に分割するとともに前記第２画像を前記第１分割画像に対応する複数の第２分割画像に分割し、前記第１分割画像およびこれに対応する前記第２分割画像の一対ごとに、該複数の一対の画像のそれぞれの画素値の分布の類似性を表す分割画像類似度に基づいて、前記第１画像を変形した画像と前記第２画像の類似性を評価するものであればいかなるものでもよい。例えば、評価関数を、両画像が類似しているほど評価値が大きくなるように定義してもよく、両画像が類似しているほど評価値が小さくなるように定義してもよい。
【００１３】
また、本発明に係る画像処理装置において、前記評価関数が、前記変形した第１画像の画素値および前記第２画像の画素値を離散確率変数とした相互相関情報量を表すものであることが好ましい。
【００１４】
なお、ここでいう図形とは、一定の決まりによって定められる形状を意味し、点、直線、平面、曲線、球などの立体またはその一部を含む。第１画像および第２画像の対応する図形は、点、直線、曲線、平面、曲面や球などの立体やこれらの一部など任意に定義してよいが、両図形は第１画像と第２画像の対応する位置に位置することが必要であり、さらに、計算の便宜上、同一の種類の図形であることが望ましい。また、第１および第２の空間パラメータは、変形した第１画像Ｖ１ａおよび第２画像間の対応する図形に対して同様の決まりにより規定される同一の種類のパラメータであればよく、１つであってもよく複数であってもよい。例えば、第１および第２の空間パラメータとして距離や角度などを用いることができる。
【００１５】
なお、対応する図形の初期位置は周知の任意の方法で特定してよい。例えば第１画像と第２画像の対応する位置をユーザが入力した位置により特定してもよく、周知の自動認識技術により得られた解剖学構造物の特徴的な位置から特定してもよい。
【００１６】
また、本発明に係る分割条件は、一定の決まりに基づいて、変形した第１画像及び第２画像を分割できるものであれば、いかなるものであってもよい。
【００１７】
一例として、本発明に係る画像処理装置において、前記分割条件が、所定の第１の空間パラメータに基づいて前記変形した第１画像を複数の第１分割画像に分割するものであり、前記第１の空間パラメータに対応する第２の空間パラメータに基づいて前記第２画像を前記第１分割画像に対応する複数の第２分割画像に分割するものであることが好ましい。
【００１８】
また、本発明に係る画像処理装置において、前記第１の空間パラメータが、前記第１画像内の所定の図形からの距離を表すものであり、前記第２の空間パラメータが、前記所定の図形に対応する第２画像内の図形からの距離を表すものであってもよい。
【００１９】
また、上記場合に、前記所定の図形が点であり、前記分割条件が、前記第１の空間パラメータに応じて、前記第１の空間パラメータの範囲ごとに前記第１画像を分割するとともに、該前記第１の空間パラメータの範囲に対応する前記第２の空間パラメータの範囲ごとに前記第２画像を分割するものであることが好ましい。
【００２０】
また、本発明に係る画像処理装置において、前記第１の空間パラメータが、前記第１画像内の所定の図形に対する角度を表すものであり、前記第２の空間パラメータが、前記所定の図形に対応する第２画像内の図形に対する角度を表すものであることが好ましい。
【００２１】
また、上記場合に、前記所定の図形が、前記第１画像内の所定の軸であり、前記第１の空間パラメータが、前記第１画像内の前記所定の軸からの角度を表すものであり、前記第２の空間パラメータが、前記第２画像内の前記所定軸に対応する軸からの角度を表すものであることが好ましい。
【００２２】
上記の所定の軸は、ユーザの要求により任意に設定でき、例えば、ｘｙｚ座標系におけるｘｙ平面やｘ、ｙ、ｚ軸などを用いることができる。
【００２３】
また、本発明に係る画像処理装置において、前記評価関数が、前記分割画像類似度の総和により、前記変形した第１画像および前記第２画像の類似度を定義するものであることが好ましい。
【発明の効果】
【００２４】
本発明による画像処理装置は、同一の被写体を異なるタイミングでそれぞれ撮像して得られた第１画像および第２画像を取得する画像取得部と、前記第１画像を変形させて、該変形した第１画像と前記第２画像との類似度を、前記変形した第１画像の画素値と前記第２画像の対応する画素値の分布の相関性を評価する評価関数により評価して前記変形した第１画像と前記第２画像が類似する前記第１画像の画像変形量を推定する画像変形量推定部と、前記推定された変形量に基づいて、第１画像を変形させた画像を生成する画像生成部とを備え、前記評価関数が、所定の分割条件に従って、前記変形した第１画像を複数の第１分割画像に分割するとともに前記第２画像を前記第１分割画像に対応する複数の第２分割画像に分割し、前記第１分割画像およびこれに対応する前記第２分割画像の一対ごとに、該複数の一対の画像のそれぞれの画素値の分布の類似性を表す分割画像類似度に基づいて、前記第１画像を変形した画像と前記第２画像の類似性を評価するものであるので、同一の画素値範囲に属する複数の被写体を含む画像間において、被写体の空間的な特徴を類似度の評価に反映させることにより、従来よりも正確に類似度を評価することができる。結果として好適に第１画像を第２画像に一致させた画像を生成することができる。
【図面の簡単な説明】
【００２５】
【図１】本発明の第１の実施形態における画像処理装置の電気的な概略ブロック図である。
【図２】本発明の第１の実施形態における画像処理装置の動作を示すフローチャートである。
【図３Ａ】本発明の第１の実施形態において、変形した第１画像と第２画像を基準点からの距離に基づいて分割した図を表す図である。
【図３Ｂ】本発明の第１の実施形態において、変形した第１画像と第２画像をｘ軸からの角度に基づいて分割した図を表す図である。
【図４Ａ】位置合わせ前の第１画像（ＭＲ画像）と第２画像（ＣＴ画像）の一例を示す図である。
【図４Ｂ】位置合わせ後の変形した第１画像（ＭＲ画像）と第２画像（ＣＴ画像）の一例を示す図である。
【図５】本発明の第２の実施形態における画像処理装置の電気的な概略ブロック図である。
【図６】本発明の第２の実施形態における画像処理装置の動作を示すフローチャートである。
【発明を実施するための形態】
【００２６】
以下、本発明の画像処理装置、画像処理プログラムおよび画像処理方法の実施形態について、図面を参照しながら詳細に説明する。本発明は、異なるモダリティにより同一の被写体を異なるタイミングで撮影して取得した２つの画像を位置合わせする処理を行う様々な分野で応用可能であるが、ここでは、医療分野における画像診断に本発明を適用した例をもとに説明を行う。
【００２７】
図１に、医師が使用するワークステーションに、画像処理プログラムをインストールすることにより実現された画像処理装置の概略構成を示す。画像処理装置１は、標準的なワークステーションの構成として、プロセッサおよびメモリ（いずれも図示せず）を備え、さらに、ＨＤＤ（Hard Disk Drive）等のストレージ２を備えている。また、画像処理装置１には、ディスプレイ３と、マウス、キーボード等の入力装置４が接続されている。
【００２８】
画像処理プログラムと画像処理プログラムが参照するデータは、インストール時にストレージ２に記憶され、起動時にメモリにロードされる。画像処理プログラムは、ＣＰＵに実行させる処理として、画像取得処理と、画像変形量推定処理と、画像生成処理と、表示制御処理とを規定している。
【００２９】
そして、プログラムの規定にしたがって、ＣＰＵが上記各処理を実行することにより、汎用のワークステーションは、後述の画像取得部１１と、画像変形量推定部１２と、画像生成部１４と、表示制御部１５として機能する。
【００３０】
ストレージ２には、撮影を担当する検査部門から転送された、第１画像および第２画像、もしくはデータベース検索により取得された第１画像および第２画像が記憶される。本実施の形態では、ある患者の検査において、同じ日の異なる時間に異なるモダリティにより一人の患者の頭部を撮影して得られた、第１画像Ｖ１（ＭＲ画像）および第２画像Ｖ２（ＣＴ画像）がそれぞれ、検査部門から転送され、両画像Ｖ１、Ｖ２のヘッダ情報などの任意の情報から周知の方法で取得した画素サイズ（Pixel Spacing）とスライス間隔(Slice Spacing)に基づいて両画像Ｖ１、Ｖ２のいずれか一方に両画像Ｖ１、Ｖ２のスケールを統一する処理が施され、ストレージ２に記憶されている。ここでは、図３Ａに第１画像Ｖ１および第２画像Ｖ２のイメージ画像を示す。第１画像Ｖ１および第２画像Ｖ２は、第１画像内の頭蓋骨などの解剖学的構造物が、やや第２画像の対応する解剖学的構造物より大きく、また、左右の大脳半球の境界部分の位置など解剖学的構造物の特徴的な位置もずれている。
【００３１】
画像取得部１１は、第１画像Ｖ１および第２画像Ｖ２をストレージ２からメモリに取得する。本実施形態では、画像処理装置１は、選択メニューにおいて所定の位置合わせ機能が選択されたことを検出すると、ユーザに、第１および第２画像の特定に必要な情報の選択または入力を促す。そして、ユーザの入力装置４の操作により、第１および第２画像が特定されると、画像取得部１１は、第１画像Ｖ１および第２画像Ｖ２をストレージ２からメモリに取得する。
【００３２】
画像変形量推定部１２は、第１画像Ｖ１を変形させて、変形した第１画像Ｖ１ａと第２画像Ｖ２との類似度を、変形した第１画像Ｖ１ａの画素値と第２画像Ｖ２の対応する画素値の分布の相関性を評価する評価関数Ｓ（μ）により評価して変形した第１画像と第２画像の被写体が一致する第１画像の画像変形量を推定する。
【００３３】
さらに詳細には、画像変形量推定部１２は、変形した第１画像Ｖ１ａおよび第２画像Ｖ２において、画像空間を所定間隔で区切る制御点ｘ１、ｘ２、…、ｘｎからなる集合Ｘをそれぞれ設定する。以下、制御点ｘ１、ｘ２、…、ｘｎの集合を制御点Ｘと記載する。なお、画像変形量推定部１２は、第１画像Ｖ１の制御点Ｘを変換関数ｇにより画像変形量μだけ変位させることにより、第１画像Ｖ１を変形する。以下、第１画像Ｖ１の制御点Ｘを変換関数ｇにより画像変形量μだけ変位させた制御点をｇ(Ｘ,μ)と記載し、第１画像Ｖ１の制御点Ｘを変換関数ｇにより画像変形量μだけ変位させることにより第１画像Ｖ１を変形した画像を、変形した第１画像Ｖ１ａと記載する。なお、本実施形態における画像変形量推定部１２は、非特許文献１に開示された方法により制御点Ｘや画像変形量μを設定し、非特許文献１に記載された変換関数ｇを変換関数ｇとして用いる。
【００３４】
次いで、画像変形量推定部１２は、変形した第１画像Ｖ１ａの制御点ｇ(Ｘ,μ)における画素値Ｍ(ｇ(Ｘ,μ))を取得するとともに第２画像Ｖ２の制御点Ｘの画素値Ｆ(Ｘ)を取得する。そして、変形した第１画像Ｖ１ａの各制御点ｇ(Ｘ,μ)における画素値Ｍ(ｇ(Ｘ,μ))と第２画像Ｖ２の各制御点Ｘの画素値Ｆ(Ｘ)の類似性を評価する評価関数Ｓ（μ）（レジストレーション関数）が最大となる制御点Ｘの画像変形量μを決定し、このときの制御点Ｘの画像変形量μに基づいて第１画像Ｖ１に対する変換関数を推定する。
【００３５】
本実施形態では、画像変形量推定部１２は、後述の式（３）で表す評価関数Ｓを用いて、変形した第１画像Ｖ１ａと第２画像Ｖ２の類似度を評価する。本実施形態における評価関数Ｓによれば、変形した第１画像Ｖ１ａと第２画像Ｖ２の画素値分布が類似するほど評価関数Ｓの値が大きくなるため、画像変形量推定部１２は、画像変形量μを変化させながら、評価関数Ｓ（μ）の変化量（または、μに対する偏導関数▽Ｓ（μ）の絶対値）が所定のしきい値以下となる場合の画像変形量μを、両画像の類似度が最大値をとる画像変形量（両画像が最も類似する画像変形量）と判断する。そして、この画像変形量μに基づいて、第１画像Ｖ１を変形する変換関数を決定する。なお、所定のしきい値は、式（３）に示す評価関数Ｓ（μ）の変化量が十分小さい値であると見なせる任意の値を設定してよい。また、この画像変形量μに基づいて、第１画像Ｖ１を変形する変換関数を決定する方法は周知の種々の方法を適用でき、ここでは非特許文献１に記載された方法を適用するものとする。
【００３６】
なお、ここでは評価関数を、類似度（評価値）が大きいほど類似度が大きくなるように定義しているが、評価値が小さいほど類似度が大きくなるように定義してもよい。また、評価関数の評価値（類似度）が最大（または最小）となる画像変形量を特定可能な方法であれば、非剛体レジストレーション技術の評価関数の最大値（または最小値）を算出するいかなる方法を用いて類似度が最大となる画像変形量μを特定してもよい。また、複数の異なる画像変形量μに対して評価関数Ｓ（μ）による評価値をそれぞれ算出して、算出した評価値のうち評価値が最大（または最小）であるものを特定し、特定された最大評価値（または最小評価値）に対応する画像変形量μを特定してもよい。
【００３７】
本実施形態の画像変形量推定部１２は、所定の分割条件に従って変形した第１画像を複数の第１分割画像に分割するとともに第２画像を第１分割画像に対応する複数の第２分割画像に分割し、第１分割画像およびこれに対応する第２分割画像の一対ごとに、一対の分割画像間の画素値の分布の類似性を表すように定義された分割画像類似度に基づいて、第１画像を変形した画像と第２画像の類似性を評価する。
【００３８】
ここで、従来の評価関数について説明し、その後本実施形態の評価関数について詳細に説明する。
【００３９】
非特許文献１に示されるように、異なるモダリティにより同一の被写体をそれぞれ撮影して取得した第１画像と第２画像の位置合わせを行う非剛体レジストレーションの技術において、相互相関情報量（mutual information、相互情報量）に基づく評価関数を用いて、第１画像の画素値と第２画像の対応する画素値の分布の類似度を評価することができる。
【００４０】
相互相関情報量とは、２つの確率変数ｆ、ｍの同時確率分布関数ｐ（ｆ，ｍ）と２つの確率変数ｆ、ｍの周辺確率分布関数ｐ（ｆ）、ｐ（ｍ）に基づいて２つの確率変数ｆ、ｍの相関性の尺度を表すもので、代表的には下記式（１）により定義されるものである。式（１）において、ｆは集合Ｆに属する離散確率変数であり、ｍは集合Ｍに属する離散確率変数である。また、相互相関情報量は、２つの確率変数ｆ、ｍが相関しているほど情報量が大きくなるものであり、言い換えると、一方の変数が与えられたときに他方の変数を高い確率で推測できるほど情報量が大きくなるものである。なお、２つの確率変数が完全に独立であれば、相互情報量は０となる。
【数１】
0005955199.tif 000002
【００４１】
異なるモダリティにより同一の被写体を撮影した画像間においては、各画像の撮影の物理的原理が異なるため、同じ種類の被写体であっても画像ごとに異なる画素値（信号値）を示す場合があり、単に画素値を比較しても両画像の類似性を判断することができない。例えば、ＣＴ画像においては、放射線の吸収率（透過度）が高いほど画素値（ＣＴ値）が大きくなり、空気、水、筋肉や肝臓や心臓などの臓器組織、骨の順に画素値が大きくなる。一方、ＭＲ画像においては、撮影対象物に含まれる水素原子の核磁気共鳴に応じて画素値が決定され、Ｔ１、Ｔ２強調画像などの撮影方法に応じて物質の信号値の大きさが変化する。例えば、ＭＲ画像におけるＴ１強調画像では、脂肪、筋肉、水は画素値の大きさがこの順に小さくなる。しかしながら、ＣＴ画像では、脂肪、筋肉、水は画素値の大きさの順はＭＲ画像のＴ１強調画像とは異なっているので、両者の画素値を単純に比較しても両画像の類似を判断することができない。
【００４２】
しかし、このような場合であっても、両画像の同一の解剖学的構造物を表している部分は、同一の解剖学的構造物に基づく共通の特徴に従って画素値が分布しているため、両者の画素値分布は互いに相関するものとなる。このことを利用して、非特許文献１では、異なるモダリティにより同一の被写体を撮影した２つの画像を、両者の画素値をそれぞれ確率変数とする相互相関情報量に基づいて、両画像の画素値の分布の相関が大きいほど、両画像が類似していると判断している。詳細には、式（２）で示すように、評価関数として、第１のモダリティにより撮影して得られた画像を変形した画像の画素値ｍ;μと、第２のモダリティにより撮影して得られた画像の各画素値ｆをそれぞれ確率変数とする相互相関情報量を用いている。（以下、第１画像を画像変形量μだけ変形した画像の各画素値を、第１画像からの画像変形量μを用いてｍ;μと表す。）式（２）において、相互相関情報量は、第２画像と第１画像を変形した画像のそれぞれの画素値ｆ、ｍ;μの分布の相関を表しているため、第２画像と第１画像を変形した画像の類似性の尺度として機能している。ここで、集合Ｆは変形第１画像の全画素値からなる集合であり、集合Ｍは第２画像の全画素値からなる集合である。
【数２】
0005955199.tif 000003
【００４３】
本発明者は、非特許文献１に記載された式（２）で表されるような従来の２つの画像の画素値の分布の相関性のみに基づいて２つの画像の類似性を判断している評価関数に関し、画像の空間的な特徴を判別できず、類似性の判定を誤る場合があるという問題（第１の問題）に着目した。
【００４４】
例えば、非特許文献１に記載された方法では、同じ画素値範囲に属する複数の被写体を撮影した２つの画像間において被写体の数や空間的位置が異なっていても、２つの画像間で同じ画素値範囲に属する画素の総数が同じ場合には、類似していると判断してしまう。具体的には、患者の胸部を撮影した一つの画像において、膵臓と肝臓を表す各画素の画素値は、所定の同じ画素値範囲に属するものとなる。もし、同一の患者の胸部を撮影した他方の画像に、この所定の範囲の画素値範囲を示す１つの陰影があり、この陰影の体積が、先述の一つの画像の肝臓と膵臓の体積の総和と等しい場合、非特許文献１の方法では両画像が類似していると判断してしまう可能性がある。
【００４５】
本発明者は上記問題に鑑み、第１変形画像と第２画像を互いに対応する領域ごとに所定の分割条件に従って分割し、この分割した領域ごとに類似度を算出し、分割した領域ごとの類似度に基づいて、第１画像と第２画像の類似度を評価することが有効であることを見出した。この場合には、離れて位置する同じ画素値範囲に属する複数の被写体が別々の分割された領域に含まれる可能性が高いため、互いに異なる分割領域に位置する複数の被写体を表す画素値間で相関性が評価されることが抑制され、上記第２の問題の発生を低減することができるからである。
【００４６】
なお、分割条件は、所定の決まりに基づいて第１変形画像と第２画像を互いに対応する領域ごとに分割できるものであれば、いかなる方法でもよい。
【００４７】
例えば、分割条件として、変形した第１画像を所定の図形に対する所定の第１の空間パラメータに基づいて複数の第１分割画像に分割し、第２画像を第１の空間パラメータに対応する第２の空間パラメータに基づいて第１分割画像に対応する複数の第２分割画像に分割することを規定することが考えられる。なお、ここでいう図形とは、一定の決まりによって定められる形状を意味し、点、直線、平面、曲線、球などの立体またはその一部を含む。第１画像および第２画像の対応する図形は、点、直線、曲線、平面、曲面や球などの立体やこれらの一部など任意に定義してよいが、両図形は第１画像と第２画像の対応する位置に位置することが必要であり、さらに、計算の便宜上、同一の種類の図形であることが望ましい。また、第１および第２の空間パラメータは、変形した第１画像Ｖ１ａおよび第２画像間の対応する図形に対して同様の決まりにより規定される同じ種類のパラメータであればよく、１つであってもよく複数であってもよい。例えば、第１および第２の空間パラメータとして距離や角度などを用いることができる。
【００４８】
なお、上記第１分割画像と第２分割画像を対応する範囲を表すものとするために、第1分割画像および第２分割画像のスケールを統一する必要がある。このために、例えば、第１分割画像と第２分割画像のスケールを統一する処理が、ヘッダ情報などの任意の情報から周知の方法で取得した画素サイズ（Pixel Spacing）とスライス間隔(Slice Spacing)に基づいて、両画像の第１分割画像と第２分割画像の設定の前に、第１画像と第２画像に対して実施することが考えられる。なお、このスケールを統一する処理は、各第１および第２分割画像の設定後に、第１分割画像と第２分割画像のそれぞれに対して実施してもよい。
【００４９】
なお、対応する図形の初期位置は周知の任意の方法で特定してよい。例えば第１画像と第２画像の対応する位置をユーザが入力した位置により特定してもよく、周知の自動認識技術により得られた解剖学構造物の特徴的な位置から特定してもよい。
【００５０】
本第１の実施形態における分割条件は、変形した第１画像を所定位置からの距離に応じて複数の第１分割画像に分割し、第２画像を第１画像の所定位置に対応する位置からの距離に基づいて第１分割画像に対応する複数の第２分割画像に分割することを規定する。そして、第２の実施形態における評価関数Ｓ（μ）は、分割条件に従って、第１分割画像とこれに対応する第２分割画像とからなる一対の分割画像ごとに、一対の分割画像のそれぞれの画素値の分布の類似性を表す分割画像類似度を定義し、各分割画像類似度に基づいて、第１画像を変形した画像Ｖ１ａと第２画像Ｖ２の類似度を評価する。
【００５１】
ここで、評価関数Ｓ（μ）は、複数の分割画像類似度に基づいて変形した第１画像Ｖ１ａと第２画像Ｖ２との類似度を算出するものであれば、種々の方法を用いることができる。一例として、評価関数Ｓ（μ）は、複数の分割画像類似度の総和により定義できる。
【００５２】
図３Ａは第２の実施形態における第１および第２画像の分割方法を説明する図である。図３Ａにおける変形した第１画像Ｖ１ａおよび第２画像Ｖ２は、共に３次元画像の中心を原点とするｘｙ平面を表しており、変形した第１画像Ｖ１ａの各制御点Ｘを説明のために、ｘｙ平面上でのみ変位させて表している。図４Ａに示すように、本実施形態における分割条件は、変形した第１画像Ｖ１ａおよび第２画像Ｖ２のそれぞれの基準点Ｐ１、Ｐ２を設定し、各基準点Ｐ１、Ｐ２から距離に応じて、各基準点Ｐ１、Ｐ２からの距離の範囲に応じてＶ１ａおよび第２画像Ｖ２をそれぞれ分割することを規定している。詳細には、基準点からの距離に応じて、０≦ｄ＜ｄ０、ｄ０≦ｄ＜ｄ１、ｄ１≦ｄ＜ｄ２、…、ｄｋ−１≦ｄ＜ｄｋを満たすｋ個の球状または中空の球状の領域に，変形した第１画像Ｖ１ａおよび第２画像Ｖ２をそれぞれ分割する。そして、評価関数Ｓ（μ）を、各第１分割画像と第２分割画像のペア（Ａ０、Ｂ０）、（Ａ１、Ｂ１）、…、（Ａｋ、Ｂｋ）ごとの分割画像類似度の総和によって規定する。
【００５３】
また、評価関数として、以下の式（３）を用いる。
【数３】
0005955199.tif 000004
式（３）において、ｄは、各距離の範囲を意味し、集合Ｄは第１画像の所定点からの距離の範囲ｄ０、ｄ１、ｄ２、…、ｄｋ（ｋは正の整数）からなる集合を表す。なお、場合により、距離の範囲０≦ｄ＜ｄ０を距離の範囲ｄ０と、記載する。また、各基準点の位置は、ユーザのマニュアル操作により入力装置４から入力されたものとする。また、式（７）に示す評価関数Ｓ（μ）においては、第１分割画像と第２分割画像のペア（Ａ０、Ｂ０）、（Ａ１、Ｂ１）、…、（Ａｋ、Ｂｋ）の類似度（分割画像類似度）を、各ペア画像の相互相関情報に先述の蓋然性情報により重みづけた項目により定義している。
【００５４】
画像生成部１４は、画像変形量推定部１２により決定された変換関数により第１画像Ｖ１を変換した画像Ｖ１Ａを生成する。
【００５５】
表示制御部１５は、画像生成部１４により生成された画像Ｖ１Ａおよび第２画像Ｖ２を比較可能にディスプレイ３に表示させる。また、ユーザの入力などの必要に応じて、取得した第１画像Ｖ１および第２画像Ｖ２および／または本実施形態の画像処理プログラムの実行過程で生成した各画像をディスプレイ３に表示させる。
【００５６】
図２は、第２の実施形態の画像処理装置の処理の流れを示すフローチャートである。図２に従って、処理の流れを説明する。まず、画像取得部１１は、被写体を撮像して得られた第１画像（第１画像データ）Ｖ１および第２画像（第２画像データ）Ｖ２を取得する（Ｓ０１）。
【００５７】
次に、画像変形量推定部１２は、予め設定された分割条件に従って、第２画像を分割した第２分割画像を設定する。具体的には、ユーザの入力装置の入力に基づいて第２画像Ｖ２の基準点Ｐ２を特定し、基準点Ｐ２と第２画像における各制御点Ｘとの距離をそれぞれ算出する。そして、基準点Ｐ２からの距離に応じて、所定の距離の範囲ｄ０、ｄ１、ｄ２、…、ｄｋごとに、第２画像Ｖ２を複数の第２分割画像として設定し、各分割画像を特定する情報をメモリに記憶する（Ｓ０２）。
【００５８】
次に、画像変形量推定部１２は、画像変形量μを設定して（Ｓ０３）、変換関数ｇを用いて第１画像Ｖ１を変形する（Ｓ０４）。
【００５９】
また、画像変形量推定部１２は、予め設定された分割条件に従って、第１画像を分割した第１分割画像を設定する。具体的には、ユーザの入力装置の入力に基づいて変形した第１画像Ｖ１ａの基準点Ｐ１を特定し、基準点Ｐ１と変形した第１画像Ｖ１ａにおける各制御点Ｘとの距離をそれぞれ算出する（Ｓ０５）。そして、基準点Ｐ１からの距離に応じて、所定の距離の範囲ｄ０、ｄ１、ｄ２、…、ｄｋごとに、変形した第１画像Ｖ１ａを複数の第１分割画像として設定し、各分割画像を特定する情報をメモリに記憶する。
【００６０】
次いで、画像変形量推定部１２は、評価値として式（３）であらわされる評価関数Ｓ（μ）の値の変化量｜Ｓ（μ）−Ｓ（μ−Δμ）｜を算出する（Ｓ０６）。そして、算出された関数Ｓ（μ）の変化量｜Ｓ（μ）−Ｓ（μ−Δμ）｜が所定のしきい値より大きい場合には（Ｓ０７のＮ）、画像変形量μを所定量Δμ増加させたμ＋Δμを新しい画像変形量μとして設定し（Ｓ１０）、再度Ｓ０３からＳ０６の処理を繰り返す。なお、画像変形量推定部１２は、評価値として、評価関数Ｓ（μ）の変化量｜Ｓ（μ）−Ｓ（μ−Δμ）｜に換えて、式（７）であらわされる評価関数Ｓ（μ）の偏導関数の絶対値｜▽Ｓ（μ）｜を用いてもよい。
【００６１】
一方、評価関数Ｓ（μ）の変化量｜Ｓ（μ）−Ｓ（μ−Δμ）｜が所定のしきい値以下である場合には（Ｓ０７のＹ）、画像変形量推定部１２は、この場合の画像変形量μを評価関数Ｓ（μ）が最大となる場合の画像変形量μとして取得し、この画像変形量μに基づいて、第１画像を変形する変換関数を決定する。画像生成部１４は、決定された変換関数を用いて第１画像を変換して再構成し、画像Ｖ１Ａを生成する（Ｓ０８）。そして、表示制御部１５は、生成された画像Ｖ１Ａと第２画像Ｖ２を比較可能に並列表示する（Ｓ０９）。
【００６２】
第１の実施形態によれば、変形した第１画像Ｖ１ａと第２画像Ｖ２をそれぞれ対応する複数の分割画像に分割し、評価関数Ｓ（μ）が複数の第１分割画像と第２分割画像のペアごとに第１および第２分割画像間の画素値の分布の相関性をそれぞれ定義した複数の分割画像類似度に基づいて類似度を評価することにより、互いに異なる分割領域に位置する複数の被写体を表す画素値間で相関性が評価されることが抑制され、上記第１の問題の発生を低減することができる。結果として、さらに精度良く類似度を評価することができるため、より好適に第１画像を第２画像と一致するように変換した画像を生成することができる。
【００６３】
また、分割条件を基準点からの距離という１つの空間パラメータに応じて変形した第１画像および第２画像を分割するものとしたため、基準点Ｐ１、Ｐ２からの距離に応じて、所定の距離の範囲ｄ０、ｄ１、ｄ２、…、ｄｋごとに、対応する各第１分割画像および各第２分割画像を特定する情報をメモリに記憶することにより、変形した第１画像および第２画像の類似度算出のための計算負荷を必要以上に増加させることがない。また、第１の空間パラメータを、基準となる点からの距離としたため、第１画像および第２画像において基準となる図形の設定が容易である。
【００６４】
第１の実施形態の分割条件のさらなる変形例として、距離に換えて角度により変形した第１画像Ｖ１ａおよび第２画像Ｖ２を分割してもよい。このように、分割条件を基準となる図形に対する角度θに応じて所定角度範囲ごとに、変形した第１画像Ｖ１ａおよび第２画像Ｖ２をそれぞれ分割することが考えられる。なお、この場合、式（３）における距離ｄを角度θに置換したものを評価関数として用いることができる。
【００６５】
また、第１分割画像および第２分割画像を分割する単位を小さくするほど、各分割画像に含まれる被写体が少なくなる可能性が高いため、第１の問題をより抑制でき、変形した第１画像および第２画像の類似度を精度よく評価できる。
【００６６】
図３Ｂに、第２の実施形態における第１および第２画像の分割方法の他の例として、変形した第１画像Ｖ１ａおよび第２画像Ｖ２をｘｙ平面（所定の平面）において、ｘ軸（所定の軸）からの角度θに基づいて分割した例を示す。図３Ｂにおける変形した第１画像Ｖ１ａおよび第２画像Ｖ２は、共に３次元画像の原点を通るｘｙ平面を表しており、変形した第１画像Ｖ１ａの各制御点Ｘを説明のために、ｘｙ平面上でのみ変位させて表している。なお、角度θは、両画像Ｖ１ａ、Ｖ２をそれぞれ円柱座標系で表した場合の角度θとして取得できる。
【００６７】
図３Ｂでは、０≦θ＜９０°、９０°≦θ＜１８０°、１８０°≦θ＜２７０°、２７０°≦θ＜３６０°の４つの範囲に、角度に応じて変形した第１画像Ｖ１ａを第１分割画像Ａ０、Ａ１、Ａ２、Ａ３に分割し、第２画像Ｖ２を第２分割画像Ｂ０、Ｂ１、Ｂ２、Ｂ３に分割している。θの範囲を任意に設定することにより、分割条件を基準となる図形に対する任意の角度θに応じて、０≦θ＜θ０、θ０≦θ＜θ１、θ１≦θ＜θ２、…、θｋ−１≦θ＜θｋのｋ個に変形した第１画像Ｖ１ａおよび第２画像Ｖ２を分割するものとしてもよい。
【００６８】
上記のように、分割条件を基準図形からの角度という１つの空間パラメータに応じて変形した第１画像および第２画像を分割するものとし、基準図形（基準となるｘ軸）からの角度に応じて、所定の角度の範囲θ０、θ１、θ２、…、θｋごとに、対応する各第１分割画像および各第２分割画像を特定する情報をメモリに記憶した場合には、変形した第１画像および第２画像の類似度算出のための計算負荷を必要以上に増加させることがない。
【００６９】
また、ｘ軸に対する任意の角度（０≦θ＜θ０、θ０≦θ＜θ１、θ１≦θ＜θ２、…、θｋ−１≦θ＜θｋ）に基づいて、変形した第１画像Ｖ１ａおよび第２画像Ｖ２をそれぞれ分割し、この分割した画像をさらに、ｚ軸との任意の角度（０≦β＜β０、β０≦β＜β１、β１≦β＜β２、…、βｍ−１≦β＜βｍ）（ｍは０以上の整数）に基づいて、さらに分割してもよい。このように２つの軸に対するそれぞれの角度に基づいて変形した第１画像Ｖ１ａおよび第２画像Ｖ２を分割した場合には、変形した第１画像Ｖ１ａおよび第２画像Ｖ２を３次元的に近接した領域ごとに分割できるため、第１の問題をより好適に抑制でき、変形した第１画像および第２画像の類似度を精度よく評価できる。
【００７０】
また、変形した第１画像Ｖ１ａを、距離に基づいて分割するとともに、角度にも基づいてさらに分割して第１分割画像を設定してもよい。
【００７１】
また、本実施形態においては、異なるモダリティにより同一の被写体を撮影して得られた２つの画像の相関性の尺度を表す項目として、式（２）で表す相互相関情報量を用いるが、２つの画像の相関性の尺度を表す項目は、相互相関情報量を含むものであれば、さらに周知の種々の変形を加えたものであってもよい。また、評価関数は、２つの画像の相関性の尺度を表す項目が蓋然性情報により重み付けされたものであれば、変形の滑らかさの限界条件を規定するなどの他の項目を更に加えたものであってもよく、周知の種々の変形を行ったものであってもよい。
【００７２】
次いで、本発明の第２の実施形態について以下に説明する。図５は、第２の実施形態の機能ブロック図を示す概略図である。図５に示すように、第２の実施形態では、画像変形量推定部１２が、第１画像および第２画像のそれぞれの対応する画素値の組合せの確からしさを蓋然性情報として取得する蓋然性情報取得部１３を備え、評価関数における類似度の尺度を表す項目に対しこの蓋然性情報による重み付けを行う点が第１の実施形態と異なる。
【００７３】
第２の実施形態においては、評価関数Ｓ（μ）を用いた第１画像を変形した画像と第２画像の類似度の評価以外の処理については、第１の実施形態と同じであり、各機能ブロックの機能も共通している。以下、第１の実施形態と異なる点を中心に説明し、第１の実施形態と同じ点については説明を省略する。
【００７４】
本発明者は、先述の式（２）のような従来の評価関数では、各モダリティの種類などにより得られる画素値の組合せの確からしさを全く評価していないため、第１画像を変形した画像の画素値に対応する第２画像の画素値がある程度決まった範囲となることが予想される場合などであっても、類似度算出の評価に反映させることができないという問題に着目した（第２の問題）。
【００７５】
本発明者は上記第２の問題に鑑み、評価関数に、第１画像を変形した画像と第２画像の対応する画素値の組合せの確からしさを表す蓋然性情報を表す項目を組み込むことが有効であることを見出した。そして、式（３）の評価関数の変形第１画像と第２画像の類似性の尺度を表す項目に式（４）のように、蓋然性情報ＰＬ（ｍ；μ｜ｆ）を適用した。つまり、第１画像を変形した画像と第２画像の類似性の尺度を表す項目を、第１画像を変形した画像と第２画像の画素値の組合せごとに、両画像の画素値の組合せの確からしさに応じて、両画像の組合せの確からしさを表す蓋然性情報ＰＬ（ｍ；μ｜ｆ）で重み付けするものとした。
【００７６】
式（４）に示すように、本実施形態における評価関数においては、第１画像を変形した画像と第２画像の類似性の尺度を表す項目を、第１画像を変形した画像と第２画像の画素値の組合せの確からしさに応じて、両画像の組合せの確からしさを表す蓋然性情報ＰＬ（ｍ；μ｜ｆ）で重み付けしている。なお、図３Ａに示すように、本第２の実施形態における分割条件は、変形した第１画像Ｖ１ａおよび第２画像Ｖ２のそれぞれの基準点Ｐ１、Ｐ２を設定し、各基準点Ｐ１、Ｐ２から距離に応じてＶ１ａおよび第２画像Ｖ２をそれぞれ分割することを規定する。
【数４】
0005955199.tif 000005
【００７７】
蓋然性情報ＰＬ（ｍ；μ｜ｆ）は、所定の被写体を第１モダリティにより撮影して得られた画像と同じ種類の所定の被写体を第２モダリティにより撮影された画像の対応する画素値の組合せの確からしさを規定するものであればいかなる情報でもよい。ここでいう所定の被写体は、第１画像Ｖ１および第２画像Ｖ２に表す被写体と同じであってもよく、異なっていてもよい。
【００７８】
本実施形態では、第１モダリティにより撮影して得られた画像と第２モダリティにより撮影して得られた画像の対応する画素値の組合せの確からしさを、第１画像の各画素の画素値がｆとなるという事象をもとに第２画像の対応する各画素の画素値がｍ；μとなる事象が生じる条件付き確率として規定する。なお、下記式（５）に表されるように、条件付き確率は、同時確率分布関数ｐＬ（ｆ，ｍ；μ）と確率変数ｆの周辺確率分布関数ｐＬ（ｆ）により算出することができる。
【数５】
0005955199.tif 000006
【００７９】
本実施形態においては、蓋然性情報取得部１３は、第１画像を画像変形量μだけ変形した画像の全ての制御点の画素値ｆ（ｆ∈Ｆ）と第２画像の対応する各制御点の画素値ｍ（ｍ∈Ｍ）を取得する。そして、周知の方法により、各制御点の画素値ｆ（ｆ∈Ｆ）ごとに対応する各制御点の各画素値ｍ（ｍ∈Ｍ）の分布をそれぞれ調べることにより同時確率分布関数ｐＬ（ｆ，ｍ；μ）を取得し、各制御点の画素値ｆについて分布を調べることにより周辺確率分布関数ｐＬ（ｆ）を取得する。そして、式（５）に基づいて、蓋然性情報ＰＬ（ｍ；μ｜ｆ）を取得する。
【００８０】
なお、上記の条件付き確率の算出には、所定の被写体を第１モダリティにより撮影して得られた第１参照画像および同じ所定の被写体を第２モダリティにより撮影して得られた第２参照画像を用いることができる。なお、所定の被写体は必ずしも第１画像および第２画像に表される被写体と同じ種類の被写体である必要はないが、蓋然性情報をより精度よく算出するために、第１参照画像及び第２参照画像は、第１画像および第２画像に表す被写体と同じ種類の被写体を撮影したものであることが好ましい。また、第１参照画像と第２参照画像のペアは、１つであってもよく、複数であってもよい。第１参照画像と第２参照画像のペアを複数用いて、条件付き確率を算出した場合には、画素値の組合せの確からしさをより精度よく推定できると考えられるため好ましい。
【００８１】
また、蓋然性情報ＰＬ（ｍ；μ｜ｆ）の算出または取得は、評価関数Ｓを最大とする画像変形量を算出する処理以前に行うものであれば任意のタイミングで行ってよい。
【００８２】
図６は、第２の実施形態の画像処理装置の処理の流れを示すフローチャートである。図６に従って、処理の流れを説明する。まず、画像取得部１１は、被写体を撮像して得られた第１画像（第１画像データ）Ｖ１および第２画像（第２画像データ）Ｖ２を取得する（Ｓ１１）。
【００８３】
次に、画像変形量推定部１２は、予め設定された分割条件に従って、第２画像を分割した第２分割画像を設定する。具体的には、ユーザの入力装置の入力に基づいて第２画像Ｖ２の基準点Ｐ２を特定し、基準点Ｐ２と第２画像における各制御点Ｘとの距離をそれぞれ算出する。そして、基準点Ｐ２からの距離に応じて、所定の距離の範囲ｄ０、ｄ１、ｄ２、…、ｄｋごとに、第２画像Ｖ２を複数の第２分割画像として設定し、各分割画像を特定する情報をメモリに記憶する（Ｓ１２）。
【００８４】
次に、画像変形量推定部１２は、画像変形量μを設定して（Ｓ１３）、変換関数ｇを用いて第１画像Ｖ１を変形する（Ｓ１４）。
【００８５】
また、画像変形量推定部１２は、予め設定された分割条件に従って、第１画像を分割した第１分割画像を設定する。具体的には、ユーザの入力装置の入力に基づいて変形した第１画像Ｖ１ａの基準点Ｐ１を特定し、基準点Ｐ１と変形した第１画像Ｖ１ａにおける各制御点Ｘとの距離をそれぞれ算出する（Ｓ１５）。そして、基準点Ｐ１からの距離に応じて、所定の距離の範囲ｄ０、ｄ１、ｄ２、…、ｄｋごとに、変形した第１画像Ｖ１ａを複数の第１分割画像として設定し、各分割画像を特定する情報をメモリに記憶する。
【００８６】
次いで、画像変形量推定部１２は、変形した第１画像Ｖ１ａと第２画像Ｖ２各画素値の各画素値を取得して同時確率分布ｐ（ｆ,ｍ；μ）や周辺確率分布ｐ（ｆ）、ｐ（ｍ；μ）を算出し、蓋然性情報取得部１３は、算出された同時確率分布ｐ（ｆ,ｍ；μ）、周辺確率分布ｐ（ｆ）、ｐ（ｍ；μ）および蓋然性情報ＰＬ（ｍ；μ｜ｆ）に基づいて、先述のように条件付き確率を算出することにより蓋然性情報ＰＬ（ｍ；μ｜ｆ）を取得する（Ｓ１６）。そして、画像変形量推定部１２は、式（７）であらわされる評価関数Ｓ（μ）の変化量｜Ｓ（μ）−Ｓ（μ−Δμ）｜を算出する（Ｓ１７）。そして、算出された評価関数Ｓ（μ）の変化量｜Ｓ（μ）−Ｓ（μ−Δμ）｜が所定のしきい値より大きい場合には（Ｓ１８のＮ）、画像変形量μを所定量Δμ増加させたμ＋Δμを新しい画像変形量μとして設定し（Ｓ２１）、再度Ｓ１３からＳ１８の処理を繰り返す。なお、画像変形量推定部１２は、評価値として、評価関数Ｓの変化量｜Ｓ（μ）−Ｓ（μ−Δμ）｜に換えて、評価関数Ｓ（μ）の偏導関数の絶対値｜▽Ｓ（μ）｜を用いてもよい。
【００８７】
一方、評価値の変化量｜Ｓ（μ）−Ｓ（μ−Δμ）｜が所定のしきい値以下である場合には（Ｓ１８のＹ）、画像変形量推定部１２は、この場合の画像変形量μを評価関数Ｓ（μ）が最大となる場合の画像変形量μとして取得し、この画像変形量μに基づいて、第１画像を変形する変換関数を決定する。画像生成部１４は、決定された変換関数を用いて第１画像を変換して再構成し、画像Ｖ１Ａを生成する（Ｓ１９）。そして、表示制御部１５は、生成された画像Ｖ１Ａと第２画像Ｖ２を比較可能に並列表示する（Ｓ２０）。
【００８８】
上記本実施形態によれば、評価関数において、変形第１画像と第２画像の各画素値の類似性を表す項目を、変形第１画像および第２画像の各制御点の画素値の組合せの確からしさを表す蓋然性情報により重み付けしたため、第１モダリティにより撮影された画像の各制御点の画素値と第２モダリティにより撮影された画像の画素値とが、実際にはあり得ないような画素値の組合せに対して、蓋然性条件に応じて類似度を低くするように重み付けすることができる。このことにより、従来よりも正確に類似度を評価することができ、結果として好適に第１画像を第２画像に一致させた画像を生成することができる。
【００８９】
また、評価関数の第１モダリティにより撮影された画像と第２モダリティにより撮影された画像の各画素値の類似性を表す項目を、第１モダリティにより撮影された画像の各制御点の画素値がｆとなるという事象をもとに第２モダリティにより撮影された画像の対応する制御点の画素値がｍ；μとなるという事象が生じる条件付き確率により重み付けしたものとしたため、両画像の画素値の組合せの確からしさに応じてより精度良く両画像の類似度を算出することができる。
【００９０】
また、第１画像及び第２画像の各画素値に基づいて条件付き確率を算出して蓋然性情報として用いる場合には、両画像の画素値分布の情報さえあれば蓋然性情報を設定することができるため、蓋然性情報の設定のための情報収集や解析が不要であり、評価関数への適用が容易である。また、第１参照画像及び第２参照画像の各画素値に基づいて条件付き確率を算出して蓋然性情報として用いる場合には、両画像の画素値分布の情報さえあれば蓋然性情報を設定することができるため、評価関数への適用が容易である。また、同時確率分布ｐ（ｆ,ｍ；μ）や周辺確率分布ｐ（ｆ）に基づいて、相互相関情報量および条件付き確率の両方を算出できるため、計算効率がよい。
【００９１】
また、第２の実施形態の変形例として、蓋然性情報が、被写体の種類ごとに、所定の種類の被写体を前記第１のモダリティで撮影した場合に得られる第１の画素値範囲と、所定の種類の被写体を前記第２のモダリティで撮影した場合に得られる第２の画素値範囲をそれぞれ対応付けたものであり、評価関数が、蓋然性情報に基づいて、変形した第１画像と第２画像の各画素値が、第１の画素値範囲と第１の画素値範囲に対応付けられた第２の画素値範囲の両方を満たすものでない場合には、変形した第１画像と第２画像の類似度を低くするように重み付けされたものであってもよい。なお、本明細書において、類似度を低くするように重み付けするとは、類似度が大きいほど第１変形画像と第２画像が類似していると評価する評価関数を用いる場合には、類似度を小さくするように重み付けすることを意味し、類似度が小さいほど第１変形画像と第２画像が類似していると評価する評価関数を用いる場合には、類似度を大きくするように重み付けすることを意味する。
【００９２】
一例として、上記第２の実施形態の変形例の場合の評価関数を、下記条件式（６）のように定義できる。以下の例では、第２モダリティにより撮影された画像の画素値ｆ（ｆ∈Ｆ）ごとに、第１モダリティにより撮影された画像の取り得る画素値の範囲ｍｍｉｎｆ≦ｍ（μ）≦ｍｍａｘｆを対応づけて記憶する。そして、第１画像の画素値ｆ（第１の画素値範囲）に対して第２画像の画素値ｍ（μ）が対応づけられた範囲ｍｍｉｎｆ≦ｍ（μ）≦ｍｍａｘｆ（第２の画素値範囲）に属する場合には、評価関数の第１画像と第２画像の相関性の尺度を示す項目に１により重み付けし、それ以外の場合には、実際にはあり得ない画素値の組合せと判断して、評価関数の第１画像と第２画像の相関性の尺度を示す項目に０により類似度が低くなるように重み付けする。
【数６】
0005955199.tif 000007
【００９３】
上記の第２の実施形態の変形例において、第１画像と第２画像の画素値の対応する画素値の範囲を、予め撮影した試験画像を解析した既知の情報により規定することができる。例えば、第１のモダリティにより撮影された画像における、水、空気、各解剖学的構造物（または解剖学的構造物を構成する各要素）などの被写体の種類ごとに、各種類の被写体の示す画素値の範囲を、第２のモダリティにより撮影された画像における、水、空気、各解剖学的構造物などのそれぞれ対応する対象物の示す画素値の範囲に対応づけて蓋然性情報として記憶する。そして、画像変形量推定部１２が、蓋然性情報に対応づけられた画素値の範囲の組合せのいずれにも該当しない画素値の組合せに対して、類似度を低くすることが考えられる。この場合には、蓋然性情報において、同一の被写体を表していることが既知である２つの画像の画素値の情報に基づいて対応する画素値範囲を設定しているため、より正確に蓋然性情報を定義することができる。
【００９４】
また、第１画像と第２画像の画素値の対応する画素値の範囲を、各モダリティの撮像原理に基づいて推定される理論的な画素値から決定してもよい。例えば、被写体となりうる対象物の種類毎に、各モダリティによりその対象物を撮影した場合の画素値を推定して対応付けて記憶しておく。例えば、水、空気、各解剖学的構造物（または解剖学的構造物を構成する各要素）などの被写体の種類ごとに、各被写体を撮影したときに得られる第１のモダリティによる画像の画素値の範囲と、第２のモダリティによる画像における画素値の範囲を、撮影原理に基づいて推定する。そして、対象物ごとに、推定された第１のモダリティによる画像の画素値の範囲および第２のモダリティによる画像の画素値の範囲を対応付けて記憶する。そして、蓋然性情報に対応づけられた画素値の範囲の組合せのいずれにも該当しない画素値の組合せに対して、類似度の重み付けを低くするようにすることが考えられる。この場合には、蓋然性情報において、撮像原理から推定される２つの画像の画素値の情報に基づいて対応する画素値範囲を設定しているため、より正確に蓋然性情報を定義することができる。
【００９５】
また、上述の各蓋然性情報を、さらに重み付けを調整するようにスケーリング（正規化）したものとしてもよい。この場合には、条件付き確率が非常に小さい値になる場合など、第１および第２参照画像の相関性を表す項目に対する重み付けが大きくなりすぎない（または小さくなりすぎない）ように調整することができる。また、重み付けの調整をする周知の種々の方法を適用することができる。
【００９６】
例えば、条件付き確率を蓋然性情報に用いる場合、画素値の分布の形状による影響を受けにくくなるように、蓋然性情報を正規化することが好ましい。条件付き確率の算出に用いる要素である、画素値を確率変数とした確率密度関数は、画素値の分布が平坦である場合には値が小さくなり、画素値の分布が集中する部分を有する急峻なものである場合には値が大きくなるように、画素値の分布の形状に応じて値が変動する性質を有するからである。このため、まず、第２モダリティにより撮影された画像の画素値ｆ（ｆ∈Ｆ）ごとに、第１モダリティにより撮影された画像の取り得る画素値の範囲ｍｍｉｎｆ≦ｍ（μ）≦ｍｍａｘｆを対応づけて記憶する。そして、第２モダリティにより撮影された画像の画素値ｆ（ｆ∈Ｆ）ごとに、この画素値ｆに対応付けられた第１参照画像の画素値分布を一様分布と近似し、第１参照画像の確率密度関数ＰＵ（ｍ）（＝ｍｍａｘｆ−ｍｍｉｎｆ）を算出する。そして、下記式（７）に示すように、蓋然性情報を確率密度関数ＰＵ（ｍ）の逆数Ａｍによりさらに重み付けしたものとする。
【数７】
0005955199.tif 000008
【００９７】
この場合には、第２モダリティにより撮影された画像の画素値ｆに対応づけられた第１モダリティにより撮影された画像の画素値ｍｍａｘｆ＜ｍ（μ）＜ｍｍｉｎｆの分布形状の影響を抑制することができ、第１モダリティにより撮影された画像の画素値と第２モダリティに寄り撮影された画像の画素値の組合せの確かしさをより好適に反映して、より適切に重み付けを行うことができる。
【００９８】
上記の各実施形態はあくまでも例示であり、上記のすべての説明が本発明の技術的範囲を限定的に解釈するために利用されるべきものではない。
【００９９】
この他、上記の実施形態におけるシステム構成、ハードウェア構成、処理フロー、モジュール構成、ユーザインターフェースや具体的処理内容等に対して、本発明の趣旨から逸脱しない範囲で様々な改変を行ったものも、本発明の技術的範囲に含まれる。
【０１００】
また、画像処理装置１は、複数台のコンピュータにより、手段としての機能を分担する構成としてもよい。また、入力装置、ディスプレイ等、システムを構成する装置としては、公知のあらゆる装置を採用することができる。
【符号の説明】
【０１０１】
１  画像処理装置、  ２  ストレージ、  ３  ディスプレイ（表示装置）
４  入力装置
１１  画像取得部
１２  画像変形量推定部
１３  蓋然性取得部
１４  画像生成部
１５  表示制御部
Ｖ１  ＭＲ画像（第１画像）
Ｖ２  ＣＴ画像（第２画像）
【図１】
0005955199.tif 000009
【図２】
0005955199.tif 000010
【図３Ａ】
0005955199.tif 000011
【図３Ｂ】
0005955199.tif 000012
【図５】
0005955199.tif 000013
【図６】
0005955199.tif 000014
【図４Ａ】
0005955199.tif 000015
【図４Ｂ】
0005955199.tif 000016