JP6467221B2 - 画像処理装置および方法 - Google Patents

Description

本発明は、画像間の位置合わせの誤差を推定する画像処理装置および方法に関する。

医用画像（被検体内部の情報を表す３次元断層画像）を用いた画像診断において、医師は、複数の撮像装置（モダリティ）、異なる体位、時刻、撮像パラメータ等で撮像した画像を対比しながら診断を行う。しかし、画像間で被検体の姿勢や形状が異なると、画像間における病変部の同定や対比を行うことが困難になる。そこで、複数の画像間の位置合わせを行うことにより、一方の画像に姿勢の変換や変形を施して他方と同一にした画像を生成することが試みられている。また、このような位置合わせを行うことにより、一方の画像上で注目した点の、他方の画像上における対応点の位置を算出して提示することが可能となる。その結果、医師は、複数の画像間における病変部の同定や対比を容易に行うことが可能となる。医療以外の分野においても、物体の内部状態を検査する目的において同様の作業が実施される場合がある。

しかしながら、一般に、位置合わせの結果は誤差を含むのであり、必ずしも正確ではないため、医師は、位置合わせの結果をどの程度信頼して良いのか判断できないという課題がある。これに対して特許文献１は、変形位置合わせの結果として推定した変形パラメータの不安定さ（解の曖昧さ）に基づいて誤差を推定する手法を開示している。特許文献１の手法では、不安定なパラメータを意図的に変動させた際に注目点の対応点の推定位置が変動する範囲に基づいて、当該位置の誤差を推定している。

特開２０１３−１９８７２２号公報

しかし、特許文献１に記載の誤差の推定方法では、推定パラメータの不安定さに基づいて任意の位置における偶然誤差（解の曖昧さ）を推定することはできるものの、それ以外の要因を考慮できないという課題があった。

本発明は、上記の課題に鑑みてなされたものであり、対応点以外の任意の位置における位置合わせ誤差をより正確に推定可能とすることを目的とする。

上記の目的を達成するための本発明の一様態による画像処理装置は、以下の構成を備える。すなわち、
第１の画像および第２の画像に対する複数の対応情報を取得する取得手段と、
前記複数の対応情報を用いて前記第１の画像と前記第２の画像を位置合わせし、位置合わせされた前記第１の画像および前記第２の画像に対する前記複数の対応情報によって定義された複数の位置に生じる位置合わせ誤差を算出する算出手段と、
前記複数の位置とは異なる更なる位置に生じる位置合わせ誤差を、前記位置に生じる位置合わせ誤差に基づいて推定する推定手段と、を備える。

本発明によれば、対応点以外の任意の位置における位置合わせ誤差をより正確に推定することが可能となる。

第１実施形態に係る画像診断システムおよび画像処理装置の機能構成を示す図。 第１実施形態に係る画像処理装置の処理手順を示すフローチャート。 第１実施形態に係る誤差推定部１７０の処理手順を示すフローチャート。 第１実施形態における推定誤差の算出処理を説明するための図。 第１実施形態で算出した推定誤差の例を示す模式図。 第１実施形態における推定誤差表示の例を示す模式図。 変形例２における推定誤差表示の例を示す模式図。 第２実施形態における推定誤差の算出処理を説明するための図。 第２実施形態における代表誤差増加率の取得処理を示すフローチャート。

以下、添付の図面を参照しながら、本発明の実施の形態について詳説する。ただし、本発明の範囲は図示例に限定されるものではない。

（第１実施形態）
第１実施形態に係る画像処理装置は、複数の３次元断層画像間の対応情報に基づいて位置合わせ処理（変形位置合わせ）を行い、一方の３次元断層画像を他方の３次元断層画像における位置や形状と一致するように変形させた変形画像を生成する。その際に、画像処理装置は、位置合わせ処理の誤差に基づいて、生成した変形画像の夫々の位置における位置合わせの推定誤差を求め、その分布を変形画像と対応付けて表示する。ここで、対応情報とは、２つの画像間で対応する点や線や面の情報のことである。本実施形態では、３次元断層画像間の複数の対応点を対応情報として用いて注目点における誤差を推定する場合について説明するが、対応する線や対応する面を用いることもできる。以下、第１実施形態に係る画像処理装置について説明する。なお、本実施形態では、位置合わせ処理の誤差、推定誤差には、各軸方向の誤差を成分とするベクトルが用いられる。

＜１．画像診断システム１の構成＞
図１は、本実施形態に係る画像診断システム１の構成例を示す。図１に示すように、第１実施形態における画像診断システム１は、画像処理装置１０と、データサーバ２０と、ディスプレイ３０とを有している。画像処理装置１０は、画像取得部１１０、対応情報取得部１２０、対応情報選択部１３０、位置合わせ部１４０、誤差算出部１５０、変形画像生成部１６０、誤差推定部１７０および表示制御部１８０を備える。

画像取得部１１０は、位置合わせの対象である被検体の複数の３次元断層画像を、データサーバ２０から取得する。以下、位置合わせの対象となる３次元断層画像が第１の画像と第２の画像であるとして説明する。なお、位置合わせの対象となる画像は、３次元断層画像に限られるものではなく、２次元断層画像等であってもかまわない。

対応情報取得部１２０は、画像取得部１１０により取得された第１の画像および第２の画像における複数の対応点の情報を、画像間の対応情報（以下、第１の対応情報）として取得する。ここで、第１の対応情報における夫々の対応点の情報とは、第１の画像と第２の画像の夫々において被検体上の同一の位置（対応している地点）を表す座標の対である。対応情報選択部１３０は、第１の対応情報（複数の対応点）の中から誤差算出の対象となる対応点を選択し、当該選択された対応点を除いた残りの対応点を、後述する第２の位置合わせに用いるための対応情報（以下、第２の対応情報）として選択する。

位置合わせ部１４０は、対応情報取得部１２０により取得された第１の対応情報に基づいて、第１の画像と第２の画像の変形推定処理（第１の位置合わせ処理）を実行する。また、位置合わせ部１４０は、第２の対応情報に基づいて、第１の画像と第２の画像の変形推定処理（第２の位置合わせ処理）を実行する。

誤差算出部１５０は、第１の位置合わせ処理の結果に対して、複数の対応点の夫々の位置における位置の残差ベクトルを算出する。なお、この残差ベクトルは、ＦＲＥ（Fiducial Registration Error）であり、以下の説明ではこれを第１の誤差と呼ぶ。また、誤差算出部１５０は、第２の位置合わせ処理の結果に対して、対応情報選択部１３０が選択した対応点（すなわち、第２の位置合わせ処理には用いていない対応点）の位置における位置の誤差ベクトルを算出する。なお、この誤差ベクトルは、ＴＲＥ（Target Registration Error）であり、以下の説明ではこれを第２の誤差と呼ぶ。なお、対応点ではなく対応する線や対応する面を用いる場合には、位置合わせ後の対応する線と線や面と面の間の距離が最大となる方向ベクトルを第１の誤差や第２の誤差とすればよい。ここで、線と線や面と面の間の距離の算出には公知の方法を用いることができる。変形画像生成部１６０は、第１の位置合わせ処理の結果に基づいて、第２の画像と一致するように第１の画像に座標変換を施した新たな３次元画像（変形画像）を生成する。

誤差推定部１７０は、誤差算出部１５０が算出した各対応点位置における第１の誤差および第２の誤差に基づいて、変形画像上の各点における位置合わせの推定誤差を計算し、その分布を示す推定誤差画像を生成する。たとえば、変形画像における各ボクセルまたは各ピクセルの値を、推定された任意の位置の位置合わせ誤差に基づいて決定することで推定誤差画像が生成される。なお、本実施形態における推定誤差画像は、各ボクセル値が、その位置における位置合わせの推定誤差を示すボリューム画像である。なお、推定誤差画像が２次元画像として生成される場合は、各ピクセル値がその位置における位置合わせの推定誤差を示す画像となる。表示制御部１８０は、第２の画像の断面画像と、変形画像（変形画像生成部１６０による変形後の第１の画像）の断面画像（対応断面画像）とを、ディスプレイ３０に表示する制御を行う。また、表示制御部１８０は、誤差推定部１７０が推定した誤差の情報（たとえば、推定誤差画像から取得される）を断面画像に関連付けてディスプレイ３０に表示する制御を行う。

データサーバ２０は、位置合わせの対象である被検体の複数の３次元断層画像を保持している。なお、夫々の３次元断層画像は、付帯情報として、画像のサイズ、解像度、モダリティの種類、撮像情報（撮像パラメータ、撮像部位、体位など）、症例情報（患者情報、検査情報、診断情報、臓器領域情報、注目領域情報など）、画像間の対応情報等を含んでいる。これらの付帯情報は、必要に応じて、画像と共に画像処理装置１０へと送信される。ディスプレイ３０は、画像処理装置１０（表示制御部１８０）の制御下で各種の表示を行う。たとえば、図６を参照して後述するような断層画像の表示が行われる。

＜２．画像処理装置１０が実施する処理＞
次に、図２のフローチャートを参照して、第１実施形態に係る画像処理装置１０が実施する処理の手順について説明する。なお、以下の説明では、同一の被検体を異なる変形状態で撮像した２つの３次元断層画像を第１の画像と第２の画像として読み込み、第１の画像に変形処理を施して、第２の画像と位置や形状が合うように変形させた変形画像を生成する場合を例にとって説明する。

本実施形態の画像処理装置１０では、位置合わせされた画像上の任意の位置における位置合わせ誤差の推定を、
・被検体の第１の画像と第２の画像を位置合わせするための複数の対応情報を取得する取得処理、
・上記複数の対応情報から少なくとも一つの対応情報を除外した残りの対応情報を用いて第１の画像と第２の画像を位置合わせした場合の、その除外された対応情報の位置に生じる位置合わせ誤差を算出する算出処理、そして、
・複数の対応情報を用いて第１の画像と第２の画像を位置合わせをした場合の任意の位置における位置合わせ誤差を、上記算出された位置合わせ誤差に基づいて推定する推定処理、
により実現する。なお、取得処理は、対応情報取得部１２０により実行される（ステップＳ２１０）。算出処理は、対応情報選択部１３０、位置合わせ部１４０、誤差算出部１５０により実行される（ステップＳ２４０〜Ｓ２７０）。また、推定処理は、誤差推定部１７０により実行される（ステップＳ２８０）。以下、図２のフローチャートに示される各ステップを詳細に説明する。

まず、ステップＳ２００において、画像取得部１１０は、位置合わせの対象である被検体の複数の３次元断層画像（すなわち、第１の画像と第２の画像）を、データサーバ２０から取得する。画像取得部１１０は、取得した画像を、変形画像生成部１６０および表示制御部１８０に送信する。

次に、ステップＳ２１０において、対応情報取得部１２０が、第１の画像と第２の画像における第１の対応情報を取得し、取得した第１の対応情報を対応情報選択部１３０および位置合わせ部１４０に送信する。ステップＳ２１０は上述の取得処理の一例である。なお、対応情報取得部１２０による第１の対応情報の取得は、例えば、ユーザが目視で同定した画像間の対応点を画像処理装置１０に入力することによって実行される。すなわち、対応情報取得部１２０は、
・表示制御部１８０を介して、夫々の３次元断層画像の断面画像（第１の画像と第２の画像）をディスプレイ３０に表示させ、
・ユーザに、ディスプレイ３０に表示された第１の画像と第２の画像を比較させ、夫々の画像上において解剖学的に同じと見做した位置を、不図示のマウスのクリック等によって対応点として入力させ、
・ユーザにより入力された対応点を第１の対応情報として取得する。

ステップＳ２２０において、位置合わせ部１４０は、ステップＳ２１０で取得した第１の対応情報に基づいて、第１の画像と第２の画像について第１の位置合わせ処理（変形推定処理）を実行する。すなわち、位置合わせ部１４０は、第１の画像を変形した場合の第２の画像との間での対応点位置の残差（あるいは、該残差を含むコスト関数）が最小となるような変形情報（変形パラメータ）を推定する。このような変形推定処理は公知の技術なので、詳細な説明は省略する。そして、位置合わせ部１４０は、得られた推定結果（以下、第１の位置合わせ結果と呼ぶ）を、誤差算出部１５０および変形画像生成部１６０へと送信する。

ステップＳ２３０において、変形画像生成部１６０は、複数の対応情報を用いて第１の画像を第２の画像に位置合わせすることにより変形画像を生成する。すなわち、変形画像生成部１６０は、ステップＳ２２０で得た第１の位置合わせ結果に基づいて、第２の画像と一致するように第１の画像に座標変換を施し、新たな３次元画像（変形画像）を生成する。

次のステップＳ２４０〜Ｓ２７０は、上述の算出処理の一例である。まず、ステップＳ２４０において、対応情報選択部１３０は、ステップＳ２１０で取得した第１の対応情報（複数の対応点）の中から、誤差算出の対象となる対応点を順次（１つずつ）選択する。そして、当該対応点を除いた（Leave-one-outした）複数の対応点（残りの対応点）を、第２の対応情報として選択する。

ステップＳ２５０において、位置合わせ部１４０は、ステップＳ２４０で選択した第２の対応情報に基づいて、第１の画像と第２の画像について第２の位置合わせ処理（変形推定処理）を実行する。すなわち、ステップＳ２２０における第１の位置合わせ処理と同様に、第１の画像を変形した場合の第２の画像との間での対応点位置の残差が（あるいは、該残差を含むコスト関数が）最小となるような変形情報（変形パラメータ）を推定する。そして、得られた推定結果（以下、第２の位置合わせ結果と呼ぶ）を、誤差算出部１５０へと送信する。なお、本実施形態では、ステップＳ２２０における第１の位置合わせ処理（変形推定処理）とステップＳ２５０における第２の位置合わせ処理（変形推定処理）は同一のアルゴリズムであり、使用される対応情報（第１の対応情報と第２の対応情報）が異なる。

ステップＳ２６０において、誤差算出部１５０は、ステップＳ２２０で得た第１の位置合わせ結果（変形推定結果）に対して、ステップＳ２４０で選定された対応点の位置における残差ベクトル（ＦＲＥ）を算出し、これを該対応点の第１の誤差とする。また、誤差算出部１５０は、ステップＳ２５０で得た第２の位置合わせ結果（変形推定結果）に対して、ステップＳ２４０で選択された対応点の位置における誤差ベクトル（ＴＲＥ）を算出し、これを該対応点の第２の誤差とする。なお、第１の誤差および第２の誤差は、画像の軸方向成分ごとの誤差を成分とする誤差ベクトルであるが、該誤差ベクトルの大きさ（ノルム）を当該対応点の誤差（スカラー値）として算出する構成であってもよい。

ステップＳ２７０において、対応情報選択部１３０は、対応情報の選択処理を終了するか否かの判定を行う。すなわち、全ての対応点を選択し終えたか否かを判定する。選択処理を終了すると判定した場合には、対応情報の選択処理を終了して、ステップＳ２８０へと処理を進める。一方、選択処理を終了すると判定されなかった場合には、処理はステップＳ２４０へ戻り、残りの対応点（未選択の対応点）の誤差を算出するべくステップＳ２４０からＳ２６０の処理を順次実行する。

ステップＳ２８０は、上述の推定処理の一例である。ステップＳ２８０において、誤差推定部１７０は、ステップＳ２６０で各対応点について算出された第１の誤差および第２の誤差に基づいて、ステップＳ２３０で得た変形画像上の各点（任意の点）における位置合わせの推定誤差を計算する。そして、誤差推定部１７０は、その推定誤差の分布を示す推定誤差画像を生成し、表示制御部１８０へ送信する。ここで、ステップＳ２６０において、各対応点の誤差として誤差ベクトルが算出される場合には、誤差推定部１７０は、任意の位置における位置合わせの推定誤差もベクトルとして軸方向成分ごとに推定する。そして、各軸方向の推定誤差画像を生成する。さらに、各軸方向の推定誤差から推定誤差の大きさ（ノルム）を算出し、各位置における推定誤差の大きさを表す推定誤差画像を生成する。一方、ステップＳ２６０において、各対応点の誤差としてスカラー値（誤差の大きさ）が算出される場合には、誤差推定部１７０は、任意の位置における位置合わせの推定誤差もスカラー値として推定する。そして、推定誤差の大きさを表す推定誤差画像を生成する。ステップＳ２８０における推定誤差画像の生成については、図３のフローチャートを用いて後述する。

ステップＳ２９０において、表示制御部１８０は、ステップＳ２００で取得された第２の画像の断面画像とステップＳ２３０で生成された変形画像の断面画像を、ユーザの操作に応じて、ディスプレイ３０に表示する制御を行う。なお、表示される画像は、第２の画像と変形画像の互いに対応する断面画像である。また、表示制御部１８０は、ステップＳ２８０で取得した推定誤差画像から、変形画像の断面画像に対応する断面を切り出して、推定誤差マップとしてディスプレイ３０に表示する制御を行う。

次に、ステップＳ２８０における推定誤差の算出処理（推定処理）を、図３のフローチャートを参照して説明する。ステップＳ３００において、誤差推定部１７０は、ステップＳ２６０で得た夫々の対応点の位置における第２の誤差（ＴＲＥ）に基づいて、それらの値を補間して任意の位置ｘにおける第２の誤差の値を算出するための、第２の誤差の補間場Ｅ２（ｘ）を導出する。ここで、任意の位置ｘとは３次元位置または２次元位置である。補間場Ｅ２（ｘ）の導出とは、より具体的には、夫々の対応点の位置における第２の誤差に基づいて、対応点以外の位置における誤差の値を補間する補間モデルのパラメータを導出するものである。補間モデルとしては、例えば、放射基底関数を夫々の対応点位置に配置した関数群を用いることができる。なお、本実施形態では、放射基底関数としてＴＰＳ（Thin Plate Spline）関数を用いるものとする。

図４（ａ）は、対応点と、夫々の対応点の位置における第２の誤差と、放射基底関数によって記述される補間場の関係を示す図である。なお、図４では、座標や誤差を１次元に簡略化して記載している。この例では、３つの対応点４０１（４０１ａ，ｂ，ｃ）が存在していて、夫々の点において第２の誤差４０２（４０２ａ，ｂ，ｃ）が得られている。このとき、空間中の各点における第２の誤差は、誤差４０２のＴＰＳによる補間場４０３として取得される。また、図５（ａ）は、３次元空間における対応点４０１（４０１ａ，ｂ，ｃ，ｄ）と、夫々の対応点の位置における第２の誤差４０２（４０２ａ，ｂ，ｃ，ｄ）と、その補間場４０３との関係を表す模式図である。

ここで、実際に変形画像の生成に用いる第１の位置合わせ処理の結果として生じる各対応点位置における誤差は、第２の誤差（ＴＲＥ）ではなく第１の誤差（ＦＲＥ）である。そこで、ステップＳ３１０において、誤差推定部１７０は、ステップＳ２６０で得た夫々の対応点の位置における第２の誤差と第１の誤差の差異に基づいて、ステップＳ３００で得た第２の誤差の補間場Ｅ２（ｘ）を補正するための補正場Ｄ（ｘ）を導出する。具体的には、夫々の対応点の位置における補正後の誤差が第１の誤差と一致するように第２の誤差を第１の誤差に補正する補正ベクトルを考え、当該補正ベクトルに基づいて、対応点以外の位置における誤差の補正値を与える補正モデルのパラメータを導出する。補正モデルとしては、コンパクトサポートを持つ放射基底関数を夫々の対応点位置に配置し、対応点の位置に近いほど当該対応点と同様な補正を行うような（距離に応じて補正値を伝搬させる）モデルを用いることができる。なお、本実施形態では、放射基底関数としてＷｅｎｄｌａｎｄ関数を用いるものとする。

図４（ｂ）は、対応点と、夫々の対応点の位置における第２の誤差及び第１の誤差と、該位置における補正ベクトルの関係を示す図である。また、図４（ｃ）は、補正ベクトルと、該補正ベクトルによって定義される補正場との関係を示す図である。図４（ｂ）は、３つの対応点４０１（４０１ａ，ｂ，ｃ）の夫々に関して第２の誤差４０２（４０２ａ，ｂ，ｃ）と第１の誤差４０４（４０４ａ，ｂ，ｃ）が得られた場合に、夫々の補正ベクトル４０５（４０５ａ，ｂ，ｃ）が定義される様子を示している。このとき、図４（ｃ）に示すように、夫々の対応点位置に設置された放射基底関数によって、夫々の補正ベクトルに基づく補正値を近傍領域に減衰させながら伝搬させるような補正場４０６（４０６ａ，ｂ，ｃ）が定義される。

次に、ステップＳ３２０において、誤差推定部１７０は、ステップＳ３００で得た第２の誤差の補間場Ｅ２（ｘ）と、ステップＳ３１０で得た補正場Ｄ（ｘ）とを用いて、変形画像上の所定のボクセル位置における推定誤差を算出する。具体的には、任意の座標ｘにおける推定誤差を算出する関数をＥ（ｘ）＝Ｅ２（ｘ）−Ｄ（ｘ）と定義し、Ｅ（ｘ）の値（ベクトル）を座標ｘにおける推定誤差とする。以上の処理を、変形画像上の所定のボクセルに対して実行することで、推定誤差画像を生成する。

以上のように、本実施形態の推定処理では、各対応情報の算出された位置合わせ誤差（第２の誤差）を補間することにより、補間された位置合わせ誤差を示す補間場（Ｅ２（ｘ））が生成される。そして、各対応情報の位置において補正量が極大値となり、対応情報の位置からの距離に応じて補正量が減少する補正場（Ｄ（ｘ））により補間場が修正され（Ｅ（ｘ）＝Ｅ２（ｘ）−Ｄ（ｘ））る。そして、修正された補間場（Ｅ（ｘ））から任意の位置における推定された位置合わせ誤差が取得される。あるいは、位置合わせ誤差が各対応情報の位置で極小値となるように補間場を修正し、修正された補間場から任意の位置における推定された位置合わせ誤差を取得するようにしてもよい。

また、上記実施形態では、補間場の各対応情報の位置の位置合わせ誤差が、複数の対応情報を用いて第１の画像と第２の画像を位置合わせした場合の各対応情報の位置に残る位置合わせ誤差と一致するように補間場を修正するようにした（図４（ｂ）、（ｃ））。しかしながら、補間場の修正はこれに限られるものではない。第２の位置合わせにより各対応情報の位置に生じる位置合わせ誤差に比べて第１の位置合わせにより各対応情報の位置に生じる位置合わせ誤差が非常に小さくなる場合、各対応情報の位置の位置合わせ誤差がゼロになるように補間場を修正するようにしてもよい。これは、図４（ｂ）において第１の誤差４０４ａ〜４０４ｃをゼロとすることと同等である。

なお、推定誤差を計算する領域は、変形画像の全体（全てのボクセル）としてもよいし、注目する臓器領域内や病変等の注目領域内としてもよい。注目する臓器領域や病変等の注目領域は、例えば、ステップＳ２００で取得された画像の臓器領域や注目領域の情報をデータサーバ２０から取得することにより得ることができる。あるいは、画像の閾値処理や既存の領域抽出の方法を用いて、画像から臓器領域や病変領域を抽出するようにしてもよい。あるいは、ユーザが指定した画像上の領域を注目領域として取得するようにしてもよい。そして注目する臓器領域内や病変等の注目領域内の点を注目点として設定し、夫々の注目点について推定誤差を取得する。このようにすることで、後段の処理に不要な計算を省略できる。また、注目点は、領域内の全ボクセルとしてもよいし、一定間隔毎（例えば、５ボクセルおきなど）としてもよい。このように計算対象のボクセルを減らすことで、推定誤差の計算時間を短縮することができる。

図４（ｄ）は、図４（ａ）で示した第２の誤差の補間場４０３と、図４（ｃ）で示した第１の誤差への補正場４０６と、推定誤差４０７の関係を示す図である。図４（ｄ）に示すように、推定誤差４０７は、対応点の位置では当該位置における第１の誤差と等しい値となり、また、対応点から離れるに応じて第２の誤差の補間値に近い値となる。また、図５（ｂ）は、３次元空間における対応点４０１（４０１ａ，ｂ，ｃ，ｄ）と、夫々の対応点の位置を中心に分布する誤差の補正場４０６（４０６ａ，ｂ，ｃ，ｄ）と、それによって補正された推定誤差４０７の分布との関係を表す模式図である。以上のようにして、誤差推定部１７０によるステップＳ２８０の処理が実行される。

次に、ステップＳ２９０における、表示制御部１８０による推定誤差マップの表示について説明する。表示制御部１８０は、推定処理（ステップＳ２８０）によって推定された位置合わせ誤差を変形画像の断面画像に重畳して表示部としてのディスプレイ３０に表示する。以下、この重畳表示について具体例を示す。

図６（ａ）に、ディスプレイ３０に表示される変形画像の断面画像６１０の例を示す。また、図６（ｂ）に、ディスプレイ３０に表示される推定誤差マップの表示の一例を示す。この例では、変形画像の断面画像６１０に重畳する形で、推定誤差マップ６２０が表示されている。ここで表示される推定誤差マップ６２０は、ステップＳ２８０で取得した推定誤差画像（ボリュームデータ）のいずれかを、変形画像の断面画像６１０と対応する断面で切った断面画像である。画像の輝度値は、推定誤差マップのボクセル値を変換した値とする。例えば、所定の推定誤差（例えば、推定誤差＝１０ｍｍ）を輝度値＝２５５、推定誤差＝０ｍｍを輝度値＝０に変換することでグレースケールの推定誤差マップを作成して、これに疑似カラーを割りあて疑似カラーマップとして表示する。なお、これらの値の範囲はユーザが設定できることが望ましい。

なお、推定誤差マップは、図６（ｂ）のように変形画像の断面画像に重畳して表示されてもよいし、変形画像の断面画像と並べて表示されてもよい。また、ユーザが重畳表示のオン／オフを操作できるようにしてもよい。例えば、ディスプレイ３０に表示されるＧＵＩ（不図示）に重畳表示ボタンを提供し、これをユーザがオン／オフすることで重畳表示のオン／オフが制御されるようにする。この場合、重畳表示ボタンがオフの場合は、推定誤差マップを非表示にして断面画像６１０のみを表示し、オンの場合は、推定誤差マップ６２０を断面画像６１０に重畳表示する。

また、ステップＳ２８０で各軸方向の推定誤差を求めている場合には、夫々の推定誤差マップを並べて表示するようにしてもよいし、何れの軸方向の推定誤差マップを表示するかを選択できるようにしてもよい。或いは、推定誤差の大きさ（ノルム）の表示を選択できるようにしてもよい。

さらに、断面画像６１０上に（あるいは、その近傍に）対応点が存在する場合には、夫々の対応点の近傍位置に第１の誤差（ＦＲＥ）および／または第２の誤差（ＴＲＥ）の値を表示するようにしてもよい。その際に、カーソル等でユーザが指示した該断面画像上の点を注目ボクセルとして、その近傍の対応点のみの誤差を表示してもよい。また、誤差を表示する位置は対応点の近傍位置に限らず、表示画面の横や上などの所定の位置に、例えばリスト形式で表示してもよい。

なお、ステップＳ２８０で表示制御部１８０が表示する推定誤差マップは、推定誤差画像に所定の加工を施したものであってもよい。例えば、推定誤差が閾値以上のボクセルの上にのみ赤色等を半透明で重畳表示することで、推定誤差の不十分な部位をより明示的に確認できるようにしてもよい。また、どのような加工が施された推定誤差マップを重畳表示するかを、ユーザが選択できるようにしてもよい。

以上説明したように、第１実施形態によれば、対応点が設定されていない、変形画像上の任意の位置における、変形位置合わせの推定誤差を取得することができる。変形画像の断面画像上に推定誤差マップが重畳表示されるため、ユーザは、表示された断面画像の各点が実際にはその位置からどの程度ずれている可能性があるのかを、容易に把握できる。

（変形例１）
第１実施形態では、第２の誤差の補間場を放射基底関数で表現する例について説明したが、補間場の表現はこれに限られるものではない。例えば、夫々の対応点の位置における第２の誤差の線形補間によって補間場を求めてもよいし、他の任意の関数による近似で補間を行ってもよい。また、夫々の対応点の位置における第２の誤差を夫々の位置における変位とみなして、その変位に基づいて全体の変位場を求める（変形を推定する）何れかの方法を用いて変位場を求め、得られた変位場を補間場としてもよい。例えば、公知の変形表現手法であるＦＦＤ（Free-Form Deformation）等を用いることができる。より具体的には、等間隔なグリッド点にＦＦＤの制御点を設定し、夫々の対応点の位置における第２の誤差に基づいてＦＦＤのパラメータ（各制御点が持つ制御量）を算出することによって、補間場を表現することができる。

（変形例２）
第１実施形態では、推定誤差の分布を可視化した推定誤差マップを変形画像の断面画像上に重畳表示する例について説明したが、これに限らず、マップ以外の形態で推定誤差の表示を行ってもよい。

例えば、ステップＳ２９０において、表示制御部１８０は、第２の画像の断面画像と、ステップＳ２３０で生成した変形画像の対応断面画像とを、ユーザの操作に応じて、ディスプレイ３０に表示する。そして、表示制御部１８０は、変形画像の該断面画像上でユーザが指示した点の推定誤差を、ステップＳ２８０で取得した推定誤差画像から取得して、ディスプレイ３０に表示する制御を行う。ここで、ステップＳ２８０で推定誤差画像を取得する代わりに、ユーザが指定した点の推定誤差のみを取得してもよい。この場合には、ステップＳ２８０における推定誤差画像の生成処理は不要となる（変形例３で詳述する）。

図７に、変形例２における推定誤差の表示の例を示す。表示制御部１８０は、ディスプレイ３０に表示した変形画像の断面画像６１０上の座標の、マウス等により指定された位置を取得すると、該座標における推定誤差を推定誤差画像から取得する。そして、例えば図７（ａ）に示すように、カーソル７１０の示す座標（表示上の指定された位置）に対して、該座標における推定誤差の大きさを数値で表す文字情報７２０を、カーソル７１０の近傍に重畳表示する。

あるいは、推定誤差の大きさに対応した大きさの図形を指定された位置に表示するようにしてもよい。そのような図形の表示の一例を図７（ｂ）に示す。図７（ｂ）では、カーソル７３０の示す座標に対して、該座標における推定誤差の断面画像上での分布（誤差推定範囲）７４０を示す楕円が重畳表示されている。ここで、表示制御部１８０は、カーソル７３０の示す座標におけるＸ軸（断面画像の横軸）方向の推定誤差とＹ軸（断面画像の縦軸）方向の推定誤差の大きさを長径または短径として用いて、当該座標を中心とした楕円を描画する。これにより、表示された断面画像６１０上でカーソルを移動させるだけで、ユーザは、当該座標における誤差の推定値を知ることができる。特に、図７（ｂ）のような楕円を表示することにより、ユーザはどの方向に誤差が大きくなっているかを把握でき、変形画像の精度を向上させるために対応点をどこに追加すればよいかを容易に把握できる。

このように、本変形例では、カーソルの示す座標に対して、推定誤差を表す文字情報や誤差推定範囲を変形画像の断面画像上に表示している。ただし、この際の変形情報の取得方法や推定誤差の取得方法は、第１実施形態で説明した方法に限定されず、種々の方法を用いてもよい。変形画像と推定誤差とが対応付けられて表示されることにより、表示された断面画像の各点がその位置からどの程度ずれているのかをユーザが把握できればよい。以上のような変形例２によれば、変形の推定誤差または信頼度をマップ以外の表示方法で提示することが可能となり、誤差推定マップ等の重畳によって断面画像が見にくくなるという事態を回避できる。

（変形例３）
上述のように、ステップＳ２８０の処理において、ボリューム画像としての推定誤差画像の生成は必ずしも必須ではない。その代わりに、ステップＳ２９０の処理において、表示することが決定された変形画像の断面画像上の各ボクセルについてのみ、ステップＳ２８０と同様な推定誤差の取得を行い、表示する推定誤差の断面画像を直接生成してもよい。また、変形例２の表示を行う場合には、ステップＳ２９０の処理においてカーソル等で注目座標が指定された後に、ステップＳ２８０と同様な処理によって当該座標の推定誤差を適宜（オンデマンドで）推定するようにしてもよい。

（変形例４）
ステップＳ２１０の処理において、対応情報取得部１２０が実行する第１の対応情報の取得は、画像解析処理によって自動的に行ってもよい。例えば、夫々の画像から画像パターンの特徴的な点や線を検出し、画像パターンの類似性に基づいて自動で取得するようにしてもよい。また、画像解析処理によって自動取得した対応点を候補として、ユーザが手動で修正した点を最終的な対応点の位置としてもよい。なお、第１の対応情報の取得は、画像の付帯情報としてデータサーバ２０が画像間の対応情報を保持している場合には、該情報を読み込むことにより行ってもよい。

（変形例５）
第１実施形態や変形例２で示したステップＳ２９０における画像表示の処理は、必ずしも必須ではない。その代わりに、ステップＳ２３０で得た変形画像とステップＳ２８０で得た推定誤差画像をデータサーバ２０等に保存して、画像処理装置１０の処理を終了してもよい。この場合、ステップＳ２９０と同様な画像表示の処理を行う画像表示装置が、データサーバ２０から変形画像とその誤差推定値（推定誤差画像）を取得し、表示することにより、ユーザは任意の位置の推定誤差を観察することができる。また、変形画像や推定誤差画像の観察は、３次元断層画像を表示する通常の医用画像ビューアで行ってもよい。なお、画像表示を行わない場合であって、かつ、データサーバ２０から第１の対応情報が取得される場合には、第１の画像と第２の画像を取得する処理も省略できる。

（変形例６）
ステップＳ２８０の処理において、変形画像上の各点における位置合わせの推定誤差を計算する際に、変形に関する事前知識に基づく補正を行ってもよい。たとえば、第１の画像と第２の画像における、特定の部位から各対応情報の位置までの距離の違いに基づいて、算出処理により算出された第２の位置合わせ誤差を補正するようにしてもよい。そのような補正の一例として、被検体が乳房の場合には、事前知識として、乳頭から乳房内の対応点までの距離が体位等によらずあまり変化しないというモデルを用いることができる。この場合、例えば、夫々の対応点について得られる、夫々の画像における乳頭−対応点間の２つの距離の不一致度に基づいて、夫々の対応点の位置における第２の誤差の補正を行う。より具体的には、２つの距離の短いほうを分母として該２つの距離の比を算出し、その値を夫々の対応点の位置における第２の誤差に掛け合わせればよい。これによると、乳頭−対応点間距離が体位等によらずあまり変化しないというモデルから乖離するほど推定誤差が大きくなるように補正されるため、推定誤差の算出精度の向上が期待される。

事前知識としては、他にも、領域毎の変形推定の困難度を用いることもできる。すなわち、算出処理により算出された第２の位置合わせ誤差を、各対応情報が属する被検体の領域に設定されている変形推定の困難度に基づいて補正する。例えば、被検体が乳房の場合には、乳房内の各領域（例えば、乳癌取扱い規約において規定されるＡ〜Ｅ領域ごと）における変形推定の困難度を１以上の係数で表現する。そして、夫々の対応点が属する乳房内の領域に応じて、当該係数を夫々の対応点の位置における第２の誤差に掛け合わせればよい。例えば、体位等の変化による被検体の変形が大きいため変形推定が困難なＣ領域やＤ領域の困難度を１．５とし、体位等の変化による被検体の変形がより小さいＡ領域やＢ領域の困難度を１とすればよい。これによると、変形推定が困難な領域の方が、変形推定が容易な領域よりも推定誤差が大きくなるように補正されるため、推定誤差の算出精度の向上が期待される。

あるいは、病変の性質（硬さ等）を事前知識として用いてもよい。この場合には、周囲組織の硬さに対する病変の硬さの比と、夫々の対応点からの病変までの距離の逆数とを、夫々の対応点の位置における第２の誤差に掛け合わせればよい。これによると、病変が周囲組織よりも硬いほど変形推定の誤差が大きくなるような場合に、推定誤差の算出精度の向上が期待される。

（第２実施形態）
第１実施形態では、各対応点における第２の誤差から補間場と補正場を生成し、これらに基づいて任意の点における誤差を推定した。第２実施形態に係る画像処理装置では、第２の誤差に基づいて対応点からの距離と誤差の関係を求め、該関係に基づいて推定誤差を取得する。以下、第２実施形態に係る画像処理装置について、主として第１実施形態と異なる部分について説明する。なお、第２実施形態では、第１の誤差、第２の誤差、推定誤差には、ベクトルではなく、その大きさ（スカラー）が用いられる。

第２実施形態における画像診断システム１の構成、画像処理装置１０の各部の動作、処理手順は、第１実施形態と概ね同様である。ただし、ステップＳ２８０における推定誤差の算出方法（注目ボクセルにおける推定誤差の取得方法）が第１実施形態と異なる。第２実施形態の画像処理装置１０では、誤差推定部１７０が、
・各対応情報について、算出された第２の位置合わせ誤差に基づいて、対応情報の位置からの距離に応じた位置合わせ誤差の増加を表す誤差増加率を算出し、
・任意の位置の位置合わせ誤差を、対応情報の位置からの距離と誤差増加率に基づいて推定する。

以下、ステップＳ２８０における誤差推定部１７０の処理について、図８を参照して説明する。特に、本実施形態では、複数の対応情報の各々について、各対応情報の最近傍の対応情報までの距離と該最近傍の対応情報について算出された位置合わせ誤差とから各対応情報の誤差増加率を取得する。そして、得られた複数の誤差増加率を統計処理して代表誤差増加率を取得し、任意の位置に最も近い対応情報からの距離と代表誤差増加率から、任意の位置における位置合わせ誤差を推定する。以下、より具体的に説明する。

まず、誤差推定部１７０は、ステップＳ２６０で得た夫々の対応点の位置における第２の誤差に基づいて、対応点４０１（４０１ａ，ｂ，ｃ）からの距離に応じた誤差モデルを生成する。第２実施形態では、夫々の対応点における誤差増加率８０１（８０１ａ，ｂ）を求め、全対応点の誤差増加率を統計処理することで代表誤差増加率８０２を得る。そして、図８（ｂ）に示されるように、対応点からの距離に応じてこの比率（代表誤差増加率）で誤差が増加するという誤差モデルを生成する。例えば、夫々の対応点位置における誤差増加率８０１（８０１ａ，ｂ）の平均値（平均誤差増加率）を代表誤差増加率８０２とする。

ここで、代表増加率は、例えば図９のフローチャートで示される手順により取得される。まず、夫々の対応点における誤差増加率が、ステップＳ９０１〜Ｓ９０５により取得される。すなわち、誤差推定部１７０は、複数の対応点のうちの一つを選択し（ステップＳ９０１）、選択した対応点について、最も近い対応点（最近傍対応点）を同定する（ステップＳ９０２）。誤差推定部１７０は、ステップＳ９０２で同定した最近傍対応点の位置における第２の誤差（ＴＲＥ）を取得する（ステップＳ９０３）。そして、誤差推定部１７０は、ステップＳ９０３で取得した第２の誤差を、選択された対応点から同定された対応点までの距離で除することにより、選択された対応点の誤差増加率を得る。図８（ａ）では、対応点４０１ａの最近傍対応点は、対応点４０１ｂであり、対応点４０１ｂの最近傍対応点は、対応点４０１ｃである。そして、該誤差の値を該対応点から該最近傍対応点までの距離で除算した値を、該対応点における誤差増加率８０１ａ、８０１ｂとする。以上のステップＳ９０１〜Ｓ９０４の処理を全ての対応点について実行し、全ての対応点について誤差増加率が取得されると、処理はステップＳ９０６へ進む（ステップＳ９０５）。ステップＳ９０６において、誤差推定部１７０は、取得された全対応点の誤差増加率を統計処理（たとえば、平均値を算出）することにより代表誤差増加率を得る。

次に、上記で生成した誤差モデルに基づいて、注目ボクセルにおける誤差推定値を算出する。具体的には、注目ボクセルの位置から最も近い対応点を同定し、該対応点までの距離と誤差モデルとに基づいて、注目ボクセルにおける誤差推定値を算出する。すなわち、前記距離に代表誤差増加率を掛け合わせた値を注目ボクセルにおける推定誤差とする。たとえば、図８（ｂ）において、注目点の位置８０３の推定誤差は、位置８０３に最も近い対応点４０１ｂからの距離と上述の代表誤差増加率との積を対応点４０１ｂにおける第１の誤差４０４ｂに加算することで得られる。あるいは、位置８０３の近傍に存在する所定数の近傍対応点の夫々からの距離と上述の代表誤差増加率との積を、対応点４０１ｂにおける第１の誤差４０４ｂに夫々加算し、それらの中で最少となる値を推定誤差としてもよい。

なお、上記の代表誤差増加率は、夫々の対応点の誤差増加率の平均値に限定されるものではなく、他の方法で求めてもよい。例えば、平均値の代わりに誤差増加率の中央値を用いてもよい。また、代表誤差増加率を用いる以外の他の方法で誤差モデルを求めてもよい。例えば、夫々の対応点について得られる最近傍対応点までの距離と該最近傍対応点における第２の誤差との関係（距離と誤差の関係）を、任意の距離の関数（例えば２次関数）で近似して、該関数を誤差モデルとしてもよい。この場合、注目ボクセルにおける誤差推定値は、該ボクセルから最近傍対応点までの距離を当該関数に代入することによって算出される。すなわち、代表誤差増加率を距離に掛け合わせる上記の実施形態は、該関数として距離の１次関数を用いる場合に相当する。

また、誤差モデルの生成方法は、夫々の対応点の最近傍対応点のみを利用する方法に限定されるものではなく、その他の近傍対応点を利用する方法であってもよい。例えば、誤差を推定すべき任意の位置の近傍に存在する所定数の近傍対応点の夫々について上記と同様に誤差と距離との関係を得て局所的な代表誤差増加率を取得するようにしてもよい。例えば、誤差を推定する任意の位置について、その任意の位置の近傍の所定数の対応情報を選択する（たとえば、任意の位置から最も近い順に所定数の対応情報を選択する）。そして、選択した所定数の対応情報の各々について、対応情報の最近傍の対応情報までの距離と該最近傍の第２の位置合わせ誤差とから各対応情報の誤差増加率を取得し、これら誤差増加率を統計処理して局所的な代表誤差増加率を取得する。そして、任意の位置に最も近い対応情報からの距離と局所的な代表誤差増加率から、当該任意の位置における位置合わせ誤差を推定する。以上のように、任意の位置における位置合わせ誤差が、任意の位置の近傍の対応情報からの距離と代表誤差増加率や、任意の位置の近傍の対応情報からの距離と局所的な代表誤差増加率等に基づいて推定される。

以上説明したように、第２実施形態によれば、対応点からの距離と誤差の関係に基づいて推定誤差を取得することが可能となる。これによると、簡単な処理で推定誤差を算出できるという効果がある。また、局所的な代表誤差増加率を用いれば、より、画像の内容に即した位置合わせ誤差の推定が可能になる。

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

１：画像診断システム １０：画像処理装置 ２０：データサーバ ３０：ディスプレイ １１０：画像取得部 １２０：対応情報取得部 １３０：対応情報選択部 １４０：位置合わせ部 １５０：誤差算出部 １６０：変形画像生成部 １７０：誤差推定部 １８０：表示制御部

Claims (16)

1. 第１の画像および第２の画像に対する複数の対応情報を取得する取得手段と、
   前記複数の対応情報を用いて前記第１の画像と前記第２の画像を位置合わせし、位置合わせされた前記第１の画像および前記第２の画像に対する前記複数の対応情報によって定義された複数の位置に生じる位置合わせ誤差を算出する算出手段と、
   前記複数の位置とは異なる更なる位置に生じる位置合わせ誤差を、前記位置に生じる位置合わせ誤差に基づいて推定する推定手段と、を備えることを特徴とする画像処理装置。
2. 前記推定手段は、
   前記複数の位置に生じる位置合わせ誤差を補間することにより、補間された位置合わせ誤差を示す補間場を生成し、
   前記補間場に基づいて前記更なる位置に生じる位置合わせ誤差を推定することを特徴とする請求項１に記載の画像処理装置。
3. 前記推定手段は、前記複数の位置からの距離に応じて補正量が変化する補正場により前記補間場を修正し、修正された前記補間場から前記更なる位置に生じる位置合わせ誤差を取得することを特徴とする請求項２に記載の画像処理装置。
4. 前記推定手段は、前記補間場の前記複数の位置における位置合わせ誤差が、前記複数の対応情報を用いて前記第１の画像と前記第２の画像を位置合わせした場合の前記複数の位置に生じる位置合わせ誤差と対応するように前記補間場を修正することを特徴とする請求項２に記載の画像処理装置。
5. 前記推定手段は、前記補間場の前記複数の位置における位置合わせ誤差がゼロになるように前記補間場を修正することを特徴とする請求項２に記載の画像処理装置。
6. 前記推定手段は、
   前記複数の位置に生じる位置合わせ誤差に基づいて、前記複数の位置からの距離に応じた位置合わせ誤差の増加を表す誤差増加率を算出し、
   前記更なる位置に生じる位置合わせ誤差を、前記複数の位置からの距離と前記誤差増加率に基づいて推定する、ことを特徴とする請求項１に記載の画像処理装置。
7. 前記複数の対応情報を用いて前記第１の画像を前記第２の画像に位置合わせすることにより変形画像を得る変形手段と、
   前記推定手段により推定された位置合わせ誤差を前記変形画像に重畳して表示部に表示する表示制御手段と、をさらに備えることを特徴とする請求項１乃至６のいずれか１項に記載の画像処理装置。
8. 前記変形画像における各ボクセルまたは各ピクセルの値を、前記推定手段により推定された当該ボクセルまたはピクセルの位置の位置合わせ誤差に基づいて決定することで推定誤差画像を生成する生成手段をさらに備え、
   前記表示制御手段は、前記推定誤差画像を前記変形画像に重畳表示することを特徴とする請求項７に記載の画像処理装置。
9. 前記複数の位置に生じる位置合わせ誤差は、２次元または３次元の誤差ベクトルであり、
   前記推定手段は、前記誤差ベクトルの軸方向成分ごとの誤差を用いて前記更なる位置における位置合わせ誤差を軸方向成分ごとに推定し、
   前記生成手段は、前記軸方向成分ごとに前記推定誤差画像を生成することを特徴とする請求項８に記載の画像処理装置。
10. 前記表示制御手段は、前記変形画像の表示上で指定された位置について推定された位置合わせ誤差の大きさを数値により表示することを特徴とする請求項７に記載の画像処理装置。
11. 前記表示制御手段は、前記変形画像の表示上で指定された位置について推定された位置合わせ誤差の大きさに対応した大きさの図形を、前記指定された位置に表示することを特徴とする請求項７に記載の画像処理装置。
12. 前記複数の位置に生じる位置合わせ誤差は、２次元または３次元の誤差ベクトルであり、
    前記推定手段は、前記誤差ベクトルの軸方向成分ごとの誤差を用いて前記更なる位置における位置合わせ誤差を軸方向成分ごとに推定し、
    前記表示制御手段は、前記変形画像の表示上で指定された位置について、前記軸方向成分ごとに推定された誤差の大きさを径とする楕円を前記指定された位置を中心として表示することを特徴とする請求項７に記載の画像処理装置。
13. 前記第１の画像と前記第２の画像における、特定の部位から前記複数の位置までの距離の違いに基づいて、前記複数の位置に生じる位置合わせ誤差を補正する第１の補正手段をさらに備えることを特徴とする請求項１乃至１２のいずれか１項に記載の画像処理装置。
14. 前記複数の位置に生じる位置合わせ誤差を、前記複数の対応情報が属する被検体の領域に設定されている変形推定の困難度に基づいて補正する第２の補正手段をさらに備えることを特徴とする請求項１乃至１３のいずれか１項に記載の画像処理装置。
15. 第１の画像および第２の画像に対する複数の対応情報を取得する取得工程と、
    前記複数の対応情報を用いて前記第１の画像と前記第２の画像を位置合わせし、位置合わせされた前記第１の画像および前記第２の画像に対する前記複数の対応情報によって定義された複数の位置に生じる位置合わせ誤差を算出する算出工程と、
    前記複数の位置とは異なる更なる位置に生じる位置合わせ誤差を、前記位置に生じる位置合わせ誤差に基づいて推定する推定工程と、を備えることを特徴とする画像処理方法。
16. 請求項１乃至１４のいずれか１項に記載された画像処理装置の各手段としてコンピュータを機能させるためのプログラム。