JP6821303B2 - 情報処理装置、情報処理方法およびプログラム - Google Patents

Description

本発明は、情報処理装置、情報処理方法およびプログラムに関する。

複数のモダリティ、異なる日時、異なる体位、異なる撮像モードで撮像された画像を用いた画像診断において、複数の画像を位置合わせして表示／解析することが提案されている。非特許文献１には、画像間における被検体の変形を推定する変形位置合わせ技術が開示されている。ここで、より高精度な位置合わせを実現するために、ユーザが目視で対応付けした対応点の情報を取得して、この対応点の情報を位置合わせ処理に用いることが行われている。したがって、対応付けの精度によって位置合わせの精度が異なってくる。

そこで、一方の画像（以下、参照画像）上で指摘した注目点を他方の画像（以下、対象画像）から探索して対応付けする作業をユーザが行う際に、この作業を支援するための技術が開示されている。特許文献１には、注目点の推定位置を対象画像上に重畳表示する技術が開示されている。また、特許文献２には、注目点の存在確率の分布を対象画像上に重畳表示する技術が開示されている。

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

Ｄ．Ｒｕｅｃｋｅｒｔ, Ｌ．Ｓｏｎｏｄａ, Ｃ．Ｈａｙｅｓ, Ｄ．Ｈｉｌｌ, Ｍ．Ｌｅａｃｈ，ａｎｄ Ｄ．Ｈａｗｋｅｓ，"Ｎｏｎｒｉｇｉｄ ｒｅｇｉｓｔｒａｔｉｏｎ ｕｓｉｎｇ ｆｒｅｅ−ｆｏｒｍ ｄｅｆｏｒｍａｔｉｏｎｓ： ａｐｐｌｉｃａｔｉｏｎ ｔｏ ｂｒｅａｓｔ ＭＲ ｉｍａｇｅｓ"，ＩＥＥＥ ｍｅｄ． ｉｍａｇ．，ｖｏｌ．１８（８），ｐｐ．７１２−７２１，１９９９． Ｍａｒｉｕｓ Ｓｔａｒｉｎｇ，Ｓｔｅｆａｎ Ｋｌｅｉｎ，ａｎｄ Ｊｏｓｉｅｎ Ｐ．Ｗ．Ｐｌｕｉｍ，"Ａ ｒｉｇｉｄｉｔｙ ｐｅｎａｌｔｙ ｔｅｒｍ ｆｏｒ ｎｏｎｒｉｇｉｄ ｒｅｇｉｓｔｒａｔｉｏｎ"，Ｍｅｄ．Ｐｈｙｓ．３４（１１），Ｎｏｖｅｍｂｅｒ ２００７．

しかしながら、従来の技術では、ユーザが入力した対応点の正しさ（整合性）が分からないという課題がある。

本発明は、上記の課題に鑑みてなされたものであり、ユーザが入力した対応点を評価する技術を提供することを目的としている。

上記の目的を達成するために、本発明の一態様による情報処理装置は以下の構成を備える。即ち、
第一の画像における複数の特徴点と該複数の特徴点に対応する第二の画像の複数の特徴点との対である対応点の候補群を取得する取得手段と、
前記対応点の候補群に基づいて複数の異なる正則化条件下で第一の変位場と第二の変位場を算出するとともに、前記第一の画像における各特徴点を前記第一の変位場と前記第二の変位場により変位させて、第一の特徴点群と第二の特徴点群を算出する算出手段と、
前記第一の画像における特徴点を前記第一の変位場によって変位させた第一の特徴点群における第一の位置と、前記第一の画像における前記特徴点を前記第二の変位場によって変位させた第二の特徴点群における第二の位置との間の距離に基づいて、前記対応点の候補群の内、前記特徴点に関する一つの対応点の候補の整合度を算出する整合度算出手段と
を備え、
前記整合度算出手段は、前記第一の特徴点群と前記第二の特徴点群の夫々の位置に基づいて、前記対応点の候補群に含まれる対応点の候補の整合度を夫々算出することを特徴とする。

本発明によれば、ユーザが入力した対応点を評価することができる。

本発明の一実施形態に係る情報処理装置の構成例を示す図である。 本発明の一実施形態に係る情報処理装置が実施する処理の手順を示すフローチャートである。 本発明の一実施形態に係る第１の画像、第２の画像および対応点候補群を説明するための図である。

以下、添付図面に従って本発明に係る情報処理装置の好適な実施形態について説明する。なお、本発明は、添付図面に例示された例に限定されるものではない。

＜第１の実施形態＞
本実施形態に係る情報処理装置は、複数の３次元画像の変形位置合わせを、画像間の対応点の情報に基づいて行う。その際に、画像上に指定された画像間の対応点の候補の集合について、その夫々が、正しい対応点であることの尤もらしさの度合いを算出する。そして、算出した度合いをユーザが視認できる形態で表示する。本実施形態では、正しい対応点であることの尤もらしさの度合いを対応点の整合度と称する。ここで、「対応点」とは、「第１の画像上の注目点」と、「第２の画像上において該注目点と対応する点（解剖学的に同一の点）」との２点の組のことである。また、対応点の候補としてユーザが指定した第１および第２の画像上の２点のことを対応点候補と称する。また、対応点候補の集合を対応点候補群と称する。

本実施形態では、同一被検体を互いに異なる変形条件で撮像した３次元画像を用いる場合を例として説明する。また複数の３次元画像は、具体的には２つの３次元画像である場合を例として説明する。また、画像間の対応点の候補は、例えばユーザが手動で設定した画像中の特徴的な箇所（特徴点）の位置座標として取得することを想定する。この場合、当該対応点の候補の多くは、正しい対応点、即ち、被検体中の略同一の位置を表す位置座標値の組であるが、候補の一部には、略同一ではなく誤った対応点（実際の対応点ではない２点の組。誤対応点）が含まれることが想定される。また被検体中の略同一の位置を表す位置座標値であっても、そこに混入する誤差は対応点の候補毎にばらつきがあることが想定される。この場合、誤った対応点の候補や大きな誤差が混入した対応点の候補は、他の対応点候補が表す大域的な変位との間に不整合が生じる。

本実施形態における情報処理装置は、対応点候補の夫々について、この不整合の度合いを算出し、それを表示することで、誤った対応点候補や、大きな誤差が混入した対応点候補の存在を、ユーザが知覚できるように提示する仕組みを提供する。これにより、当該対応点候補の修正をユーザに促すことができる。また、ユーザによる対応点候補の修正を取得する構成や、修正後の対応点候補について再び整合度を算出して表示する構成を備える。すなわち、ユーザが対応点候補群の整合性を視認しながら、対応点候補群を修正できる構成を備える。

以下、図１乃至図３を用いて本実施形態に係る情報処理装置の構成および処理を説明する。

＜１．情報処理装置の構成＞
図１は、本実施形態に係る情報処理装置１００の機能構成例を示すブロック図である。図１に示すように、本実施形態における情報処理装置１００は、データサーバ１７０、操作部１８０、及び表示部１９０と接続されており、全体として情報処理システムを構成している。

データサーバ１７０は、ネットワークを介して情報処理装置１００と接続されており、以下に説明する第１の画像及び第２の画像を保持している。データサーバ１７０が保持する第１の画像及び第２の画像は、それぞれ異なる条件（異なるモダリティ、撮像モード、日時、体位等）で被検体を予め撮像して得られた３次元断層画像である。３次元断層画像を生成するモダリティは、ＭＲＩ装置、Ｘ線ＣＴ装置、３次元超音波撮影装置、光音響トモグラフィ装置、ＰＥＴ／ＳＰＥＣＴ、ＯＣＴ装置などである。第１の画像及び第２の画像は、例えば、異なるモダリティや異なる撮像モードで同時期に撮像されたものであってもよく、経過観察のために同一患者を同一モダリティ、同一体位で異なる日時に撮像した画像であってもよい。第１の画像及び第２の画像は、データ取得部１０２を介して情報処理装置１００に入力される。

さらに、データサーバ１７０は、第１の画像及び第２の画像における対応点候補群の情報として、夫々の対応点の組の３次元座標値の情報を保持している。これらのデータも、データ取得部１０２を介して情報処理装置１００に入力される。対応点候補群の中で、互いに対応付けられている第１の画像上の点と第２の画像上の点には、共通するＩＤ番号が付与されているものとする。すなわち、対応点候補群の各点は、点のＩＤ番号と、夫々の画像座標系における３次元座標（３次元空間におけるＸ，Ｙ，Ｚの３つの位置成分）とによって表されているものとする。

操作部１８０は、ユーザによるマウスやキーボードの操作を受け付ける。この操作により、ユーザが指定する対応点候補の修正情報が対応点修正情報取得部１０４を介して情報処理装置１００に入力される。また、ユーザは、操作部１８０を介して、対応点の修正処理を終了するか否かを指示することができる。

表示部１９０は、情報処理装置１００が生成する表示画像や整合度情報等の各種の情報を表示する。また、表示部１９０には、ユーザからの指示を取得するためのＧＵＩも配置されている。

情報処理装置１００は、データ取得部１０２と、対応点修正情報取得部１０４と、正則化条件生成部１０５と、変位場算出部１０６と、整合度算出部１０８と、表示制御部１１０とを備えている。データ取得部１０２は、データサーバ１７０から、第１の画像及び第２の画像と、対応点候補群の情報とを取得する。対応点修正情報取得部１０４は、ユーザによる操作部１８０の操作に従って、第１の画像上と第２の画像上とにおける対応点候補群に関する修正情報を取得する。正則化条件生成部１０５は、後述する変位場算出部１０６が変位場を算出する際の算出条件である正則化条件を複数生成する。

変位場算出部１０６は、対応点候補群を用いた第１の画像と第２の画像との間の変形位置合わせ処理を、複数の正則化条件の夫々に基づいて行い、第１の画像を変形させて第２の画像と位置合わせするための変位場を生成する。整合度算出部１０８は、変位場算出部１０６から複数の変位場の情報を取得し、これらの情報に基づいて対応点候補群の夫々に関する整合度を算出する。表示制御部１１０は、第１の画像及び第２の画像と、整合度に関する情報とを、表示部１９０に表示させる制御を行う。

＜２．情報処理装置が実施する処理＞
続いて、図２は、本実施形態に係る情報処理装置１００が実施する処理の手順を示すフローチャートである。以下、本実施形態における処理の詳細について説明する。

（ステップＳ３００：データの取得）
ステップＳ３００において、データ取得部１０２は、データサーバ１７０から、第１の画像、第２の画像、および対応点候補群の情報を取得する。そして、対応点候補群の情報を変位場算出部１０６へ出力し、第１の画像及び第２の画像を表示制御部１１０へ出力する。

ここで、図３（ａ）及び図３（ｂ）は、本処理ステップで取得する第１の画像、第２の画像、および対応点候補群の情報を説明する図である。図３（ａ）において、第１の画像２００は、データサーバ１７０から取得した画像であり、本実施形態の被検体を撮像した画像である。なお、書面による説明の都合上、同図の第１の画像２００として２次元の画像を表記しているが、これは、データ取得部１０２が取得する第１の画像２００の一断面を表すものであり、第１の画像２００は３次元画像である。本実施形態における第１の画像２００は、被検体をＭＲＩ装置で撮像した３次元画像である場合を例として説明する。図３（ａ）の第２の画像２１０は、第１の画像２００と同じ被検体を撮像した画像である。本実施形態における第２の画像２１０は、被検体をＣＴ装置で撮像した３次元画像である場合を例として説明する。第１の画像２００および第２の画像２１０は同一の被検体の略同一部位を撮像した画像であり、例えば血管などの被検体の解剖学的な構造が、双方の画像に類似した形態で描出される。ただし、両画像の撮像時における被検体の体位などの撮像条件の違いにより、両画像間に描出される構造には相違が生じうる。

図３（ｂ）は、第１の画像２００および第２の画像２１０に、取得した対応点候補（１番から５番）を付記した画像である。書面による説明の都合上、全ての対応点が同一平面上に含まれているように図示しているが、実際の対応点は同一平面外にも分布している。対応点候補は、例えば第１の画像２００および第２の画像２１０をユーザが観察し、両画像間で被検体の同一部分が描出されていると考えた位置を、夫々の画像上においてユーザが指定し、それを記録したものをデータ取得部１０２が取得する。本実施形態では、図３（ｂ）に示すように、被検体の血管構造の分岐に相当する位置等にユーザが対応点候補を設定した場合を例として説明する。

図３（ｂ）において、対応点候補１、２、５は、被検体の真に正しく対応する位置が指定されており、対応点候補３，４は、ユーザの誤りにより、正しく対応していない位置が指定されている状態を表している。具体的には、第１の画像２００上の対応点候補３，４は、第２の画像２１０上の対応点候補３，４と順番が入れ替わっている。対応点候補の誤りは、上記の例に限らず、例えば、任意の対応点の位置に誤差を含むような場合であってもよい。

本実施形態では、第１の画像２００における対応点候補の３次元の位置座標値をp_i(1≦i≦N)と表記する。また、第２の画像２１０における対応点候補の３次元の位置座標値をq_i(1≦i≦N)と表記する。ここで、Nは対応点候補の総数であり、図３に例示した場合であれば、N=5である。対応点候補のp_i、q_iの添え字iには、ユーザが同一の位置であるとして指定した点同士に、同一の数字を付する。つまり、p₁と q₁、p₂とq₂、...、はユーザが対応する位置として第１の画像２００および第２の画像２１０に夫々設定した点の位置座標値である。

次に、ステップＳ３０１以降の処理について説明する。本実施形態において、ステップＳ３０１からステップＳ３０６の処理は、複数回繰り返して実行されることがある。この場合、二回目以降に実行される場合には、初回の実行時と処理の一部が異なる。まず、ステップＳ３０１からステップＳ３０６が初回に実行する処理について説明し、二回目以降の実行時の処理については後述する。

（ステップＳ３０１：第１の変位場を算出）
ステップＳ３０１において、正則化条件生成部１０５は、第１の正則化条件を生成する。そして、変位場算出部１０６は、現在の対応点候補群の情報に基づき、第１の正則化条件の下で第１の変位場を算出する。そして、算出した第１の変位場を整合度算出部１０８へ出力する。

第１の変位場は、第１の画像上の対応点候補群の位置を、対応する第２の画像上の対応点候補群の位置へと近づけるように変位させ、かつ空間的に滑らかに変位を補間した場である。変位場の算出は、第１および第２の画像上の対応点候補群の位置の一致度合いと、変位場の空間的な滑らかさの度合いとの両方を加味し、これらをバランスさせるようにして実行する。この時、上記の二つの度合いの重みとして第１の正則化条件が適用される。変位場の算出に際して、正則化を強く働かせる正則化条件が適用される場合には、変位場は空間的な滑らかさ度合いが相対的に重視され、対応点候補群の位置の一致度合いは相対的に軽視される。

上記に説明した正則化条件に基づく変位場の算出方法の具体的な実現手段として、例えばThin Plate Spline(TPS)法、またはそれに類する方法を用いることができる。

TPS法では、対応点候補群の一致度合いを評価するデータ項と、変位場の滑らかさを曲げエネルギーの大きさにより評価する正則化項との加重和による、式（１）に示すコスト関数を最小化する変位場Φを算出する。

ここで、式（１）におけるスカラー値λは、データ項と正則化項の重みを司るパラメータであり、本実施形態では正則化パラメータと称する。変位場Φは関数として作用し、例えば式（２）に示す計算をすることで、変位前の座標値xを変位後の座標値x'へと座標変換することができる。

なお、変位場Φを用いた座標変換は、具体的には式（３）に示す計算により実行される。

ここで、x_iは制御点位置、w_iは制御点毎に定める非線形成分のパラメータ、a_v(1≦v≦4)は空間全体に定める線形成分のパラメータであり、これらの情報により変位場Φが特徴づけられる。本実施形態においてx_i 、w_i、 a_vの夫々は３次元ベクトルである。また、関数U(a,b)は放射基底関数の一種であり、引数のベクトルaおよびベクトルbの距離に基づいて算出されるスカラー値の関数である。本実施形態では、関数U(a,b)として、放射基底関数の一種であるＴＰＳ関数を用いる。またL_vは引数のベクトル要素値等を返す関数であり、x={x,y,z}^Tとした場合には、L₀(x)=x、L₁(x)=y、L₂(x)=z、L₃(x)=1である。本実施形態では、式（３）における制御点位置x_iとして第１の画像における対応点候補の位置座標値p_iを用い、w_i、 a_vを式（１）のコスト関数に基づいて最適化する。すなわち、コスト関数E(Φ)を最小化する変位場Φを与えるようなw_iとa_vを算出する。

ここで、式（１）におけるE_distance(Φ)はコスト関数におけるデータ項であり、具体的には式（４）に示す計算により座標変換による対応点の位置の距離の誤差が計算される。

また式（１）におけるE_smoothness(Φ)は正則化項であり、式（５）に示す計算により、変位場Φの曲げエネルギーが計算される。

本処理ステップにおける具体的な処理手順は以下のとおりである。まず正則化条件生成部１０５が、第１の正則化条件として、第１の正則化パラメータλ₁を生成する。ここで、正則化パラメータλ₁は所定の任意の正のスカラー値であってよいが、後述するステップＳ３０２で正則化条件生成部１０５が第２の正則化条件として生成するλ₂とは異なる値をとる。このとき、変形後の対応点の位置の残差が０になるように（すなわち、Ｅ_distance（Φ）＝０となるように）、非常に小さな値をλ₁の値として設定できる。次に、変位場算出部１０６が、式（１）に示したコスト関数のλとしてλ₁を用いて、コスト関数E(Φ)を最小化する変位場Φを算出する。なお、非常に小さな値をλ₁の値として設定して、Ｅ_distance（Φ）＝０となる解を得る場合には、正則化を含まない一般的なTPSの解法（線形方程式による求解）により変位場Φを算出してもよい。本実施形態では、算出した変位場を第１の変位場Φ₁と表記する。

（ステップＳ３０２：第２の変位場を算出）
ステップＳ３０２において、正則化条件生成部１０５は、第２の正則化条件を生成する。そして、変位場算出部１０６は、現在の対応点候補群の情報に基づき、第２の正則化条件の下で第２の変位場を算出する。そして、算出した第２の変位場を整合度算出部１０８へ出力する。

第２の変位場の算出は、ステップＳ３０１で実行した第１の変位場の算出と同様の処理により実行される。ただし、第２の変位場の算出では、第１の正則化条件とは異なる第２の正則化条件が適用される。具体例として、正則化条件生成部１０５は、第２の正則化条件である正則化パラメータλ₂を、λ₂=λ₁+Δλとして生成する。ここでΔλは、第１の正則化条件と第２の正則化条件の差異を表す所定の微小な正のスカラー値である。なお、第２の正則化条件に基づいて第２の変位場を生成する変位場算出部１０６の具体的な処理はステップＳ３０１と同様であるため、詳しい説明は省略する。本処理ステップを実行した結果として第２の変位場Φ_２が算出される。

（ステップＳ３０３：整合度を算出）
ステップＳ３０３において、整合度算出部１０８は、ステップＳ３０１で算出した第１の変位場Φ₁とステップＳ３０２で算出した第２の変位場Φ₂とを用いて、対応点候補群の夫々に関する整合度を算出する。そして、算出した整合度を表示制御部１１０へ出力する。

ここで、整合度とは、対応付けの尤もらしさを指標化した値であり、例えば現実の被検体に生じえないような変位を表す対応点候補が存在する場合には、当該対応点候補に対する整合度は低い数値とする。

まず、整合度算出部１０８は、第１の画像における対応点候補群の位置座標p_iの夫々を、ステップＳ３０１で算出した変位場Φ₁によって変位させた位置情報を算出する。つまり、式（６）の計算処理を実行する。

この位置情報p'_iを第１の変位後対応点候補位置（第１の変位後の特徴点群）と称する。

次に、整合度算出部１０８は、第１の画像における対応点候補群の位置座標p_iの夫々を、ステップＳ３０２で算出した変位場Φ₂によって変位させた位置情報を算出する。つまり、式（７）の計算を実行する。

この位置情報p''_iを第２の変位後対応点候補位置と称する（第２の変位後の特徴点群）。

さらに、整合度算出部１０８は、第１の変位後対応点位置および第２の変位後対応点位置の夫々について、式（８）に示す計算処理を実行して、整合度c _iを算出する。

ここで整合度c_iは、対応点候補群の夫々に関する整合度であり、より正確には、全ての対応点候補群の中における当該対応点候補の整合度である。

（ステップＳ３０４：整合度の可視化）
ステップＳ３０４において、表示制御部１１０は、第１の画像、第２の画像、および対応点候補群の整合度の情報をユーザに提示する表示制御処理を実行する。

この処理は例えば、第１の画像および第２の画像の夫々をボリュームレンダリングした画像を生成し、その画像に対して対応点候補群の位置にマークを重畳する。この時、対応点候補群の夫々の位置に重畳するマークの色や大きさ、点滅表示の有無などの表示態様を、当該対応点候補の整合度に基づいて変えることで、ユーザが整合度の違いを視認できるようにできる。これ以外にも、対応点候補位置に整合度の数値を文字として提示してもよい。

これ以外にも、例えば対応点候補の整合度について任意の閾値を設定し、閾値よりも整合度が低い（絶対値が大きい）対応点候補だけを表示するようにして、修正の必要性が高い対応点候補をユーザが視認しやすくしてもよい。また、整合度の低い対応点候補の位置をレンダリング画像の中心に設定してもよいし、また、当該位置をズーム表示するなどしてもよい。また、第１の画像や第２の画像は表示せずに、３次元画像空間中における対応点候補の夫々の位置にマークを表示し、夫々の整合度の情報のみを表示してもよい。

表示画像の生成方法は上記の方法に限らず、例えば第１の画像および第２の画像を任意の平面で切断した断面画像を夫々生成し、その断面画像上の対応点候補位置に上述のような様態でマーク等を表示するようにしてもよい。この時、ユーザによる操作に基づいて、第１の画像および第２の画像を切断する平面の位置を切り替え、画像中の様々な位置の断面を表示できるようにしてもよい。

また、整合度の低い対応点候補の位置を切断した画像が表示されるように、平面の位置を対応点候補の整合度に基づいて変更するようにしてもよい。例えば、整合度の値で対応点候補をソートしておいて、ユーザの指示に応じて、最も整合度の低い対応点候補を含む第１の画像と第２の画像の夫々の所定の断面（例えば、アキシャル断面）の画像を並べて表示するようにできる。さらに、ユーザの指示に応じて表示を切り替え、整合度が２番目に低い対応点候補について同様の表示を行うようにできる。また同様に、整合度が低い順に対応点候補を順次切り替えて表示することができる。

また、整合度の表示は、上記のように対応点候補の位置に応じた分布を示すようなグラフィカルな表示ではなく、対応点候補のＩＤと整合度の値の対をテキスト情報としてコンソール等に表示するような形態であってもよい。このとき、情報は対応点候補のＩＤの順に表示してもよいし、整合度の値でソートして、整合度の低い対応点候補から順に表示してもよい。また、整合度が閾値以下の対応点候補だけを表示するようにしてもよい。

以上に説明した方法により、表示制御部１１０は、対応点候補の整合度をユーザが視認できる表示画像を生成し、表示制御部１１０の制御により表示画像を表示部１９０に表示する。

（ステップＳ３０５：ユーザによる対応点の修正）
ステップＳ３０５において、対応点修正情報取得部１０４は、ユーザによる操作部１８０の操作情報に基づいて、対応点候補群に関する修正情報を取得する処理を実行する。そして、取得した修正情報に基づいて情報処理装置１００に記録されている対応点候補群の情報を変更する処理を実行する。

例えば、対応点修正情報取得部１０４は、操作部１８０であるマウスやキーボードに対するユーザの操作を取得して、それに基づいて対応点候補群の位置情報の修正を行ってもよい。この時、ステップＳ３０４で表示部１９０に表示した表示画像上にＧＵＩを配置し、そのＧＵＩに対するユーザの操作に基づいて対応点候補群の位置情報の修正を行ってもよい。例えば、ユーザのマウス操作に基づいて表示画像上の対応点候補の位置を移動し、対応点候補の位置情報を移動後の位置に基づいて変更してもよい。

また、対応点修正情報取得部１０４が取得する情報は、位置情報の修正に限らず、例えば、対応点候補群の中から任意の対応点を削除したり、追加したりできるようにしてもよい。例えば、表示画像上のＧＵＩにより、対応点候補群の中の任意の対応点候補をユーザが選択できるようにし、さらにＧＵＩとして「削除」ボタンなどを配置する。そして、ユーザが「削除」ボタンを押した場合に、選択中の対応点候補を削除できるようにしてもよい。また、表示画像上のＧＵＩとして修正操作を終了させるボタンを配置し、当該ボタンをユーザが押した場合に本処理ステップの処理を終了できるようにしてもよい。

（ステップＳ３０６：終了判定）
ステップＳ３０６において、情報処理装置１００は、本実施形態におけるステップＳ３０１からステップＳ３０６までの繰り返し実行される処理のループを終了させるか否かの判定を行う。ユーザ操作に基づいて、終了又は継続の指示が入力されたことに応じて判定を行ってもよい。そして、ループを終了すると判定した場合には処理ステップをＳ３０７へと進め、そうでない場合には処理ステップをＳ３０１へと戻す。

以上の処理により、本実施形態におけるステップＳ３０１からステップＳ３０６の繰り返し処理の初回の処理が実行される。

ステップＳ３０１からステップＳ３０６の二回目以降の処理は、初回の処理に対して以下の点が異なる。ステップＳ３０１およびステップＳ３０２の処理は、前回のループ処理におけるステップＳ３０５の処理により修正が施された対応点候補群に基づいて処理が実行される。それ以外の処理は初回の実行と同様の処理を実行する。

以上により、ステップＳ３０１からステップＳ３０６の処理が実行される。

（ステップＳ３０７：修正データによる変形）
ステップＳ３０７において、変位場算出部１０６は、以上の処理によって更新された対応点候補群の情報を用いて、第１の画像と第２の画像の変形位置合わせ処理を行う。すなわち、第１の画像を変形させて第２の画像に位置合わせするための変位場を生成する。また、必要に応じて、上記の変位場に基づいて第１の画像を変形させた変形画像を生成する。そして、生成した変位場や生成した変形画像を、不図示の通信部を介してデータサーバ１７０に格納する。また、情報処理装置１００は、以上の処理により更新された対応点候補群の情報を、データサーバ１７０に不図示の通信部を介して格納する。

ここで、変位場の生成は、ステップＳ３０１と同様な処理によって行うことができる。このとき、正則化パラメータλの値は、ステップＳ３０１と同じ値（すなわちλ₁）を用いることができる。なお、この場合の変位場はステップＳ３０１で最後に求めた変位場と同じなので、変位場の算出は省略してもよい。また、正則化パラメータλとして別の値を用いてもよい。例えば、変形後の対応点の位置の残差が０になるように、非常に小さな値をλに設定できる。また、対応点の情報に基づく任意の（ステップＳ３０１とは異なる）位置合わせ方法を用いてもよい。本ステップの処理はユーザによるインタラクティブな作業が不要であるので、ＴＰＳよりも計算時間を要する、非線形最適化計算を要するＦＦＤなどのより高精度な位置合わせ手法を用いることができる。また、輝度値の一致度等の対応点以外の情報をコスト関数として併用した位置合わせを行ってもよい。

以上のようにして、本実施形態における情報処理装置の処理が実行される。本実施形態によれば、対応点候補群の夫々に関する整合度をユーザが容易に視認できる。また、整合度に基づいてユーザが対応点候補の修正を行うことができる。すなわち、整合性の高い（誤対応の少ない）対応点候補群をユーザが容易に設定できる。その結果、精度の高い位置合わせを行うことができる。

（変形例１−１：整合度算出方法のバリエーション）
本実施形態では、ステップＳ３０３における整合度c_iの算出方法として、第１および第２の変位場に基づく変位後対応点位置の間の距離に基づき、式（８）に示した計算を実行する場合を例として説明した。しかし、本発明の実施形態はこれに限られるものではなく、ステップＳ３０３において、整合度算出部１０８は、第１および第２の変位場に基づく他の方法によって、対応点候補の整合度を算出してもよい。例えば、整合度算出部１０８は、第１の変位場による変位量Φ₁(p_i)と第２の変位場による変位量Φ₂(p_i)との差異に基づいて、c_i＝−||Φ₁(p_i) −Φ₂(p_i)||²として整合度を算出してもよい。また、整合度算出部１０８は、式（９）に示すように第１および第２の変位場に基づく変位後対応点位置の間の距離に関する所定の関数に基づいて整合度を算出してもよい。

ここでf_distributionは任意の分布関数であり、例えばガウス分布などを用いることができる。またf_distributionは対象物体の変形の複雑さ等の事前知識を考慮して、関数の特性を適切に設定することが望ましい。この方法によれば、より高い精度で整合性を算出できる効果がある。

（変形例１−２：整合度の低い点を自動的に削除する）
本実施形態では、ステップＳ３０３の処理で算出した整合度を、ステップＳ３０４の処理でユーザに提示し、ステップＳ３０５の処理でユーザによる修正を取得する場合を例として説明した。しかし、本発明の実施形態はこれに限定されない。例えば、情報処理装置１００は、ステップＳ３０３の処理で算出した整合度に基づいて対応点候補群の夫々に関する有効性を判定し、その判定に基づいて対応点候補群を修正する不図示の対応点修正部を備えるようにしてもよい。より具体的には、対応点修正部は、ステップＳ３０３の次の処理として、ステップＳ３０３で算出した対応点候補群の夫々に関する整合度が所定の閾値以下であるか否かを判定し、閾値以下の整合度となる対応点候補を対応点候補群から削除する処理を行う。これによれば、ユーザによる修正処理を要することなく、整合度の低い対応点候補を削除する修正を施すことができる。なお、表示制御部１１０が「削除候補」を提示して、ユーザが承認した「削除候補」だけを削除する構成としてもよい。

整合度の値により整合性の低い対応点候補を削除する方法は上記の方法に限らず、例えば、整合度が所定の閾値より小さく、かつ対応点候補群の中で最も整合性の低い対応点候補を削除するようにしてもよい。さらには、上記の削除の処理と、ステップＳ３０１からステップＳ３０３で説明した整合度の算出処理とを繰り返し実行するようにしてもよい。この方法によれば、対応点候補群の夫々の対応点候補の整合度が所定の閾値以上となるように、対応点候補群を修正することができる。

（変形例１−３：差分のベクトルを保持して、不整合点の修正の方向を提示する）
本実施形態では、対応点候補の夫々に関する整合度を提示することで、それらの対応候補に関する修正の必要性をユーザに示す場合を例として説明した。しかし、本発明の実施形態はこれに限定されない。例えば、対応点候補の夫々に関する修正の必要性に加え、当該対応点候補に関する位置の修正方向（向き）のガイドや参考情報を提示できるようにしてもよい。より具体的には、ステップＳ３０３において、整合度算出部１０８は、整合度c_iの算出処理に加えて、修正方向の候補ベクトルv_iを式（１０）に示す計算処理で算出してもよい。

そして、ステップＳ３０４において、表示制御部１１０は、生成する表示画像上に、算出したベクトルv_iに基づいて向きを変えた矢印状のガイド情報を、当該対応点候補上や、その近傍に表示してもよい。この方法によれば、当該対応点候補の位置の修正に関して、整合度の大小に影響の大きい修正方向をユーザに提示できる。そのため、当該対応点候補の修正の手掛かりをユーザに提供できるため、修正情報を与えるユーザの負荷を軽減できる効果がある。

この時、算出したベクトルv_iに基づいて、表示画像のレンダリングの視点位置や姿勢を変えたり、断面画像を表示する際の断面位置や姿勢を変えたりしてもよい。表示画像としてボリュームレンダリング画像を生成する場合には、ボリュームレンダリングの視線方向をベクトルv_iに直行するようにしてもよい。また、表示画像として断面画像を生成する場合には、ベクトルv_iが断面画像の面内に含むように、断面画像の切り出しの位置や姿勢を調整するようにしてもよい。これらの方法によれば、対応点候補の修正方向をユーザが視認しやすくなり、ユーザによる修正の負荷を、より一層軽減できる効果がある。

＜第２の実施形態：FFDで体積保存正則化を施す場合の例＞
本発明の第２の実施形態について、第１の実施形態との差異を中心に説明する。本実施形態における情報処理システムの構成例は、図１に示した第１の実施形態における情報処理システムの構成例と同様である。また、本実施形態における情報処理装置１００の処理手順は、図２のフローチャートで示した第１の実施形態における情報処理装置１００の処理手順と同様である。ただし、ステップＳ３０１およびステップ３０２で正則化条件生成部１０５および変位場算出部１０６が実行する処理が、第１の実施形態とは異なる。以下、本実施形態における情報処理装置１００が実行するステップＳ３０１およびステップ３０２の処理について説明する。

（ステップＳ３０１：第１の変位場を算出）
ステップＳ３０１において、正則化条件生成部１０５は、第１の正則化条件を生成する。そして、変位場算出部１０６は、現在の対応点候補群の情報に基づき、第１の正則化条件の下で第１の変位場を算出する。そして、算出した第１の変位場を整合度算出部１０８へ出力する。

本実施形態における変位場は、Free Form Deformation(FFD)法により表す。FFD法はグリッド状に制御点を配置し、その各制御点における制御量に基づいて、空間中の変位量を算出する。この算出処理は公知の方法で実現可能であるため、詳細な説明は省略する。ここではFFD法において変位場の特性を定める制御量の算出方法について説明する。

本実施形態において、変位場算出部１０６は、式（１１）に示すコスト関数を最小化する変位場Φを算出する。

式（１１）においてE_distance(Φ)は、コスト関数におけるデータ項であり、第１の実施形態で説明した式（４）と同様であるため、ここでは説明を省略する。また、式（１１）においてE_bending(Φ)は、曲げエネルギーに関する正則化項であり、第１の実施形態において式（５）を用いて説明したE_smoothness(Φ)と同様であるため、ここでは説明を省略する。式（１１）においてE_compression(Φ)は、変位場Φがもたらす空間の変換における局所的な体積変化に関する正則化項である。また、式（１１）においてλ_bendingおよびλ_compressionは、それぞれ曲げエネルギーに関する正則化パラメータ、体積変化に関する正則化パラメータである。

本実施形態において、E_compression(Φ)について、例えば式（１２）のような計算を行う。

ここでJ(Φ)は変位場Φに関する空間微分である3×3のヤコビ行列を返す関数である。また、det()は引数の行列の行列式を求める関数である。式（１２）は変位場Φがもたらす空間変換において、局所的な体積変化の大きさを算出し、それを空間的に積分する演算を表している。

本処理ステップにおいて、正則化条件生成部１０５は、第１の正則化条件として、式（１１）に示したλ_bendingとλ_compressionを生成する。本実施形態ではこれらを、λ_bending,1、λ_{compression,1}と表記する。ここで、λ_bendingとλ_compressionの値としては、予め定めた正のスカラー値を用いることができる。このとき、λ_bendingとλ_compressionの少なくとも何れかの値は０であってもよい。すなわち、第１の正則化条件として、正則化を施さないという正則化条件を設定してもよい。

次に、変位場算出部１０６は、λ_bending,1とλ_{compression,1}とを用いて、式（１１）に示したコスト関数の最小化を評価基準として、それを最適化することで第１の変位場Φ₁を算出する。なお、変位場Φ₁の算出にために実行される最適化演算には、例えば最急降下法など任意の非線形最適化手法が適用できる。

（ステップＳ３０２：第２の変位場を算出）
ステップＳ３０２において、正則化条件生成部１０５は、第２の正則化条件を生成する。そして、変位場算出部１０６は、現在の対応点候補群の情報に基づき、第２の正則化条件の下で第２の変位場を算出する。そして、算出した第２の変位場を整合度算出部１０８へ出力する。

第２の変位場の算出は、ステップＳ３０１で実行した第１の変位場の算出と同様の処理により実行される。ただし、第２の変位場の算出では、第１の正則化条件とは異なる第２の正則化条件が適用される。すなわち、正則化条件生成部１０５が、第１の正則化条件であるλ_bending,1、λ_{compression,1}とは異なるλ_bending,2、λ_{compression,2}を生成する。ここで、λ_bending,2、λ_{compression,2}の具体的な生成方法は、第１の実施形態と同様である。

すなわち、第１の正則化条件である各正則化パラメータλ_bending,1、λ_{compression,1}の夫々に予め定めた正の微小なスカラー値を加算した値として、λ_bending,2、λ_{compression,2}を生成する。なお、第２の正則化条件の下で変位場を算出する変位場算出部１０６の具体的な処理は、ステップＳ３０１と同様であるため詳しい説明は省略する。本処理ステップを実行した結果として第２の変位場Φ_２が算出される。

以上により、本発明の第２の実施形態における情報処理装置の処理が実行される。本実施形態によれば、変位場の体積変化に関する正則化条件が異なる、複数の変位場が算出される。そして、それらの変位場に基づいて対応点候補群の夫々に関する整合度が算出される。これにより、局所的な体積変化が小さいことが想定される被検体を対象とした対応点候補群について、局所的な体積変化の観点から整合性を算出することができる。

（変形例２−１：ヤコビアンではなく剛体性を用いる）
本実施形態における、ステップＳ３０１およびステップＳ３０２の処理では、曲げエネルギーおよび局所的な体積変化に関する正則化条件の下で変位場を算出する場合を例として説明した。しかし、本発明の実施形態はこれに限定されない。例えば、正則化条件として、剛体性からの逸脱度を用いてもよい。例えば、非特許文献２に記載されているように、変位場の線形性、正規直交性、体積保存性の夫々に関する正則化を含むコスト関数を用いて変位場を算出するようにしてもよい。この方法によれば、局所的に剛体変換から大きく相違しない変形をする被検体を対象とした場合に、上記の被検体の変形特性を考慮した整合性を算出する仕組みを提供できる。

（変形例２−２：ヤコビアンと剛体性の両方を用いる。条件1, 条件2, 条件3 でFFD近似して２次元の整合度を算出（例：ヤコビアンと剛体性の両方のスコアが悪ければ不整合にする等））
本実施形態では、異なる２種類の正則化条件下で生成した変位場を用いて整合度を算出する場合を例として説明したが、本発明の実施形態はこれに限定されない。例えば、３種類以上の異なる正則化条件下で生成した変位場を用いて整合度を算出するようにしてもよい。例えば、以下の処理を行うことができる。

まず、ステップＳ３０１の処理において、正則化条件生成部１０５は、上記の実施形態と同様に、式（１１）のコスト関数に関する第１の正則化条件として、予め定めた正の値としてλ_bending,1、λ_{compression,1}を生成する。そして、変位場算出部１０６は、第１の正則化条件に基づき第１の変位場Φ₁を算出する。

次に、ステップＳ３０２の処理において、正則化条件生成部１０５は、第２の正則化条件を、λ_bending,2=λ_bending,1+Δλ_bending、λ_{compression,2}=λ_{compression,1}の計算により生成する。ここで、Δλ_bendingは、λ_bendingを微小変動させるための正の規定値である。

そして、変位場算出部１０６は、第２の正則化条件に基づいて第２の変位場Φ₂を算出する。さらに、正則化条件生成部１０５は、第３の正則化条件を、λ_bending,3=λ_bending,1、λ_{compression,3}=λ_{compression,1}+Δλ_compressionの計算により生成する。ここで、Δλ_compressionは、λ_compressionを微小変動させるための正の規定値である。そして、変位場算出部１０６は、第３の正則化条件に基づいて第３の変位場Φ₃を算出する。

そして、ステップＳ３０３の処理において、整合度算出部１０８は、上記の第１、第２、第３の変位場に基づいて整合度を算出する。例えば、第１の変位場と第２の変位場に基づいて本実施形態のステップＳ３０３と同様の処理を実行して第１の整合度を算出する。さらに、第１の変位場と第３の変位場に基づいて本実施形態のステップＳ３０３と同様の処理を実行して第２の整合度を算出する。

そして、第１の整合度と第２の整合度に基づき、これらを統合した整合度を算出する。例えば、対応点候補の夫々について、第１の整合度と第２の整合度の最小値を統合した整合度としてもよい。この方法によれば、曲げエネルギーまたは体積保存性の両方の観点に基づき、いずれかの観点で不整合となる対応点候補をユーザが視認することができる。上記の例では３種類の正則化条件に基づく場合を例として説明したが、本発明の実施形態はこれに限らず、４種類以上の正則化条件を用いてもよい。

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

１００：情報処理装置、１０２：データ取得部、１０４：対応点修正情報取得部、１０５：正則化条件生成部、１０６：変位場算出部、１０８：整合度算出部、１１０：表示制御部、１７０：データサーバ、１８０：操作部、１９０：表示部、２００，２１０：画像

Claims (11)

1. 第一の画像における複数の特徴点と該複数の特徴点に対応する第二の画像の複数の特徴点との対である対応点の候補群を取得する取得手段と、
   前記対応点の候補群に基づいて複数の異なる正則化条件下で第一の変位場と第二の変位場を算出するとともに、前記第一の画像における各特徴点を前記第一の変位場と前記第二の変位場により変位させて、第一の特徴点群と第二の特徴点群を算出する算出手段と、
   前記第一の画像における特徴点を前記第一の変位場によって変位させた第一の特徴点群における第一の位置と、前記第一の画像における前記特徴点を前記第二の変位場によって変位させた第二の特徴点群における第二の位置との間の距離に基づいて、前記対応点の候補群の内、前記特徴点に関する一つの対応点の候補の整合度を算出する整合度算出手段と
   を備え、
   前記整合度算出手段は、前記第一の特徴点群と前記第二の特徴点群の夫々の位置に基づいて、前記対応点の候補群に含まれる対応点の候補の整合度を夫々算出することを特徴とする情報処理装置。
2. 前記整合度算出手段は、前記第一の特徴点群と前記第二の特徴点群の夫々の位置の間の距離に基づいて前記候補群に含まれる対応点の候補の整合度を夫々算出することを特徴とする請求項１に記載の情報処理装置。
3. 前記複数の異なる正則化条件は、変位場の曲げエネルギーに関する正則化条件を含むことを特徴とする請求項１又は２に記載の情報処理装置。
4. 前記複数の異なる正則化条件は、変位場の体積保存性に関する正則化条件を含むことを特徴とする請求項１又は２に記載の情報処理装置。
5. 前記複数の異なる正則化条件のうちの１つは、正則化を施さない正則化条件であることを特徴とする請求項１又は２に記載の情報処理装置。
6. 前記整合度の情報を表示部に表示させる表示制御手段をさらに備えることを特徴とする請求項１乃至５の何れか１項に記載の情報処理装置。
7. 前記算出手段は、前記候補群に含まれる対応点の修正方向を示す情報をさらに算出し、
   前記表示制御手段は、前記修正方向を示す情報を前記表示部に表示させることを特徴とする請求項６に記載の情報処理装置。
8. ユーザ操作に基づいて前記候補群に関する修正情報を取得するとともに、前記修正情報に基づいて前記候補群の情報を修正する処理を実行する修正情報取得手段をさらに備えることを特徴とする請求項１乃至７の何れか１項に記載の情報処理装置。
9. 前記整合度が閾値以下の対応点を前記候補群から削除する削除手段をさらに備えることを特徴とする請求項１乃至８の何れか１項に記載の情報処理装置。
10. 第一の画像における複数の特徴点と該複数の特徴点に対応する第二の画像の複数の特徴点との対である対応点の候補群を取得する取得工程と、
    前記対応点の候補群に基づいて複数の異なる正則化条件下で第一の変位場と第二の変位場を算出するとともに、前記第一の画像における各特徴点を前記第一の変位場と前記第二の変位場により変位させて、第一の特徴点群と第二の特徴点群を算出する算出工程と、
    前記第一の画像における特徴点を前記第一の変位場によって変位させた第一の特徴点群における第一の位置と、前記第一の画像における前記特徴点を前記第二の変位場によって変位させた第二の特徴点群における第二の位置との間の距離に基づいて、前記対応点の候補群の内、前記特徴点に関する一つの対応点の候補の整合度を算出する整合度算出工程と
    を有し、
    前記整合度算出工程では、前記第一の特徴点群と前記第二の特徴点群の夫々の位置に基づいて、前記対応点の候補群に含まれる対応点の候補の整合度を夫々算出することを特徴とする情報処理方法。
11. 請求項１０に記載の情報処理方法の各工程をコンピュータに実行させるためのプログラム。

Images