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

処理装置、処理方法、およびプログラム Download PDF

Info

Publication number
JP6660428B2
JP6660428B2 JP2018145407A JP2018145407A JP6660428B2 JP 6660428 B2 JP6660428 B2 JP 6660428B2 JP 2018145407 A JP2018145407 A JP 2018145407A JP 2018145407 A JP2018145407 A JP 2018145407A JP 6660428 B2 JP6660428 B2 JP 6660428B2
Authority
JP
Japan
Prior art keywords
processing
anatomical feature
shape
cross
image data
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Active
Application number
JP2018145407A
Other languages
English (en)
Other versions
JP2018167077A (ja
Inventor
亮 石川
亮 石川
佐藤 清秀
清秀 佐藤
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Canon Inc
Original Assignee
Canon Inc
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Canon Inc filed Critical Canon Inc
Priority to JP2018145407A priority Critical patent/JP6660428B2/ja
Publication of JP2018167077A publication Critical patent/JP2018167077A/ja
Application granted granted Critical
Publication of JP6660428B2 publication Critical patent/JP6660428B2/ja
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Description

本発明は、核磁気共鳴映像装置(MRI)、X線コンピュータ断層撮影装置(X線CT)、超音波画像診断装置(US)など、種々の医用画像収集装置(モダリティ)で撮像した医用画像を処理する処理装置、処理方法、およびプログラムに関する。
医療の分野において、あるモダリティの画像上に注目部位があった場合に、その対応部位を他のモダリティの画像上で同定し、その対比によって診断を行う場合がある。モダリティ間で撮像体位が異なる場合には、撮像時の被検体の形状が異なるため、その同定が難しくなるという課題がある。そこで、双方の間における被検体の変形を推定すること(すなわち、変形を伴う画像間の位置合わせを行うこと)が試みられている。これにより、注目部位の位置情報に基づいて対応部位の位置を推定することや、一方の画像に変形を施して形状を他方と同一にした画像を生成することが可能となる。
非特許文献1には、変形を伴う物体の形状を正規化することで、当該物体の変形前後の形状の比較を容易にする技術が開示されている。具体的には、物体の表面形状の測地線距離行列を算出し、その行列を用いた多次元尺度構成法によって正規化を行う方法が開示されている。この方法によれば、物体表面の測地線距離が変わらないような変形に関して、その変形前後の形状間を直接比較可能な形態に正規化することができるようになる。これにより、変形を伴う物体の変形形状間の比較や、形状に基づく物体認識などが容易に行える。
A. Elad and R. Kimmel, "On bending invariant signatures for surfaces," IEEE Trans. PAMI, 25(10), 2003 Daniel Rueckert, Alejandro F. Frangi, and Julia A. Shnabel, "Automatic construction of 3-D statistical deformation models of the brain using Nonrigid registration," IEEE Transaction on medical imaging, Vol. 22, No. 8, 2003
非特許文献1に記載の方法を用いると、変形を伴う複雑な形状に関して、その変形の前後における形状間の正規化が行えるため、位置合わせを比較的に容易に行えることが期待できる。しかし、物体の形状が比較的に単調で、変形形状間で対応付けできるランドマーク等が少ない場合には、正規化後の物体の姿勢に不安定さが残るという課題があった。
本発明は、上記の課題に鑑みてなされたものであり、異なる形状間の変形を高精度に表現可能な変形モデルの構築を行う仕組みを提供することを目的とする。
上記の目的を達成するための、本発明の一態様によるに処理装置は以下の構成を備える。すなわち、対象物の三次元のボリュームデータである三次元画像データを取得する画像取得手段と、前記対象物における少なくとも第一の解剖学特徴と第二の解剖学特徴を取得する特徴取得手段と、前記三次元画像データから、前記第一の解剖学特徴における第一の面と前記第二の解剖学特徴における第二の面との距離が所定の条件を満たす位置の集合で構成される断面に対応する断面画像を生成する生成手段と、有する。
本発明によれば、多数の変形に関する事例データについて、変形前および変形後の形状に基づいて変形を正規化することで、異なる形状間の変形を高精度に表現可能な変形モデルの構築を行う仕組みを提供できる。
第1実施形態による処理システムの機能構成を示す図。 第1実施形態による処理システムの装置構成を示す図。 第1実施形態による処理装置の処理手順を示すフローチャート。 伏臥位MRI画像に描出される被検体の模式図。 第1実施形態によるステップS320の処理手順を示すフローチャート。 第1実施形態による正規化座標系を説明する図。 第1実施形態による画像変形処理を説明する図。 第3実施形態による処理システムの機能構成を示す図。 第3実施形態による処理装置の学習フェーズの処理手順を示すフローチャート。 第3実施形態による処理装置のステップS580の処理手順を示すフローチャート。 第3実施形態による処理装置の変形推定フェーズの処理手順を示すフローチャート。 第4実施形態による処理システムの機能構成を示す図。 第4実施形態による処理装置の処理手順を示すフローチャート。 第5実施形態による処理システムの機能構成を示す図。 第5実施形態による処理装置の学習フェーズの処理手順を示すフローチャート。 第5実施形態による処理装置の変形推定フェーズの処理手順を示すフローチャート。 第6実施形態による処理システムの機能構成を示す図。 第6実施形態による処理装置の学習フェーズの処理手順を示すフローチャート。 第6実施形態による処理装置の変形推定フェーズの処理手順を示すフローチャート。
以下、添付図面に従って本発明による処理装置及び方法の好ましい実施形態について詳説する。ただし、発明の範囲は図示例に限定されるものではない。
[第1実施形態]
(第1実施形態の概要)
本実施形態による処理装置は、注目部位としての被検体の乳房を異なる二種類の体位で夫々撮像した医用画像が取得された場合に、夫々の医用画像を基準形状に変換する正規化変換(正規化変換情報)を求めた上で、それを介して画像間の変形位置合わせを行う。この正規化変換は、体位の違いにより、異なる変形状態で撮像された被検体の乳房を、解剖学的に略一致した空間に座標変換する変換である。
ここで、伏臥位と仰臥位の間の乳房の変形に関する特性として、生体力学的に以下の事が知られている。第一に、冠状面において、乳頭を基準とした方位は、概ね保持されること。第二に、体表面において、乳頭を基準とした測地線距離は概ね保存されること。これらの特性により、乳頭を基準とした方位と測地線距離とを基準化した空間を考え、変形状態の異なる乳房の夫々を、その空間に座標変換する。これにより、伏臥位と仰臥位との間で生じる変形に伴う位置の変動を概ね吸収し、解剖学的に共通する空間へと変換できる。この変換を介して画像間の変形位置合わせを行うことにより、元画像を直接的に変形位置合わせするよりも高い精度で位置合わせすることができるようになる。
本実施形態では、取得された医用画像の夫々から被検体の乳房の輪郭を抽出し、輪郭上の基準点に基づいて、その形状を基準形状である矩形形状に座標変換する正規化変換を算出する。なお、以下の説明では、被検体として仰臥位および伏臥位における人体の乳房を対象とする場合について説明するが、被検体の体位や部位は特に限定されるものではない。また本実施形態では、医用画像の一例として三次元のMRI画像を用いる場合を例として説明するが、MRI画像に限定されず、他の三次元画像であっても良い。例えば、X線CT画像や3D超音波画像、PET画像などに適用できる。
(機能構成)
図1は、本実施形態による処理システムの構成を示す図である。同図に示すように、本実施形態における処理装置100は、画像取得部1000、解剖学特徴抽出部1020、正規化部1040、画像変形部1060、観察画像生成部1080によって構成される。そして、処理装置100は、データサーバ120、モニタ160に接続される。MRI画像撮像装置110は、人体である被検体の内部の三次元領域に関する情報を核磁気共鳴法により取得した画像、すなわちMRI画像を取得する装置である。MRI画像撮像装置110はデータサーバ120と接続され、取得したMRI画像をデータサーバ120へ送信する。データサーバ120は、MRI画像撮像装置110が撮像したMRI画像を保持する装置であり、処理装置100からの命令により保持したMRI画像を処理装置100へ転送する。
次に、処理装置100を構成する各要素について説明する。画像取得部1000は、MRI画像撮像装置110によって撮像された被検体(対象物体)のMRI画像を、データサーバ120を介して処理装置100に取り込む。解剖学特徴抽出部1020は、画像取得部1000が取り込んだMRI画像に画像処理を施し、被検体の解剖学特徴を抽出する。正規化部1040は、解剖学特徴抽出部1020が抽出した被検体の解剖学特徴に基づいて、被検体の形状を基準形状に変換する(正規化する)ための変換を算出する。正規化に関する詳細は後述する。画像変形部1060は、正規化部1040が算出した変換に基づいて伏臥位-仰臥位間の位置合わせを行い、伏臥位MRI画像を変形させて仰臥位MRI画像に合わせた変形画像を生成する。観察画像生成部1080は、画像取得部1000が取り込んだMRI画像と、画像変形部1060が生成した変形画像の夫々から、ユーザに提示する観察画像を生成する。そして、その観察画像をモニタ160へと出力する。モニタ160は観察画像生成部1080が生成した観察画像を表示する。
(装置構成)
図2は、本実施形態による処理システムの装置構成を示す図である。本実施形態の処理システムは処理装置100、MRI画像撮像装置110、データサーバ120、モニタ160、マウス170、キーボード180により構成される。処理装置100は、例えばパーソナルコンピュータ(PC)などで実現することができる。処理装置100は、中央演算処理装置(CPU)211、主メモリ212、磁気ディスク213、表示メモリ214を有する。
CPU211は、主として処理装置100の各構成要素の動作を制御する。主メモリ212は、CPU211が実行する制御プログラムを格納したり、CPU211によるプログラム実行時の作業領域を提供したりする。磁気ディスク213は、オペレーティングシステム(OS)、周辺機器のデバイスドライブ、後述する処理等を行うためのプログラムを含む各種アプリケーションソフト等を格納する。表示メモリ214は、モニタ160のための表示用データを一時記憶する。モニタ160は、例えばCRTモニタや液晶モニタ等であり、表示メモリ214からのデータに基づいて画像を表示する。マウス170及びキーボード180は、ユーザによるポインティング入力及び文字やコマンド等の入力をそれぞれ行う。上記各構成要素は共通バス218により互いに通信可能に接続されている。
(処理フロー)
次に、処理装置100が行う処理に関して、図3のフローチャートを用いて詳しく説明する。図3は本実施形態において処理装置100が実行する処理のフローチャートである。本実施形態では、主メモリ212に格納されている各部の機能を実現するプログラムをCPU211が実行することにより実現される。また以下に説明する処理装置100が行う各処理の結果は、主メモリ212に格納することにより記録される。まず、図3に示す各処理ステップについて、その手順を追って詳しく説明する。
(ステップS300)伏臥位MRI画像を取得
ステップS300において、画像取得部1000は、MRI画像撮像装置110が被検体の乳房を伏臥位の体位で撮像したMRI画像(伏臥位MRI画像)を、データサーバ120を介して処理装置100へ取り込む処理を実行する。ここで伏臥位MRI画像は、三次元のボリュームデータであり、被検体の足側から頭側に向かう方向をZ軸、腹側から背側に向かう方向をY軸、被検体の左方向をX軸とする三次元の直交座標系を持つ(そのような座標変換が予め施されている)ものとする。本実施形態では、この座標系を伏臥位MRI画像座標系と称する。また、伏臥位MRI画像の輝度値を、伏臥位MRI画像座標系における三次元位置xを引数としたスカラ関数Ip(x)と表記する。
(ステップS310)伏臥位MRI画像から解剖学特徴を抽出
ステップS310において、解剖学特徴抽出部1020は、ステップS300で取得した伏臥位MRI画像を処理することにより、被検体の伏臥位における解剖学特徴を抽出する処理を実行する。本実施形態において解剖学特徴とは、被検体の乳頭位置、体表面形状、大胸筋面形状、大胸筋面上の基準位置である。
図4は、伏臥位MRI画像上における解剖学特徴を説明する図である。実際の伏臥位MRI画像は三次元画像であるが、ここでは紙面による説明の都合上、三次元画像上のある2次元断面を示して説明する。伏臥位MRI画像400には、空気領域403、乳房領域402、内部領域405が存在する。体表面401は、空気領域403と乳房領域402との境界の位置の集合であり、三次元の曲面である。また、大胸筋面404は、乳房領域402と内部領域405との境界の位置の集合であり、三次元の曲面である。本処理ステップにおいて、解剖学特徴抽出部1020は、伏臥位MRI画像400を閾値処理やエッジ検出など周知の方法で画像処理し、体表面401を検出する。ただし、体表面401の検出は伏臥位MRI画像に描出される被検体の体表面の全てを検出する必要はなく、乳房領域やその周辺の領域に関係する体表面だけを検出すれば良い。本実施形態では、伏臥位MRI画像における乳房領域の中心位置は、マウス170やキーボード180を用いたユーザ入力により取得され、その中心位置より所定の範囲内を処理対象とする。
上記の方法で体表面401の検出は実行される。そして、その体表面である空気領域403と乳房領域402との境界の位置の集合である体表面形状を、本実施形態ではsp,surface,i (1≦i≦Np,surface)と表記する。ここで、Np,surfaceは体表面形状を構成する位置(点)の数である。同様にして、大胸筋面形状も画像処理によって検出する。本実施形態ではこれをsp,pectral,i (1≦i≦Np,pectral)と表記する。ただし、Np,pectralは大胸筋面形状を構成する位置(点)の数である。ここで、体表面形状および大胸筋面形状には夫々、これらを構成する点の間の接続情報が付随しているものとする。つまり、体表面形状および大胸筋面形状とは、その位置を表す複数の点(点群)の情報に加え、その点群が作る面に関する情報も同時に持つものとする。
次に解剖学特徴抽出部1020は、乳頭位置の検出を行う。乳頭位置は、上記の方法で検出した体表面形状をさらに処理することによって検出することができる。例えば、体表面形状の局所的な曲率を算出し、その曲率が最大となる位置を乳頭位置として検出するようにできる。または、体表面形状を構成する全ての位置のうち、MRI画像座標系においてY軸の座標値が最も小さい(最も腹側方向)位置を選択し、これを乳頭位置とするようにできる。また、伏臥位MRI画像を画像処理して乳頭位置を検出することもできる。本実施形態では、検出した乳頭位置を三次元座標値xp,surfaceと表記する。次に解剖学特徴抽出部1020は、大胸筋面上の基準位置の検出を行う。この処理は、例えば大胸筋面形状を構成する全ての位置の中で、乳頭位置に最も近い位置を選択することで実行される。本実施形態では、検出した大胸筋面上の基準位置の三次元座標値をxp,pectralと表記する。
次に解剖学特徴抽出部1020は、上記のように検出した乳頭位置xp,surfaceが原点となるように、伏臥位MRI画像座標系を座標変換する処理を実行する。具体的には、上記の処理で取得したsp,surface,i、sp,pectral,i、xp,surface、xp,pectralを-xp,surfaceだけ並進させる処理を実行する。以上の処理によりステップS310における解剖学特徴の抽出が行われる。
なお、上記の説明では解剖学特徴抽出部1020が伏臥位MRI画像を処理して解剖学特徴を抽出する方法について説明したが、この方法に限定されない。例えば、処理装置100は、伏臥位MRI画像をモニタ160に表示し、ユーザがマウス170やキーボード180によって解剖学特徴に関する情報を処理装置100に入力できるようにしても良い。また、画像処理によって抽出した解剖学特徴を、ユーザがマウス170やキーボード180を用いて修正・変更できるようにしても良い。また、処理装置100は、解剖学特徴のうち一部を画像処理によって抽出し、それ以外をユーザによる入力によって取得するようにしても良い。その時、画像処理によって抽出した解剖学特徴の一部をモニタ160に表示するようにしても良い。例えば、処理装置100は、体表面形状と大胸筋面形状を画像処理によって抽出し、その結果をモニタ160に表示する。そして、表示された体表面形状や大胸筋面形状を参照しながら乳頭位置と大胸筋面上の基準位置をマウス170やキーボード180を使ってユーザが処理装置100に入力するようにしても良い。
(ステップS320)伏臥位正規化座標系への変換を算出
ステップS320において、正規化部1040は、ステップS310で抽出した被検体の伏臥位における解剖学特徴に基づき、伏臥位における被検体の形状を基準形状に変換する正規化変換を導出する。具体的には、正規化部1040は、伏臥位MRI画像座標系から伏臥位正規化座標系への変換を表す情報として、これらの座標系間の座標変換関数を算出する処理を実行する。この変換は、MRI画像座標系における体表面および大胸筋面の夫々が、伏臥位正規化座標系において予め定める面上に位置するような変換である。また、この変換は、変換の前後において乳房領域における任意の構造がトポロジーの観点で、なるべく損なわれない変換である。上記の座標変換関数を算出するためにステップS320が実行する具体的な処理の手順について、図5のフローチャートを用いて詳しく説明する。
〈ステップS3200〉体表面の測地線距離を算出
ステップS3200において、正規化部1040は、ステップS310で抽出した解剖学特徴に基づき、伏臥位における被検体の体表面形状を構成する各位置における、乳頭位置を基準とした測地線距離を算出する処理を実行する。すなわち、正規化部1040は、体表面形状を構成する各位置のうち、乳頭位置については測地線距離を0とし、それ以外の任意の位置における乳頭からの測地線距離を算出する。なお、測地線距離を算出する方法は周知のいかなる方法を用いても良い。本実施形態における体表面形状は、体表面を構成する位置の集合と共に、それらの接続に関する情報も付随しているため、測地線距離を算出する方法として例えばダイクストラ法などを用いることができる。以上の処理により、体表面形状を構成する各位置の測地線距離dp,surface,i (1≦i≦Np,surface)を算出する。ここで、添え字のiは体表面形状sp,surface,i(1≦i≦Np,surface)の添え字iと共通であり、i番目の体表面形状の位置sp,surface,iにおける測地線距離をdp,surface,iとする。
〈ステップS3210〉体表面の方位を算出
ステップS3210において、正規化部1040は、ステップS310で抽出した解剖学特徴に基づき、伏臥位における被検体の体表面形状を構成する各位置における、乳頭位置を基準とした方位を算出する処理を実行する。ここで方位とは、例えばMRI画像座標系におけるX-Z平面上の方位とすることができる。この場合、体表面形状を構成する各位置の座標値sp,surface,i (1≦i≦Np,surface)のうち、X座標値xiとZ座標値ziとを用いて、数1に示す計算により方位ap,surface,i [rad]を算出するようにできる。
(数1)
ap,surface,i =tan-1(zi/ xi)
ここで、添え字のiは体表面形状sp,surface,i(1≦i≦Np,surface)の添え字iと共通であり、i番目の体表面形状の位置sp,surface,iにおける方位をap,surface,iとする。
なお、方位の算出方法は上記の方法に限らず、例えば以下の方法により算出できる。すなわち、乳頭位置と大胸筋面上の基準位置とを結ぶベクトルをY軸とし、被検体の体軸方向(足側から頭側へ向かう方向)をZ軸とする。ただし、Z軸がY軸と直交しない場合には補正を行う必要がある。そして、Y軸とZ軸との外積方向をX軸と設定したうえで、上記数1の算出を行うようにしてもよい。以上の方法によれば、MRI画像座標系において、MRI画像に描出されている被検体の姿勢が斜めに傾いているような場合であっても、被検体の姿勢を基準とした座標軸で方位を算出できる効果がある。
〈ステップS3220〉大胸筋面の測地線距離を算出〉
ステップS3220において、正規化部1040は、ステップS310で抽出した解剖学特徴に基づき、伏臥位における被検体の大胸筋面形状を構成する各位置sp,pectral,iにおける、大胸筋面上の基準位置xp,pectralを基準とした測地線距離dp,pectral,i (1≦i≦Np,pectral)を算出する処理を実行する。本処理ステップは、体表面を対象としたステップS3200と同様の処理を大胸筋面に適用することで実行される。
〈ステップS3230〉大胸筋面の方位を算出
ステップS3230において、正規化部1040は、ステップS310で抽出した解剖学特徴に基づき、伏臥位における被検体の大胸筋面形状を構成する各位置sp,pectral,iにおける、大胸筋面上の基準位置xp,pectralを基準とした方位ap,pectral,i (1≦i≦Np,pectral)を算出する処理を実行する。本処理ステップは、体表面を対象としたステップS3210と同様の処理を大胸筋面に適用することで実行される。
〈ステップS3240〉体表面を正規化座標系に座標変換
ステップS3240において、正規化部1040は、ステップS3200およびステップS3210で算出した体表面上の測地線距離および方位に基づいて、伏臥位における被検体の体表面形状を、伏臥位正規化座標系における所定の面へ座標変換する変換を求める処理を実行する。具体的には、正規化部1040は、伏臥位MRI画像座標系における体表面形状を構成する夫々の位置sp,surface,iに対応する、伏臥位正規化座標系における位置s’p,surface,iを算出する。
この処理について図6を用いて説明する。図6(a)は、伏臥位MRI画像400における乳房の模式図である。図6(b)は、伏臥位正規化座標系で表わされる伏臥位正規化空間における乳房の模式図である。図6(a)、図6(b)共に、紙面による説明の都合上、2次元の画像として例示するが、実際の処理では三次元の空間を持つ。例えば、体表面401は図中では曲線であるが、実際の処理では曲面である。同様に正規化体表面411は図中では直線であるが、実際の処理では平面である。本処理ステップでは、正規化部1040は、図6(a)のように、輪郭形状が曲面となる体表面401を、矩形形状の上側の平面である正規化体表面411に座標変換する処理を実行する。ここで、伏臥位正規化座標系には、正規化乳頭413の位置が予め所定の位置に定義されているものとする。本実施形態では一例として、正規化乳頭413の位置を伏臥位正規化座標系の原点に定義する場合について説明する。
体表面401から正規化体表面411への座標変換に関する具体的な処理手順について説明する。体表面401は、本実施形態でsp,surface,iと表記するNp,surface個の点の集合である。これら個々の点については、ステップS3200およびステップS3210の処理により測地線距離dp,surface,iと方位ap,surface,iが算出されている。本処理ステップでは、正規化部1040は、これらの算出結果に基づいて、伏臥位正規化座標系における対応する位置を算出する。具体的には、正規化部1040は、数2から数4の計算により座標値を算出する。
(数2)
xi= dp,surface,i・cos(ap,surface,i)
(数3)
yi=0
(数4)
zi= dp,surface,i・sin(ap,surface,i)
この計算は、正規化部1040が体表面401上の全ての点に以下の座標変換を施すことを意味する。すなわち、正規化部1040は、伏臥位MRI画像座標系における体表面401上の全ての点を、伏臥位正規化座標系においてy=0のx-z平面上に座標変換する。また、正規化部1040は、全ての点について、伏臥位正規化座標系における正規化乳頭413を基準とした距離及び方位が、伏臥位MRI画像座標系における乳頭406を基準とした測地線距離及び方位と一致するようにする。以上の処理により算出した、伏臥位正規化座標系における正規化体表面の位置をs’p,surface,i (1≦i≦Np,surface)と表記する。
〈ステップS3250〉大胸筋面を正規化座標系に座標変換
ステップS3250において、正規化部1040は、ステップS3220およびステップS3230で算出した大胸筋面上の測地線距離および方位に基づいて、伏臥位における被検体の大胸筋面形状を、伏臥位正規化座標系における所定の面へ座標変換する変換を求める処理を実行する。具体的には、正規化部1040は、伏臥位MRI画像座標系における大胸筋面形状を構成する夫々の位置sp,pectral,iに対応する、伏臥位正規化座標系における位置s’p,pectral,iを算出する。
この処理について図6を用いて説明する。本処理ステップでは、正規化部1040は、体表面を対象としたステップS3240の処理と同様にして、大胸筋面に対して処理を実行する。すなわち、正規化部1040は、図6(a)のように輪郭形状が曲面となる大胸筋面404を、矩形形状の下側の平面である正規化大胸筋面412に座標変換する処理を実行する。ここで、伏臥位正規化座標系には、正規化大胸筋面の基準点414の位置が予め所定の位置に定義されているものとする。本実施形態では一例として、正規化大胸筋面の基準点414を伏臥位正規化座標系の座標値(0,100,0)に定義する場合について説明する。
具体的な処理は、正規化部1040が、大胸筋面404上の全ての点sp,pectral,iについて、数5から数7の計算を行う事で実行される。すなわち、正規化部1040は、正規化大胸筋面の基準点414と同一のx-z平面上に全ての点を座標変換する。このとき、正規化部1040は、全ての点について、正規化大胸筋面の基準点414を基準とした距離及び方位が、伏臥位MRI画像座標系における大胸筋面の基準点からの測地線距離dp,pectral,i及び方位ap,pectral,iと一致するようにする。
(数5)
xi= dp,pectral,i・cos(ap,pectral,i)
(数6)
yi=100
(数7)
zi= dp,pectral,i・sin(ap,pectral,i)
以上の処理により算出した、伏臥位正規化座標系における正規化大胸筋面の位置をs’p,pectral,i (1≦i≦Np,pectral)と表記する。
〈ステップS3260〉正規化変形を算出
ステップS3260において、正規化部1040は、伏臥位MRI画像座標系から伏臥位正規化座標系への変換を表す情報として、当該座標系間の座標変換関数(変形場)を算出する処理を実行する。すなわち、正規化部1040は、ステップS3240およびステップS3250で求めた体表面および大胸筋面の伏臥位正規化座標系への離散的な座標変換の結果群を空間的に補間して、伏臥位MRI画像座標系から伏臥位正規化座標系への密な変換を算出する。この処理は、具体的には放射基底関数やB-スプラインなどを用いた周知の補間方法により実現できる。本処理ステップで算出した伏臥位MRI画像座標系から伏臥位正規化座標系への変換関数を、本実施形態ではφp(x)と表記する。ただし、φp(x)は伏臥位MRI画像座標系における位置座標値を引数とし、それに対応する伏臥位正規化座標系における位置座標値を返す関数である。
また正規化部1040は、φp(x)とは逆に、伏臥位正規化座標系における位置座標値を引数として、それに対応する伏臥位MRI画像座標系の位置を返す関数も同様に算出する。本実施形態では、これをφp -1 (x)と表記する。ここでφp -1 (x)は、伏臥位正規化座標系における所定の矩形領域において定義されるものとする。前記矩形領域は、例えば、s’p,surface,iおよびs’p,pectral,iの全てを内包する矩形領域である。
なお、上記の方法で算出した変換関数φp(x)およびφp -1 (x)は、数8から数11に示す性質を備えている。
(数8)
s’p,surface,i≒φp(sp,surface,i)
(数9)
s’p,pectral,i≒φp(sp,pectral,i)
(数10)
sp,surface,i≒φp -1 (s’p,surface,i)
(数11)
sp,pectral,i≒φp -1 (s’p,pectral,i)
以上に説明したステップS3200からステップS3260の処理により、本実施形態におけるステップS320の処理が実行される。
(ステップS340)仰臥位MRI画像取得
ステップS340において、画像取得部1000は、被検体の乳房を仰臥位の体位で撮像したMRI画像(仰臥位MRI画像)を、データサーバ120から処理装置100へ取り込む処理を実行する。この処理は伏臥位のMRI画像を対象としたステップS300と同様の手順で実行できるため、詳細な説明は省略する。取得した仰臥位MRI画像は、仰臥位MRI画像座標系における三次元のボリュームデータである。
(ステップS350)仰臥位MRI画像から解剖学特徴を抽出
ステップS350において、解剖学特徴抽出部1020は、ステップS340で取得した仰臥位MRI画像を処理することにより、被検体の仰臥位における解剖学特徴を抽出する処理を実行する。この処理は伏臥位MRI画像を対象としたステップS310と同様の処理を仰臥位MRI画像に適用して実行できるため、詳細な説明は省略する。なお、以下では、本処理ステップにおいて抽出した、仰臥位における体表面形状をss,surface,i (1≦i≦Ns,surface)、大胸筋面形状をss,pectral,i (1≦i≦Ns,pectral)、乳頭位置をxs,surface、大胸筋面の基準点をxs,pectralと表記する。
(ステップS360)仰臥位正規化座標系への変換の算出
ステップS360において、正規化部1040は、ステップS350で抽出した被検体の仰臥位における解剖学特徴に基づき、仰臥位における被検体の形状を基準形状に変換する正規化変換を導出する。具体的には、正規化部1040は、仰臥位MRI画像座標系から仰臥位正規化座標系への変換を表す情報として、これらの座標系間の座標変換関数を算出する処理を実行する。この処理は伏臥位の解剖学特徴を対象としたステップS320と同様の手順を仰臥位の解剖学特徴に適用して実行できるため、詳細な説明は省略する。なお、以下では、本処理ステップにおいて取得した仰臥位MRI画像座標系から仰臥位正規化座標系への変換関数をφs(x)と表記する。また、仰臥位正規化座標系から仰臥位MRI画像座標系への変換関数をφs -1 (x) と表記する。
(ステップS380)MRI画像を変形
ステップS380において、画像変形部1060は、ステップS320およびステップS360の処理結果に基づいて、伏臥位MRI画像を仰臥位に変形した変形画像を生成する処理を実行する。具体的には、画像変形部1060は、伏臥位MRI画像を構成する全てのボクセル(ボリュームデータを構成する画素)に対して、そのボクセルの位置xpを数12の計算により座標変換し、変換後の位置xdを算出する。
(数12)
xds -1psp(xp)}]
そして、画像変形部1060は、その変換後の位置xdに輝度値Ip(xp)を持つボリュームデータを生成する。ただし関数φps(x)は伏臥位正規化座標系から仰臥位正規化座標系の間の任意の変換関数であり、ここでは数13に示す恒等関数とする。また、変換関数φps(x)の定義域はφp -1 (x)と同じ矩形領域である。
(数13)
x=φps(x)
上記の処理で生成したボリュームデータを変形MRI画像Id(x)と表記する。つまり、変形MRI画像Id(x)は、数14に示す計算により伏臥位MRI画像Ip(x)から算出される。
(数14)
Id(x)= Ipp -1ps -1s(x))}]
本処理ステップにおいて実行される数12の計算は、具体的には以下の事を意味する。すなわち、伏臥位MRI画像座標系における位置xpは伏臥位正規化座標系へと変換され、伏臥位正規化座標系と仰臥位正規化座標系とを同一視した上で、その座標値が仰臥位MRI画像座標系に変換されている。つまり、伏臥位と仰臥位の間の形状の差異が夫々の正規化によって打ち消される(解剖学的に同一の点が正規化によって正規化座標系の略一致した座標に写像される)と仮定し、その仮定に基づき、伏臥位MRI画像座標系と仰臥位MRI画像座標系との間の変換を求められている(変形位置合わせをしている)ことを意味している。
なお、上記の説明では、伏臥位正規化座標系と仰臥位正規化座標系との間の座標変換φps(x)として恒等関数を用いる場合を例として説明したが、任意の座標変換であって良い。例えば、φps(x)は変形MRI画像と仰臥位MRI画像との間の画像間類似度を高める変形関数であって良い。φps(x)の一例としては、非線形座標変換の代表的な手法の一つであるFFD(Free Form Deformation)を用いて表現するようにできる。この場合、数15に示す画像間類似度を最大化するようにFFDの変形パラメータの最適化する処理を実行する。
(数15)
ただし、Is(x)は仰臥位MRI画像であり、Ωは仰臥位MRI画像座標系における乳房領域である。FFDの変形パラメータの最適化は最急勾配法など周知の非線形最適化方法によって実行される。また、画像間類似度は、数15に示す計算方法以外にも相互相関や相互情報量を用いる方法など、周知のいかなる画像間類似度の計算方法であって良い。上記の方法によれば、MRI画像に描出されている被検体の体表面や大胸筋面などの形状情報以外に、伏臥位MRI画像および仰臥位MRI画像の輝度値の類似度が高い変形を生成することができる。このため、伏臥位MRI画像と仰臥位MRI画像との位置合わせを、より高精度に実行できる効果がある。
(ステップS390)変形画像を表示
ステップS390において、観察画像生成部1080は、ステップS380で生成した変形画像Id(x)と、仰臥位MRI画像Is(x)とを並べた観察画像を生成する。この時、観察画像生成部1080は、変形画像Id(x)と仰臥位MRI画像Is(x)の夫々をユーザの操作に応じて任意の平面で切断した画像(断面画像)を切り出し、これらを並べて観察画像を構成するようにできる。また、観察画像生成部1080は、変形画像Id(x)と仰臥位MRI画像Is(x)の夫々をボリュームレンダリングした画像を並べて観察画像を構成するようにしても良い。また、観察画像生成部1080は、変形画像Id(x)と仰臥位MRI画像Is(x)の夫々から生成した断面画像を重畳または融合して観察画像を構成するようにしても良い。そして、観察画像生成部1080は、以上の処理で生成した観察画像をモニタ160に表示する処理を実行する。
以上に説明した方法により、第1実施形態による処理装置100の処理が実施される。本実施形態によれば、乳房の変形に関する生体力学的な特性を考慮した正規化を行っているため、変形に伴う位置の変動を概ね吸収し、解剖学的に共通する空間へと変換できる。そのため、伏臥位および仰臥位で撮像された被検体の乳房の夫々を解剖学的に略一致した空間に写像する正規化変換を施す事ができる。また、伏臥位および仰臥位で撮像された被検体の乳房MRI画像を簡易な方法で変形位置合わせできる仕組みを提供できる。
(第1実施形態の変形例)矩形形状以外の正規化空間
第1実施形態の説明では、ステップS320の処理において、正規化部1040が図6(a)に示す乳房領域402が、図6(b)に示す矩形形状となる正規化変形を生成する場合を例として説明したが、形状はこの場合に限定されない。例えば、正規化座標系における乳房領域は矩形形状以外の形状でも良い。一例としては、二次曲面などの任意の幾何曲面で囲まれた形状であっても良い。また、処理装置100に予め複数の形状情報を備え、その中からユーザが任意に選択できるようにしても良い。この場合、処理装置100は、例えば伏臥位MRI画像と仰臥位MRI画像とをモニタ160に表示するなどしてユーザに提示し、ユーザは、その画像を観察しながら適切な形状情報を選択できるようにしても良い。この方法によれば、被検体毎に多様な形状を持ちうる乳房の特性に対して、正規化変形の方法を適応的に選択でき、より精度の良い位置合わせが行える効果がある。
[第2実施形態]曲断面表示
第1実施形態の説明では、ステップS380の処理において、画像変形部1060が伏臥位MRI画像を変形した変形MRI画像を生成する場合を例に説明したが、生成するMRI画像はこの例に限定されない。本発明の第2実施形態では、画像変形部1060が、伏臥位MRI画像を伏臥位正規化座標系に変換した伏臥位正規化画像と、仰臥位MRI画像を仰臥位正規化座標系に変換した仰臥位正規化画像とを生成し、これらの画像を並べて表示する場合について説明する。
以下、第1実施形態との相違点についてのみ説明する。なお、本実施形態における処理装置の機能構成は図1と同様であり、画像変形部1060の機能のみが異なっている。本実施形態における画像変形部1060は、伏臥位正規化画像Ipdと仰臥位正規化画像Isdを生成する。また本実施形態における処理装置の処理手順は図3と同様であり、ステップS380およびステップS390の処理内容のみが異なっている。以下、本実施形態におけるステップS380およびステップS390の処理について説明する。
ステップS380において、画像変形部1060は、伏臥位MRI画像Ipを変換関数φp -1 (x)に基づいて伏臥位正規化座標系に変換した伏臥位正規化画像Ipdを生成する。図7は、この処理によって生成される画像の具体例を示す図である。図7(a)において、伏臥位MRI画像400はIpを模式的に表した例である。画像変形部1060は、この伏臥位MRI画像400を変換関数φp -1 (x)に基づいて変形させることで、図7(b)に示すように、伏臥位正規化座標系における画像として伏臥位正規化画像410、すなわちIpdを生成する。またステップS380において画像変形部1060は、仰臥位MRI画像Isについても同様の処理を行い、図7(c)の仰臥位MRI画像420をφs -1 (x)に基づいて変形させ、図7(d)の仰臥位正規化画像430、すなわちIsdを生成する。以上の処理により、三次元のボリュームデータである伏臥位正規化画像Ipdおよび仰臥位正規化画像Isdが生成される。
ステップS390において、観察画像生成部1080は、ステップS380で生成した伏臥位正規化画像Ipdおよび仰臥位正規化画像Isdに基づいて、観察画像を生成する処理を実行する。具体的には、観察画像生成部1080は、伏臥位正規化画像Ipdと仰臥位正規化画像Isdの双方を y=任意の定数 となるx-z平面で切断した断面画像I slice,pdおよびI slice,sdを生成し、これらの断面画像を並べた画像を観察画像として生成する。ここで、伏臥位正規化座標系および仰臥位正規化座標系において、Y座標の値は概ね体表面および大胸筋面との間の距離の比率を表している。従って、断面画像I slice,pdおよびI slice,sdには、伏臥位MRI画像座標系および仰臥位MRI画像座標系において、体表面と大胸筋面との距離が一定の割合となる曲面形状上の断面(曲断面)が切り出される。そのため、両画像の医学的・解剖学的な観点での高度な比較が容易に行える効果がある。例えば、乳房の体表面付近に存在する表在血管の走行や、乳房領域における乳腺の広がり方などを、単一の断面画像で広い範囲に渡って比較観察できる効果がある。
(第2実施形態の変形例)
上記のようにMRI画像の曲断面を生成して表示する処理は、必ずしも伏臥位MRI画像と仰臥位MRI画像との両方に対して実行する必要は無い。例えば、処理装置100は、伏臥位MRI画像または仰臥位MRI画像の何れか一方を入力とし、第1実施形態のステップS300からステップS320の処理を実行し、上記の曲断面を生成する処理を実行して表示するようにしても良い。この方法によれば、表在血管の走行や乳腺の広がり方など医学的・解剖学的な観点での高度な画像観察を行うことができる効果がある。
[第3実施形態]統計変形モデル生成、位置合わせ
本実施形態による処理装置200は、第1実施形態と同様の正規化処理を多症例(学習症例)に対して適用し、正規化された各症例の変形を統計処理することにより、変形に関するモデル(統計変形モデル)を構築する。ここで、異なる症例間における正規化について説明する。少なくとも正常人体においては、個体の違いによらず、トポロジー的な観点でほぼ同一の解剖学的構造を有する。また、個体の違いは概ねスケール変換(相似変換)によって吸収できる。ただし、乳房については発達のメカニズムに起因し、前述のスケール変換による個体差の吸収には限りがある。一方、正常人体において、乳頭、体表面、大胸筋面、正中線、頭尾方向(体軸)などは解剖学的に全ての個体で共通する特徴的な幾何構造である。
本実施形態による処理装置200は、上記の特性を考慮して、個体間のスケール変換に加え、前記特徴的な幾何構造で基準化した空間を考え、異なる個体の人体の乳房を、その空間に座標変換する。これにより、個体間の違いを概ね吸収し、解剖学的に共通する空間へと変換できる。本実施形態では、Nsamples個の症例の伏臥位および仰臥位のMRI画像を学習症例とする場合を例として説明する。さらに処理装置200は、学習症例とは異なる未知の症例(対象症例)に対して統計変形モデルの当てはめを行うことにより、変形の推定を行う。これにより、対象症例の伏臥位および仰臥位のMRI画像間の変形位置合わせを高精度に行う。
(機能構成)
図8は、本実施形態による処理システムの構成を示す図である。本実施形態において、第1実施形態と同様の機能を持つ構成要素には図1と同一の番号を付しており、説明は省略する。同図に示すように、本実施形態における処理装置200は、画像取得部1000、解剖学特徴抽出部1020、正規化部1040、学習症例変形生成部1220、スケール算出部1230、統計変形モデル生成部1240、対象症例変形生成部1420、画像変形部1060、観察画像生成部1080を持つ。
学習症例変形生成部1220は、学習症例の夫々のMRI画像に基づき、伏臥位から仰臥位への変形を生成する。スケール算出部1230は、学習症例の夫々に関するスケールの算出を行う。統計変形モデル生成部1240は、学習症例の夫々に関する変形・スケールに基づいて統計変形モデルを生成する。対象症例変形生成部1420は、対象症例の伏臥位・仰臥位間の変形を生成する。
(処理フロー)
次に、処理装置200が行う全体の動作について説明する。本実施形態では、主メモリ212に格納されている各部の機能を実現するプログラムをCPU211が実行することにより実現される。また、以下に説明する処理装置200が行う各処理の結果は、主メモリ212に格納することにより記録される。
本実施形態における処理装置200の処理は、学習フェーズの処理と変形推定フェーズの処理から成り、まず学習フェーズの処理を実行し、その後に変形推定フェーズの処理を実行する。学習フェーズの処理では、多症例の伏臥位および仰臥位のMRI画像間の変形を学習し、統計変形モデルを生成する処理が実行される。変形推定フェーズの処理では、学習フェーズで算出した統計変形モデルを用いて対象症例の伏臥位、仰臥位間の変形位置合わせが実行される。
以下に記載する本実施形態の説明では、処理装置200が前記学習フェーズの処理と変形推定フェーズの処理の両方を実行する場合を例として説明するが、学習フェーズの処理と、変形推定フェーズの処理を異なる処理装置で実行するようにしても良い。また、本実施形態による処理装置200は、学習フェーズの処理と変形推定フェーズの処理の両方の処理を実行する場合に限らず、例えば、学習フェーズの処理だけを実行するようにしても良い。また、学習フェーズの処理の結果として得られる統計変形モデル自体の提供も、本実施形態に含まれる。
(学習フェーズの処理)
図9は本実施形態における処理装置200が行う、学習フェーズの処理の手順を説明するフローチャートである。以下、このフローチャートに示す処理手順に従って、本実施形態の学習フェーズの処理について詳しく説明する。
(ステップS500)学習症例の画像を算出
ステップS500において、画像取得部1000は、Nsamples個の学習症例の伏臥位MRI画像と仰臥位MRI画像を取得する処理を実行する。この処理は第1実施形態のステップS300およびS340と同様の処理を、Nsamples個の学習症例の夫々に対して適用することで実行できる。詳細な説明は省略する。
(ステップS510)学習症例の解剖学特徴を抽出
ステップS510において、解剖学特徴抽出部1020は、ステップS500で取得した学習症例の伏臥位MRI画像と仰臥位MRI画像の夫々を処理することにより、解剖学特徴を抽出する処理を実行する。この処理は第1実施形態のステップS310およびS350と同様の処理を、Nsamples個の学習症例の夫々に対して適用することで実行できる。詳細な説明は省略する。
(ステップS520)学習症例の正規化座標系への変形関数を算出
ステップS520において、正規化部1040は、ステップS510で抽出した学習症例の解剖学特徴に基づき、学習症例の夫々に関して、被検体の形状を基準形状に変換する正規化変換を導出する。具体的には、正規化部1040は、伏臥位MRI画像の夫々について、伏臥位MRI画像座標系から伏臥位正規化座標系への変形関数を算出する処理を実行する。また、正規化部1040は、仰臥位MRI画像の夫々に関して、仰臥位MRI画像座標系から仰臥位正規化座標系への変形関数を算出する処理を実行する。これらの処理は、第1実施形態で説明したステップS320およびS360と同様の処理で実行できる。詳細な説明は省略する。
以上の処理により算出された学習症例の伏臥位MRI画像座標系から伏臥位正規化座標系の変形関数をφp,j(x)と表記する。また、仰臥位MRI画像座標系から伏臥位正規化座標系の変形関数をφs,j(x)と表記する。ここで、jは学習症例の症例番号を意味し、1≦j≦Nsamplesである。すなわち、Nsamples個の学習症例の夫々に関して変形関数を求める処理を実行する。
(ステップS540)学習症例の伏臥位正規化座標系から仰臥位正規化座標系への変換関数を算出
ステップS540において、学習症例変形生成部1220は、学習症例の夫々に関して、伏臥位正規化座標系から仰臥位正規化座標系への変換関数を算出する処理を実行する。具体的には、学習症例変形生成部1220は、以下に説明する処理をNsamples個の学習症例の夫々について実行する。
まず、学習症例変形生成部1220は、伏臥位MRI画像座標系および仰臥位MRI画像座標系の双方で対応する位置を、マウス170やキーボード180を用いたユーザ入力等により取得する。ここでは、この対応する位置が、ユーザから入力されたものとする。入力された伏臥位MRI画像座標系および仰臥位MRI画像座標系における対応する位置を夫々、xp,corres,k、xs,corres,kとする。ただしkは対応点のインデックスであり、Ncorresを対応する位置の数としたときに1≦k≦Ncorresである。ここで、対応する位置は、例えばMRI画像における血管の分岐点や、乳腺の特徴的な構造を持つ部位などであり、ユーザが目視により対応が与えられる位置の組である。
次に、学習症例変形生成部1220は、取得したxp,corres,k、xs,corres,kの夫々を、ステップS520で取得した正規化座標系への変換を用いて変換する。具体的には、学習症例変形生成部1220は、数16の処理を実行し、伏臥位正規化座標系および仰臥位正規化座標系における位置x’p,corres,k、x’s,corres,kを算出する。
(数16)
x’p,corres,kp,j(xp,corres,k) (1≦k≦Ncorres)
(数17)
x’s,corres,ks,j(xs,corres,k) (1≦k≦Ncorres)
さらに本処理ステップでは、学習症例変形生成部1220は、算出したx’p,corres,k、x’s,corres,kに基づいて、伏臥位正規化座標系から仰臥位正規化座標系への変換関数φps,j(x)を算出する。具体的には、学習症例変形生成部1220は、数18の関係を最小の誤差で近似するように変換関数φps,j(x)を算出する。
(数18)
x’s,corres,kps,j(x’p,corres,k) (1≦k≦Ncorres)
なお、数18の誤差は、例えばNcorres個の対応する位置に関する二乗和誤差を用いることができる。ここで、変換関数φps,j(x)は、伏臥位正規化座標系において定義される連続関数であり、具体的にはFFD(Free Form Deformation)やRBF(Radial Basis Function)等を用いて表現するようにできる。
この時、学習症例変形生成部1220は、伏臥位正規化座標系における体表面の位置(本実施形態ではy=0の平面)における変換関数φps,j(x)の値が、仰臥位正規化座標系における体表面の位置(同じくy=0の平面)を出力するように所定の拘束を設けることが望ましい。同様に、学習症例変形生成部1220は、伏臥位正規化座標系における大胸筋面の位置(本実施形態ではy=100の平面)における変換関数φps,j(x)の値が、仰臥位正規化座標系における大胸筋面の位置(同じくy=100の平面)を出力するように所定の拘束を設けることが望ましい。これにより、学習症例変形生成部1220は、伏臥位と仰臥位との間の位置の関係において、体表面同士および大胸筋面同士とが合致するという条件と、ユーザが入力した乳房内部などの対応する位置の情報との両方を考慮した変換関数を求めることができる。本処理ステップでは、学習症例変形生成部1220は、以上に説明した処理を、Nsamples個の学習症例の夫々について実行し、各症例について伏臥位正規化座標系から仰臥位正規化座標系への変換関数φps,j(x)を算出する。
なお、上記の説明では、学習症例変形生成部1220は、伏臥位MRI画像座標系および仰臥位MRI画像座標系の双方で対応する位置の情報を取得し、それに基づいて変形関数φps,j(x)を算出する場合を例として説明したが、変形関数φps,j(x)の算出方法はこの例に限定されない。例えば、第1実施形態のステップS380の処理として説明したことと同様に、学習症例変形生成部1220は、伏臥位MRI画像と仰臥位MRI画像の画像間類似度に基づいて変形関数φps,j(x)を算出するようにしても良い。また、学習症例変形生成部1220は、伏臥位MRI画像座標系および仰臥位MRI画像座標系の双方で対応する位置の情報に加えて、伏臥位MRI画像と仰臥位MRI画像の画像間類似度に基づいて変形関数φps,j(x)を算出するようにしても良い。上記の方法によれば、より正確に学習症例の変形を取得できる効果がある。
(ステップS550)学習症例のスケールの算出
ステップS550において、スケール算出部1230は、学習症例の夫々に関して、症例のスケールを算出する処理を実行する。ここで症例のスケールとは、症例毎に異なる乳房領域の大きさを表す数値である。スケールの算出方法は、例えば、伏臥位における被検体の、乳頭位置と、乳頭位置から最も近い正中線直上の体表位置との間の距離の値を計測する事で算出する。この処理は、ユーザが、伏臥位における被検体に対して直接的に測定器具を用いて計測した数値を処理装置200に入力するようにしても良い。または、処理装置200が、伏臥位MRI画像の上で上記の距離の値を計測できるようにしても良い。この時、処理装置200は、伏臥位MRI画像を自動処理することで上記計測値を算出するようにしても良い。または、処理装置200は、モニタ160などを用いてユーザに伏臥位MRI画像を提示し、ユーザによるマウス170やキーボード180の操作により計測値を取得できるようしても良い。一例としては、処理装置200は、モニタ160上に伏臥位MRI画像を表示し、画像上に描出されている被検体の乳頭位置と、乳頭位置から最も近い正中線直上の体表位置とをユーザに指定させ、その間の距離を算出することで実現できる。
上記の例では、乳頭位置と、乳頭位置から最も近い正中線上の体表位置との間の距離(ユークリッド距離)を用いて症例のスケールを算出する方法について説明したが、スケールを算出する方法はこの方法に限定されない。例えば、上記二点の間の測地線距離を用いてスケールを算出するようにしても良い。これによれば、症例毎の乳房の形状の違いも考慮したスケール値を算出できる効果がある。また、スケールの算出は上記二点間の距離または測地線距離に限らず、例えば乳房領域の体積や、乳房の外端までの距離、被検体の胸囲などに基づいて算出するようにしても良い。また、スケールを算出する方法は一つとは限らず、例えば複数の種類の方法で算出した値に基づいて算出しても良い。この場合、複数の種類で算出した値をベクトル化した多次元のスケール値としても良いし、複数の種類で算出した値の平均演算や線形結合演算などによってスカラ値としてスケール値を算出するようにしても良い。何れにしても、Nsamples個の複数症例に対して同一の方法・基準でスケールを算出することが望ましい。
以上に説明した方法により算出したスケールを vj (1≦j≦Nsamples)と表記する。本実施形態では、伏臥位における乳頭位置と乳頭位置から最も近い正中線上の体表位置との間の距離(ユークリッド距離)であるスカラ値と所定の基準値との比率をスケール値とする。ここで、所定の基準値は、例えば標準的な乳房における乳頭位置と乳頭位置から最も近い正中線上の体表位置との間の距離の値である。
(ステップS580)統計変形モデルの生成
ステップS580において、統計変形モデル生成部1240は、ステップS500からステップS570の処理によって算出したNsamples個の被検体に関するφps,j(x)、vjに基づいて統計変形モデルを生成する処理を実行する。図10は、ステップS580の処理をさらに詳しく説明するフローチャートである。以下、図10のフローチャートに沿って説明する。
〈ステップS5800〉φps,j(x)のスケーリング
ステップS5800において、統計変形モデル生成部1240は、Nsamples個の症例に関する変換関数φps,j(x)と、スケールvjとに基づいて、変換関数φps,j(x)をスケーリング処理した変換関数φ’ps,j(x)を算出する処理を実行する。具体的には数19のようにして変換関数φ’ps,j(x)を算出する。
(数19)
φ’ps,j(x’)= φps,j(x)/ vj
ただし、ここでx=( x , y , z )Tとした場合に、x’=( x×vj , y , z×vj )Tである。つまり、φ’ps,j(x)はφps,j(x)を、X座標およびZ座標に関して当該症例のスケール値vjでスケーリングした関数である。なお、φ’ps,j(x)の定義域はφps,j(x)の定義域と比較して、X座標およびZ座標に関して当該症例のスケール値vj分だけ縮小した領域である。統計変形モデル生成部1240は、以上の処理をNsamples個の全症例の変換関数φps,j(x) (1≦j≦Nsamples)について実行し、各症例のスケーリング処理した変換関数φ’ps,j(x) (1≦j≦Nsamples)を算出する。
〈ステップS5820〉φ’ps,j(x)のベクトル化
ステップS5820において、統計変形モデル生成部1240は、ステップS5800で算出したスケーリング処理後の変換関数φ’ps,j(x)を離散化する処理を実行する。この離散化の処理は以下の手順により実行される。
まず、統計変形モデル生成部1240は、Nsamples個の症例に関する変換関数φ’ps,j(x)の共通の定義域を求める。Nsamples個の変換関数φ’ps,j(x)の夫々は、各症例の体表面や大胸筋面の形状、スケール値などによって異なる定義域を持つが、本実施形態では、Nsamples個の全症例においてφ’ps,j(x)が定義されている領域を共通の定義域とする。従って、この定義域内ではNsamples個の全ての変換関数φ’ps,j(x)が値を持つ。
次に、統計変形モデル生成部1240は、変換関数φ’ps,j(x)の値を上記の共通の定義域に渡ってサンプリングし、そのサンプリング結果を縦に並べた離散化ベクトルを生成する。ここで、離散化ベクトルの生成は、上記の定義域内を所定の間隔でラスタスキャン状にサンプリングした値を順に並べたベクトルとする。なお、変換関数φ’ps,j(x)はx、y、zの三次元の値を返す関数であるため、統計変形モデル生成部1240は、離散化ベクトルを各座標軸毎に生成する。ここでは、x、y、zの三次元の各座標軸毎の離散化ベクトルをpx,j、py,j、pz,j、とする。離散化ベクトルは上記のラスタスキャンによるサンプリングの回数分の次元をもつ実ベクトルである。
統計変形モデル生成部1240は、以上に説明した処理をNsamples個の全ての変換関数φ’ps,j(x)について適用する。これによりNsamples個の離散化ベクトルpx,j、py,j、pz,jを取得する。なお、上記のラスタスキャン状のサンプリングにより得た離散化ベクトルは、任意の補間関数などによって、実空間の変換関数に逆変換できる。本実施形態では、離散化ベクトルpx,j、py,j、pz,jと実空間の位置xとを引数とする補間関数finterp(px,j , py,j , pz,j , x)により、φ’ps,j(x)が近似されるものとする。
〈ステップS5840〉主成分分析による統計変形モデル生成
ステップS5840において、統計変形モデル生成部1240は、ステップS5820で算出した離散化ベクトルpx,j、py,j、pz,j (1≦j≦Nsamples) を主成分分析することにより、統計変形モデルを生成する処理を実行する。主成分分析は周知の方法により実行できるため、ここでは詳細な説明は省略する。主成分分析の結果、統計変形モデル生成部1240は、各離散化ベクトルに関する平均ベクトルeave,x、eave,y、eave,zと、固有ベクトル ex,k、ey,k、ez,k (1≦k≦Nmode)を得る。ここで、Nmodeは主成分分析で算出する固有ベクトルの総数であり、例えば主成分分析により算出される累積寄与率に所定の閾値を設定することにより設定できる。本実施形態では、上記の処理により算出した平均ベクトルおよび固有ベクトルを統計変形モデルと呼ぶ。
なお、この平均ベクトルおよび固有ベクトルは、これらを適切な加重で線形和演算を行う数20の計算を行うことで、離散化ベクトルを近似することができる。
(数20)
ここで、Exはex,kを横に並べた行列である。また、bはbkを縦に並べたベクトルであり、本実施形態ではこれを係数ベクトルと呼ぶ。一般的に、症例数Nsamplesが十分に大きい場合には、それらの変形を表す離散化ベクトルは、Nsamplesより少数の固有ベクトルおよび平均ベクトルの線形和によって高い精度で近似できることが既に知られている。
以上に説明したステップS500からステップS580の処理により、本実施形態の学習フェーズの処理が実行される。この処理の結果、統計変形モデルが生成される。
(変形推定フェーズの処理)
図11は、本実施形態における処理装置200が行う、変形推定フェーズの処理の手順を説明するフローチャートである。以下、このフローチャートに示す処理手順に従って、本実施形態の変形推定フェーズの処理について詳しく説明する。
(ステップS600)統計変形モデル読み込み
ステップS600において、処理装置200は、学習フェーズの処理によって生成した統計変形モデルを、処理装置200の主メモリ212に読みだす処理を実行する。
(ステップS602)対象症例の画像取得
ステップS602において、画像取得部1000は、対象症例の伏臥位および仰臥位のMRI画像を取得する処理を実行する。この処理は第1実施形態のステップS300およびステップS340と同様の処理で実行できる。詳細な説明は省略する。
(ステップS604)対象症例の解剖学特徴を抽出
ステップS604において、解剖学特徴抽出部1020は、ステップS602で取得した対象症例の伏臥位および仰臥位のMRI画像を処理して、対象症例の解剖学特徴を抽出する処理を実行する。この処理は第1実施形態のステップS310およびステップS350の処理と同様の処理で実行できる。詳細な説明は省略する。
(ステップS610)対象症例の正規化座標系への変換を算出
ステップS610において、正規化部1040は、ステップS604で抽出した対象症例の解剖学特徴に基づき、対象症例の伏臥位MRI画像座標系から伏臥位正規化座標系への変換、および仰臥位MRI画像座標系から仰臥位正規化座標系への変換を算出する処理を実行する。この処理は、学習フェーズの処理として説明したステップS520の処理と同様の処理を対象症例に適用して実行される。詳細な説明は省略する。この処理により算出した伏臥位MRI画像座標系から伏臥位正規化座標系への変換をφp,target(x)とする。また、仰臥位MRI画像座標系から仰臥位正規化座標系への変換をφs,target(x)とする。
(ステップS630)対象症例のスケール算出
ステップS600において、スケール算出部1230は、対象症例のスケールを算出する処理を実行する。この処理は学習フェーズの処理として説明したステップS550の処理と同様の処理を対象症例に適用して実行する。詳細な説明は省略する。この処理により算出した対象症例のスケールをvtargetとする。
(ステップS640)統計変形モデルの係数の最適化
ステップS600において、対象症例変形生成部1420は、対象症例の伏臥位MRI画像座標系と仰臥位MRI画像座標系との間の変換を算出する処理を実行する。すなわち、対象症例変形生成部1420は、伏臥位MRI画像と仰臥位MRI画像の間の変形位置合わせ処理を実行する。この処理は学習フェーズの処理で取得した統計変形モデルと、対象症例の伏臥位および仰臥位のMRI画像に基づいて実行される。具体的には、数21で示す計算によって算出される評価関数G(b)を最大化するような係数ベクトルbを算出する。
(数21)
G(b)=Gsimil{ D(Ip,b) , Is }
ここで、Ipは対象症例の伏臥位MRI画像、Isは対象症例の仰臥位MRI画像である。また、関数Gsimil(I1,I2)は、引数として与えられた2つの画像間の類似度を評価する関数であり、例えばSSDやSAD、相互相関、相互情報量などの周知の画像間類似度評価方法により実現できる。
また、関数D(I,b)は、統計変形モデルの係数ベクトルbに基づいて画像Iを変形させる関数である。関数D(I,b)は、より具体的には以下の処理を行う。すなわち、係数ベクトルbに基づき、数22で示す計算により離散化ベクトルpx、py、pzを算出する。
(数22)
そして、算出された離散化ベクトルpx、py、pzに基づいて数23によって変形関数φtarget(x)を求める。
(数23)
φtarget(x)= φs -1[ finterp{ px , py , pz , φp(x) } ]
さらに、ステップS630で算出した対象症例のスケールvtargetを用いて変換関数φtarget(x)をスケーリングした変換関数φ’target(x)を数24により算出する
(数24)
φ’target(x’)= φtarget(x)×vtarget
ただし、ここで、x=( x , y , z )Tとした場合に、x’=( x×vtarget , y , z×vtarget )Tである。つまり、φ’target(x)はφtarget(x)を、X座標およびZ座標に関して対象症例のスケール値vjでスケーリングし、さらに関数の値についてスケール値vtargetでスケーリングした関数である。
そして、関数D(I,b)は、この変形関数φ’target(x)に基づいて画像Iを変形させる。つまり、数21におけるD(Ip,b)は、以下の数25の計算を行うことになる。
(数25)
D(Ip,b)= Ip{φ’target(x) }
以上に説明したように、数21に示した評価関数G(b)は、係数ベクトルbに基づいて、伏臥位MRI画像を変形させた画像と、仰臥位MRI画像との類似度を評価する。本処理ステップでは、対象症例変形生成部1420は、最急降下法や準ニュートン法、共役勾配法などの非線形最適化の方法を用いて、数21に示す評価関数の値を最大化する係数ベクトルbを算出する処理を実行する。本処理により得た係数ベクトルをboptと表記する。
上記の説明では、伏臥位MRI画像を変形した変形MRI画像と、仰臥位MRI画像との間の画像間類似度に基づいて係数ベクトルを算出する場合を例として説明したが、係数ベクトルの算出方法はこの例に限定されない。例えば、対象症例の伏臥位MRI画像と仰臥位MRI画像の間で対応する位置が同定できる場合には、これらの位置情報を処理装置200が取得し、対象症例変形生成部1420は、その対応の関係を近似するように係数ベクトルを算出するようにしても良い。例えば、処理装置200が、この対応する位置をユーザの入力によって取得する場合には、ユーザが両画像間で合致して欲しいと考える位置が合致させるように変形を推定するようにできる効果がある。
また、数21に示した画像間類似度の評価関数に、さらに対象症例の伏臥位MRI画像と仰臥位MRI画像の間で対応する位置の誤差の評価関数を加えた関数を新たな評価関数としても良い。これによれば、画像間類似度とユーザが入力した対応する位置との両方の情報に基づいた、より精度の高い変形推定が行える効果がある。なお、上記の対応する位置の情報は、必ずしもユーザによる入力によって取得する場合に限らず、伏臥位MRI画像と仰臥位MRI画像の双方の画像からnSIFTなどによる特徴点検出・特徴点対応付けの方法などを用いて、自動的に対応する位置の情報を取得するようにしても良い。これによれば、より効率的に変形推定が実行できる効果がある。
また、本処理ステップにおいて算出する係数ベクトルは、必ずしも評価関数の最大化等によって得る必要はなく、例えば、係数ベクトルは0ベクトルであっても良い。この場合、本処理ステップで算出される変形は、本実施形態における学習フェーズの処理であるS580で算出した平均的な変形となるが、前述の第1実施形態においてφps(x)を恒等関数とする場合との比較によれば、より高精度な変形推定が実行できる。
(ステップS650)MRI画像を変形
ステップS600において、画像変形部1060は、ステップS640で算出した変換に基づいて伏臥位MRI画像を変形した変形MRI画像を生成する処理を実行する。具体的には、画像変形部1060は、数24に示した、画像の変形関数D(I,b)を用い、対象症例の伏臥位MRI画像Ipを係数ベクトルboptに基づいて変形した変形MRI画像 Idを算出する。
(ステップS660)変形画像を表示
ステップS600において、観察画像生成部1080は、ステップS650で生成した変形MRI画像Idおよび対象症例の仰臥位MRI画像Isとを並べた観察画像を生成する。本処理ステップの具体的な処理は第1実施形態におけるステップS390の処理と同様であるため、詳細な説明は省略する。
以上に説明したステップS600からステップS660の処理により、本実施形態の変形推定フェーズの処理が実行される。この処理の結果、対象症例に関する伏臥位MRI画像と仰臥位MRI画像の間の変形推定処理が実行される。そして、仰臥位MRIと対応するように伏臥位MRI画像を変形した変形MRI画像が生成され、仰臥位MRI画像と共に表示されることで、比較容易な形態で入力画像を提示される。
(第3実施形態の変形例1)標準形状(お椀型,円錐などを含む)に正規化する,伏臥位・仰臥位で異なる正規化をする場合も含む
本実施形態におけるステップS540およびステップS610の処理は、基準点である乳頭位置からの体表面および大胸筋面を構成する点群の距離と方位に基づいて正規化座標系への変換φpおよびφsを算出する場合を例として説明した。しかし、変換φpおよびφsの算出方法はこの例に限定されない。例えば、ステップS540およびステップS610における正規化座標系への変換φpおよびφsは,伏臥位および仰臥位に関する標準形状への変換関数とすることができる。ここで標準形状とは、例えば、特定の1症例の乳房形状とすることができる。
具体的には、本実施形態において学習するNsamples個の症例の中の1症例をユーザが任意に選択し、選択した症例の体表形状、大胸筋形状を標準形状とすることができる。また、標準形状は学習するNsamples個の症例の平均形状とするようにしても良い。また、標準形状は必ずしも学習症例の中から選択する場合に限らず、学習症例とは異なる症例の体表形状、大胸筋形状を用いるようにしても良い。また、実在する特定の症例の形状を用いる場合に限らず、例えば、人工的に構築した乳房の模擬形状であっても良い。例えば、球や楕円体を切断した形状や円錐,お椀のような形状等を標準形状としても良い。また、伏臥位と仰臥位とで同一の形状を標準形状とするようにしても良いし、夫々異なる形状を標準形状とするようにしても良い。また、予め複数の標準形状に関する情報を処理装置200が持ち、学習症例の解剖学特徴に基づいて適切な標準形状を選択するようにしても良い。何れの場合も、本実施形態における統計変形モデルは、伏臥位の標準形状を基準とした座標系と仰臥位の標準形状を基準とした座標系との間の変換関数を近似して表現するモデルである。
(第3実施形態の変形例2)学習症例の正規化座標系の間の変形はFFDでも良い。その場合、統計モデルはSDMを用いても良い。
本実施形態では、ステップS580の処理として、ステップS540で算出した変形φps,j(x)を,正規化座標系の空間における変形の場(離散化ベクトル)として展開する場合を例として説明したが、別の形態によってこの展開が行われても良い。例えば,ステップS540では変形φps,j(x)をFFDにより表現し、FFDのパラメータ(制御点が持つ制御量)をベクトル化して離散化ベクトルpx,j、py,j、pz,jを算出するようにしても良い。この場合、統計変形モデルは、非特許文献2に開示のある手法(Statistical Deformation Model法)により構築できる。この方法によれば,ステップS580において、変形φps,j(x)を変形場の離散化ベクトルとして展開する処理を実行する必要がなく、計算処理量とメモリ容量の消費を低減できる効果がある.
[第4実施形態]伏臥位乳房の統計アトラスを構築する例
本実施形態による処理装置800は、多症例の乳房を伏臥位の体位で撮像したMRI画像を学習データとして、各症例間で異なる乳房の形状を正規化した上で、その形状を効率的に表現する統計形状モデルを構築する。なお、本実施形態では伏臥位の体位で撮像した乳房を学習データとして用いる場合を例として説明するが、他の体位における乳房であっても良い。例えば仰臥位や立位等であっても良い。
(機能構成)
図12は、本実施形態による処理システムの構成を示す図である。本実施形態において、第1実施形態と同様の機能を持つ構成要素には図1と同一の番号を付しており、説明は省略する。同図に示すように、本実施形態における処理装置800は、画像取得部1000、解剖学特徴抽出部1020、正規化部1040、スケール算出部1620、統計形状モデル生成部1640を持つ。スケール算出部1620は学習症例の夫々に関するスケールの算出を行う。統計形状モデル生成部1640は学習症例の夫々に関する変形・スケールに基づいて統計形状モデルを生成する。
(処理フロー)
次に、処理装置800が行う全体の動作について説明する。本実施形態では、主メモリ212に格納されている各部の機能を実現するプログラムをCPU211が実行することにより実現される。また以下に説明する処理装置200が行う各処理の結果は、主メモリ212に格納することにより記録される。図13は、本実施形態における処理装置800の処理手順を説明するフローチャートである。以下、この図に沿って説明する。
(ステップS700)学習症例の画像取得
ステップS700において、画像取得部1000は、学習症例の伏臥位MRI画像を取得する処理を実行する。本実施形態では、画像取得部1000は、Nsamples個の学習症例の伏臥位MRI画像を取得する。この処理は、第1実施形態におけるステップS300と同様の処理をNsamples個の学習症例の夫々について適用して実行する。詳細な説明は省略する。
(ステップS705)学習症例の解剖学特徴を抽出
ステップS705において解剖学特徴抽出部1020は、ステップS700で取得した学習症例の伏臥位MRI画像の夫々を処理して、学習症例の解剖学特徴を抽出する処理を実行する。この処理は、第1実施形態におけるステップS310と同様の処理をNsamples個の学習症例の夫々について適用して実行する。詳細な説明は省略する。
(ステップS710)学習症例の正規化座標系への変換を算出
ステップS710において、正規化部1040は、複数の学習症例の夫々に関して、被検体の形状を基準形状に変換する正規化変換を導出する。具体的には、正規化部1040は、MRI画像座標系から正規化座標系への座標変換関数を算出する処理を実行する。この処理は、第1実施形態におけるステップS320と同様の処理を、Nsamples個の学習症例の夫々について適用して実行する。詳細な説明は省略する。
本実施形態では算出した正規化座標系への変換を関数φp,j(x)と表記する。ここで関数φp,j(x)は各症例毎に算出され、各症例における伏臥位MRI画像座標系から伏臥位正規化座標系への変換関数である。伏臥位正規化座標系から伏臥位MRI画像座標系への変換についても第1実施形態と同様に算出する。これを、φ-1 p,j(x)と表記する。ここで、jは学習症例の症例番号のインデックスであり、本実施形態では1≦j≦Nsamplesである。なお、本実施形態では学習症例の夫々に関して、伏臥位MRI画像における乳頭位置が、伏臥位MRI画像座標系の原点となるように、画像の並進、もしくは座標系の並進が予め施されているものとする。
(ステップS720)学習症例のスケールの算出
ステップS720において、スケール算出部1620は、学習症例の夫々に関して、症例のスケールを算出する処理を実行する。この処理は、第3実施形態におけるステップS550の処理と同様にして実行する。詳細な説明は省略する。算出したスケールをvj (1≦j≦Nsamples)と表記する。
(ステップS740)変換関数のスケーリング
ステップS740において、統計形状モデル生成部1640は、Nsamples個の症例に関する変換関数φ-1 p,j(x)と、スケールvjとに基づいて、変換関数φ-1 p,j(x)をスケーリング処理した変換関数φ’-1 p,j(x)を算出する処理を実行する。具体的には数26のようにして変換関数φ’-1 p,j(x)を算出する。
(数26)
φ’-1 p,j (x’)= φ-1 p,j (x) / vj
ただし、ここでx=( x , y , z )Tとした場合に、x’=( x×vj , y , z×vj )Tである。つまり、φ’-1 p,j (x)はφ-1 p,j(x)を、定義域についてX座標およびZ座標の方向に学習症例のスケール値vjでスケーリングし、さらに関数の値についてスケール値vjでスケーリングした関数である。統計形状モデル生成部1640は、以上の処理をNsamples個の全症例の変換関数φ-1 p,j(x) (1≦j≦Nsamples)について実行し、各症例のスケーリング処理した変換関数φ’-1 p,j(x) (1≦j≦Nsamples)を算出する。
(ステップS760)統計形状モデル生成
ステップS760において、統計形状モデル生成部1640は、ステップS740で算出した変換関数φ’-1 p,j(x) (1≦j≦Nsamples)を統計処理することで、統計形状モデルを生成する処理を実行する。具体的な処理手順について説明する。
まず、統計形状モデル生成部1640は、変換関数φ’-1 p,j(x)を離散化する処理を実行する。この処理は第3実施形態のステップS5820においてφ’ps,j(x)を対象とした処理と同様の処理を、本実施形態における変換関数φ’-1 p,j(x)を対象として実行する。詳細な説明は省略する。これによって得た離散化ベクトルを、qx,j、qy,j、qz,j、とする。統計形状モデル生成部1640は、この処理をNsamples個の全ての変換関数φ’-1 p,j(x)について適用する。ここで、算出した離散化ベクトルqx,j、qy,j、qz,jの夫々は正規化座標系における変換関数をサンプリングしたベクトルである。正規化座標系においては、学習症例の夫々が解剖学的に略一致した位置関係を持つことから、これら離散ベクトルの同一の次元の値は、学習症例の夫々において解剖学的に略同一の位置における変換関数の値を持つ。
次に、統計形状モデル生成部1640は、算出した離散化ベクトルqx,j、qy,j、qz,j(1≦j≦Nsamples)について主成分分析による統計処理によって統計形状モデルを生成する処理を実行する。この処理は第3実施形態のステップS5840の処理と同様にして実行する。詳細な説明は省略する。以上の処理によって得た平均ベクトルを、e’ave,x、e’ave,y、e’ave,zと表記する。また得られた固有ベクトルを e’x,k、e’y,k、e’z,k と表記する。ここでkは主成分分析によって得た固有ベクトルのインデックスであり、固有ベクトルをN’mode個取得する場合には、1≦k≦N’modeである。
以上に説明したステップS700からステップS760の処理により、本実施形態における処理装置800は、学習症例に基づいて統計形状モデルを生成する。この統計形状モデルは、第3実施形態のステップS5840における説明として数20に示したように、学習症例の夫々の変換関数φ’-1 p,j(x)を少ない数の基底で効率的に近似するモデルである。このモデルは学習症例に共通する統計的な特性を考慮したものであり、例えば学習症例には含まれない未知の症例が与えられた場合でも、その症例における変換関数を精度よく近似できることが期待できる。ここで変換関数は正規化座標系を定義域とした関数であり、正規化座標系における体表面および大胸筋面の位置は既知であるから、以下のようにして統計形状モデルを利用できる。
すなわち、未知症例の伏臥位MRI画像から抽出される体表面や大胸筋面の形状と、統計形状モデルが表す体表面や大胸筋面の形状とが略合致するように統計形状モデルの係数ベクトルを求めることができる。ここで、未知症例の伏臥位MRI画像から抽出する体表面や大胸筋面の形状は一部に情報の欠損がある場合や、雑音を含むような場合であっても良い。つまり、未知症例の伏臥位MRI画像から体表面や大胸筋面の形状に関する限られた観測情報から係数ベクトルを求めることができる。係数ベクトルが算出されれば、未知症例に対する変換関数が推定できるため、限られた観測情報を補う情報を生成することができる。つまり、統計形状モデルは乳房領域のセグメンテーションに利用できる。
[第5実施形態]
本実施形態では、第3実施形態で説明した統計変形モデルおよび、第4実施形態で説明した統計形状モデルの両方を構築し、これらのモデルを用いることで対象症例に対する正規化処理を簡便かつロバストに実行する処理装置900を例示する。
(機能構成)
図14は、本実施形態による処理システムの構成を示す図である。本実施形態において、第1乃至第4実施形態と同様の機能を持つ構成要素には同一の番号を付しており、説明は省略する。同図に示すように、本実施形態における処理装置900は、画像取得部1000、解剖学特徴抽出部1020、正規化部1040、スケール算出部1230、統計形状モデル生成部1640、統計変形モデル生成部1240を持つ。また処理装置900は、対象症例形状抽出部1800、対象症例正規化部1820、対象症例スケール算出部1830、対象症例変形生成部1420、画像変形部1060、観察画像生成部1080を持つ。
対象症例形状抽出部1800は、画像取得部1000が取得した対象症例のMRI画像から、対象症例の体表および大胸筋の形状の抽出を行う。対象症例スケール算出部1830は、画像取得部が取得した対象症例のMRI画像、および、対象症例形状抽出部1800が抽出した対象症例の体表および大胸筋の形状に基づいて、対象症例のスケールを算出する。対象症例正規化部1820は、統計形状モデル生成部1640が生成した統計形状モデルと、対象症例形状抽出部1800が抽出した対象症例の体表および大胸筋の形状、および、対象症例スケール算出部1830が算出した対象症例のスケールに基づき、対処症例の正規化座標系への変換を算出する。
(処理フロー)
次に、処理装置900が行う全体の動作について説明する。本実施形態では、主メモリ212に格納されている各部の機能を実現するプログラムをCPU211が実行することにより実現される。また以下に説明する処理装置900が行う各処理の結果は、主メモリ212に格納することにより記録される。本実施形態における処理装置900の処理は、学習フェーズの処理と変形推定フェーズの処理から成り、まず学習フェーズの処理を実行し、その後に変形推定フェーズの処理を実行する。学習フェーズの処理では、多症例の伏臥位および仰臥位のMRI画像間の変形を学習し、統計変形モデルを生成する処理を実行する。また学習フェーズの処理では、多症例の伏臥位および仰臥位の形状を学習し統計形状モデルを生成する処理を実行する。変形推定フェーズの処理では学習フェーズで生成した統計変形モデルおよび統計形状モデルを用いて対象症例の伏臥位、仰臥位間の変形位置合わせを実行する。
以下に記載する本実施形態の説明では、処理装置900が前記学習フェーズの処理と変形推定フェーズの処理の両方を実行する場合を例として説明するが、学習フェーズの処理と、変形推定フェーズの処理を異なる処理装置で実行するようにしても良い。また、本実施形態による処理装置900は、学習フェーズの処理と変形推定フェーズの処理の両方の処理を実行する場合に限らず、例えば、学習フェーズの処理だけを実行するようにしても良い。また、学習フェーズの処理の結果として得られる統計変形モデルおよび統計形状モデルの提供も、本実施形態に含まれる。
(学習フェーズの処理)
図15は、本実施形態による処理装置900の学習フェーズの処理手順を説明するフローチャートである。以下、この図に沿って説明する。
(S6000)〜(S6050)
ステップS6000からステップS6050において処理装置900は、第3実施形態における処理装置200が実行するステップS500からステップS580と同様の処理を実行する。ステップS6040で算出した学習症例のスケール値をvjと表記する。ここで、jは学習症例の症例番号のインデックスであり、本実施形態では1≦j≦Nsamplesである。また、ステップS6050で生成した統計変形モデルの平均ベクトルをe deform_ave,x、e deform_ave,y、e deform_ave,z、固有ベクトルをe deform_x,k、e deform_y,k、e deform_z,kと表記する。kは主成分分析によって得た複数の固有ベクトルのインデックスである。本実施形態では固有ベクトルをN deform_mode個取得することとする。すなわち、1≦k≦N deform_modeである。詳細な説明は省略する。
(S6070)伏臥位統計形状モデル生成
ステップS6070において、統計形状モデル生成部1640は、学習症例の伏臥位の形状に対して第4実施形態におけるステップS740およびステップS760と同様の処理を実行し、学習症例の伏臥位の形状に関する統計形状モデルを生成する。この処理において、学習症例のスケーリングに用いるスケール値は、ステップS6040で算出したスケール値vjである。以上の処理により伏臥位の統計形状モデルが生成される。ここで生成した伏臥位の統計形状モデルの平均ベクトルをe p_shape,ave,x、e p_shape,ave,y、e p_shape,ave,z、固有ベクトルをe p_shape,x,k、e p_shape,y,k、e p_shape,z,kと表記する。ここでkは主成分分析によって得た複数の固有ベクトルのインデックスである。本実施形態では固有ベクトルをN p_shape_mode個取得することとする。すなわち、1≦k≦N p_shape_modeである。本処理ステップの詳細な説明は、第4実施形態のステップS740およびステップS760の処理の説明と重複するため省略する。
(S6080)仰臥位統計形状モデル生成
ステップS6080において、統計形状モデル生成部1640は、学習症例の仰臥位の形状に対して第4実施形態におけるステップS740およびステップS760と同様の処理を実行し、学習症例の仰臥位の形状に関する統計形状モデルを生成する。この処理において、学習症例のスケーリングに用いるスケール値は、ステップS6040で算出したスケール値vjである。
以上の処理により伏臥位の統計形状モデルが生成される。ここで生成した伏臥位の統計形状モデルの平均ベクトルをe s_shape,ave,x、e s_shape,ave,y、e s_shape,ave,z、固有ベクトルをe s_shape,x,k、e s_shape,y,k、e s_shape,z,kと表記する。ここでkは主成分分析によって得た固有ベクトルのインデックスである。本実施形態では固有ベクトルをN s_shape_mode個取得することとする。すなわち、1≦k≦N s_shape_modeである。本処理ステップの詳細な説明は、第4実施形態のステップS740およびステップS760の処理の説明と重複するため省略する。以上に説明した処理により、処理装置900は、学習症例に関する統計変形モデルと伏臥位および仰臥位の統計形状モデルを生成する。
(変形推定フェーズの処理)
次に図16のフローチャートを用いて処理装置900が実行する変形推定フェーズの処理手順について説明する。
(S6100)
ステップS6100において処理装置900は、学習フェーズの処理で生成した伏臥位および仰臥位の統計形状モデルを、処理装置900の主メモリ212に読みだす処理を実行する。
(S6110)から(S6120):統計変形モデル読み込み、対象症例の画像取得
ステップS6110からステップS6120において処理装置900は、第3実施形態における処理装置200が変形推定フェーズの処理として実行するステップS600からステップS602と同様の処理を実行する。説明は省略する。
(S6130)対象症例の形状を抽出
ステップS6130において、対象症例形状抽出部1800は、ステップS6120で取得した対象症例の伏臥位MRI画像および仰臥位MRI画像を処理することにより、対象症例の伏臥位および仰臥位における体表面および大胸筋面の形状を抽出する処理を実行する。具体的には、対象症例形状抽出部1800は、夫々のMRI画像から、体表面と大胸筋面の位置を表す複数の点群を抽出する。この処理は、第3実施形態の処理装置200が実行するステップS604の処理の一部と同様の処理である。ただし、上記の形状の抽出はMRI画像に描出される被検体の乳房の全領域に渡る領域を対象とする必要はなく、乳房領域の一部やその周辺の領域における形状を検出すれば良い。また、抽出する形状は、第3実施形態のステップS604の処理として記載した解剖学特徴の抽出の結果のように密な点群である必要はなく、比較的に疎な点の集合であって良い。本実施形態では、抽出した伏臥位の形状をsp,surface,i (1≦i≦Np,surface)、仰臥位の形状をss,surface,i (1≦i≦Ns,surface)とする。ここで、N p,surfaceは伏臥位の形状を表す点の数であり、同様にN s,surfaceは伏臥位の形状を表す点の数である。
また、ステップS6130において、対象症例形状抽出部1800は、前記処理を実行するとともに、対象症例の伏臥位および仰臥位における乳頭位置および大胸筋面上の基準位置を取得する処理を実行する。この処理は例えば、処理装置900がモニタ160に伏臥位MRI画像および仰臥位MRI画像を表示し、表示された画面上をユーザがマウス170やキーボード180によって位置を指定し、その結果を対象症例形状抽出部1800が取得することにより実行される。
(S6140)対象症例のスケール算出
ステップS6140において、対象症例スケール算出部1830は、本実施形態の学習フェーズにおけるステップS6040で学習症例に対して実行する処理と同様の処理を、対象症例に対して実行する。これにより、対象症例のスケールvtargetを算出する。
(S6150)統計形状モデルを用いた正規化
ステップS6150において、対象症例正規化部1820は、対象症例の伏臥位の形状sp,surface,i、仰臥位の形状ss,surface,i、スケールvtargetに基づき、対象症例に対して統計形状モデルのパラメータを最適化する処理を実行し、正規化座標系への変換を算出する。
具体的には、対象症例正規化部1820は、以下の処理を実行する。まず、対象症例正規化部1820は、ステップS6070で取得した伏臥位の統計形状モデルに対するパラメータの最適化を行う。ここでパラメータとは、統計形状モデルの複数の固有ベクトルに対する重み係数であり、固有ベクトルの数と同数の次元をもつベクトルである。このベクトルの値の最適化は、以下の基準により実行される。すなわち、伏臥位MRI画像座標系における伏臥位の形状sp,surface,iを正規化座標系における伏臥位の形状sp,surface,i’に変換した際の、sp,surface,i’が示す位置と、正規化座標系における基準形状との差異を最小化する基準である。つまり、対象症例正規化部1820は、伏臥位の形状sp,surface,iが基準形状へと適切にマッピングされるように統計形状モデルのパラメータを最適化する。
ところで、本実施形態における伏臥位の統計形状モデルは伏臥位の正規化座標系から伏臥位MRI画像座標系への変換を表現するモデルである。そのため、上記の最適化の処理は、具体的には以下の手順で実行する。まず、対象症例正規化部1820は、伏臥位の正規化座標系における体表面と大胸筋面の形状を点群や多角形などの任意の形状表現で表す。なお、本実施形態では伏臥位の正規化座標系における体表面と大胸筋面の形状は共に平面である。次に、対象症例正規化部1820は、伏臥位の統計形状モデルのパラメータに任意の初期値を設定し、このパラメータで表現される伏臥位の統計形状モデルを用いて、体表面と大胸筋面の形状を、伏臥位のMRI画像座標系に変換する。そして、対象症例正規化部1820は、変換した体表面と大胸筋面の形状と伏臥位のMRI画像座標系の形状sp,surface,iとの差異を評価する。そして、伏臥位の統計形状モデルのパラメータを様々に変更して上記の差異の評価が最小となるパラメータを探索する。すなわち、対象症例正規化部1820は、伏臥位の統計形状モデルのパラメータを上記の差異の評価に基づいて、繰り返し処理によって最適化する。
なお、上記の処理において対象症例の伏臥位MRI画像座標系における形状sp,surface,iの座標値と、統計形状モデルとの間のスケーリングには、ステップS6140で算出したスケールvtargetが考慮される。仰臥位の統計形状モデルに関しても、同様に、ステップS6080で取得した仰臥位の統計形状モデルのパラメータに関して、仰臥位の形状ss,surface,iに基づいて最適化する。
以上の処理により、伏臥位および仰臥位に関する統計形状モデルのパラメータが対象症例の形状に対して最適化される。これにより、対象症例の伏臥位に関する正規化座標系への変換φp,target(x)および仰臥位に関する正規化座標系への変換φs,target(x)が算出される。
(S6160)から(S6180)統計変形モデルの係数の最適化、MRI画像変形、表示
ステップS6160からステップS6180において、処理装置900は、第3実施形態における処理装置200が実行するステップS640からステップS660と同様の処理を実行する。詳細な説明は省略する。
以上に説明したステップS6100からステップS6180の処理により、本実施形態の変形推定フェーズの処理が実行される。この処理の結果、対象症例に関する伏臥位MRIと仰臥位MRI画像の間の変形推定処理が実行される。そして、仰臥位MRIと対応するように伏臥位MRI画像を変形した変形MRI画像が生成され、仰臥位MRI画像と共に表示されることで、容易に比較できる形態で入力画像を提示できる。
本実施形態の処理装置によれば、第3実施形態の処理装置と比較して、対象症例の体表面形状や大胸筋面形状を表す密な点群に関する情報を含む解剖学特徴を抽出する必要がなく、比較的に疎な体表面形状、大胸筋面形状を抽出するだけで良い。なぜならば、体表面形状と大胸筋面形状を表す点群が空間的に疎な場合であっても、統計形状モデルのパラメータを推定することは可能であり、そのパラメータに基づいてMRI画像座標系から各正規化座標系への変換が算出できるからである。そのため、対象症例のMRI画像の画質等の影響により体表面形状および大胸筋面形状を表す点群を密に抽出することが困難な場合であっても、両画像間の変形位置合わせを実行することが可能となる。
(第5実施形態の変形例1)統計形状モデルの使用・不使用・切り替え
本実施形態における変形推定フェーズでは、伏臥位・仰臥位の両方の統計形状モデルを用いる場合を例として説明したが、何れか一方であっても良い。例えば、ステップS6130の処理のうち、伏臥位または仰臥位の何れか一方の体位に関する処理は、第3実施形態におけるステップS604と同様の処理を実行するようにしても良い。この場合、後段の処理であるステップS6150の処理は、上記の伏臥位または仰臥位の何れか一方の体位に関する正規化座標系への変換の算出処理を、第3実施形態におけるステップS610と同様の処理としてもよい。これによれば、変形推定フェーズにおいて、対象症例の伏臥位または仰臥位のいずれか一方の体位における密な体表面形状および大胸筋面形状が得られる場合にはMRI画像座標系から正規化座標系への変換を、より高い精度で算出できる効果がある。
また、本実施形態や上記の変形例いずれかの方法を、対象症例の解剖学特徴の抽出結果に基づいて切り替えて実行するようにしてもよい。例えば、第3実施形態に記載の処理方法を標準的な処理として実行しつつ、対象症例のMRI画像から密な体表面形状および大胸筋面形状の抽出が困難または困難と予測される場合に、本実施形態および上記変形例の処理に切り替えて実行するようにしてもよい。これによれば、対象症例のMRI画像の画質等の影響を考慮して適切な処理方法を切り替えて実行することができるため、対象症例の変形位置合わせを、より頑健に行える効果がある。
また、本実施形態は、対象症例の伏臥位および仰臥位の両方の体位におけるMRI画像が取得される場合に限定されない。例えば、対象症例の伏臥位または仰臥位のいずれか一方の体位においては、MRI画像が取得されず、比較的に疎な体表面の形状だけが取得されるような場合であってもよい。一例として、対象症例の仰臥位におけるMRI画像が取得されず、対象症例の仰臥位における体表面の形状を、位置計測が可能なスタイラス等を用いて計測される場合について説明する。この場合、変形推定フェーズの処理のステップS6120では対象症例の伏臥位に関するMRI画像を取得するまたステップS6130の処理では、伏臥位に関しては第3実施形態におけるステップS604と同様の処理を実行するとともに、仰臥位に関しては前記スタイラス等を用いて計測した対象症例の疎な体表面の形状を取得する。ステップS6150では、伏臥位に関する正規化座標系への変換の算出処理を、第3実施形態におけるステップS610と同様の処理を実行し、仰臥位に関しては本実施形態のS6150の処理を前記疎な体表面の形状に基づいて実行する。ステップS6160以降の処理は本実施形態に記載の通りの処理を実行する。
以上の方法により、対象症例の仰臥位における体表面の形状および、統計形状モデルに基づいて推定される仰臥位における大胸筋面の形状に合致するように対象症例の伏臥位のMRI画像を変形して表示できる。上記の方法によれば、例えば対象症例の手術等を効果的に支援することができる。具体的には、手術に先立って撮影した対象症例の伏臥位におけるMRI画像を、仰臥位の体位において実施される手術時に計測した体表面の形状に基づいて変形して表示することができる。これにより、伏臥位におけるMRI画像においてユーザが注目している部位等が、手術時の仰臥位において、どこに位置するかを手術中にユーザに提示できる。
(第5実施形態の変形例2)伏臥位・仰臥位の統合モデル
本実施形態では、統計形状モデル生成部1640が伏臥位および仰臥位の夫々に関する統計形状モデルを個別に構築する場合を例として説明したが、このように個別に構築することに限定されない。例えば、統計形状モデル生成部1640は、伏臥位の形状に関する情報と仰臥位の形状とを統合した統計形状モデルを構築するようにしてもよい。具体例としては、学習フェーズの処理では、ステップS6070およびステップS6080の処理に代えて、以下の処理を実行する。すなわち、統計形状モデル生成部1640は、第4実施形態に記載した処理装置800がステップS760で実行する離散化ベクトルqx,j、qy,j、qz,j、の算出処理を、学習症例の伏臥位および仰臥位の夫々について実行する。ここでは、算出した伏臥位に関する離散化ベクトルをq p_x,j、q p_y,j、q p_z,j、仰臥位に関する離散化ベクトルをq s_x,j、q s_y,j、q s_z,j、と表記する。
そして、統計形状モデル生成部1640は、これら6個のベクトルを結合したベクトルを学習症例の夫々に関して算出し、そのベクトルの群に対して主成分分析を実行することで、平均ベクトルをおよび固有ベクトルを算出する。これらの情報を本実施形態では伏・仰臥位統計形状モデルと称する。伏・仰臥位統計形状モデルは、学習症例に関する伏臥位のMRI画像座標系から伏臥位の正規化座標系への変換と、仰臥位のMRI座標系から仰臥位の正規化座標系への変換の両方の変換を表すモデルとなる。また、このモデルはこの2つの変換の間の統計的な特性を記述したモデルとなる。変形推定フェーズでは、ステップS6150の処理に代えて以下の処理を実行する。すなわち、対象症例正規化部1820は、伏臥位および仰臥位のMRI画像座標系から夫々の正規化座標系への変換を、上記のモデルのパラメータの最適化処理によって算出する。これにより、対象症例の伏臥位および仰臥位の夫々の正規化座標系への変換を、両変換の間の統計的な特性を考慮して算出できるため、正規化座標系への変換を、より精度良く実行できる効果がある。
また、対象症例における伏臥位または仰臥位の何れか一方の正規化座標系への変換を別途(例えば第4実施形態に記載の方法等により)算出できる場合に、他方の変換を上記の伏・仰臥位統計形状モデルを用いて推定することができる。例えば、対象症例の伏臥位または仰臥位の何れか一方の形状の抽出が実行できないような場合であっても、上記の統計形状モデルを用いることで、他方の正規化座標系への変換を推定できる。これにより、変形推定フェーズの処理のステップS6160以降の処理を実行できる。この方法によれば、対象症例の伏臥位または仰臥位の何れか一方の体位におけるMRI画像が取得されないような場合においても、MRI画像が取得できる体位におけるMRI画像を、他方の体位へ統計的に尤もらしく変形させることができる。
[第6実施形態]
第5実施形態の説明では、統計変形モデルと、伏臥位および仰臥位に関する統計形状モデルとを生成し、各モデルを用いて対象症例の伏臥位MRI画像と仰臥位MRI画像の間の変形を推定する場合を例に説明した。第6実施形態では、伏臥位および仰臥位のMRI画像座標系と各正規化座標系への変換と、各正規化座標系の間の変形とを統合したモデルを生成し、生成したモデルを用いて対象症例の伏臥位および伏臥位のMRI画像間の変形を推定する方法について説明する。本実施形態では、生成するモデルを統計モデルと称する。
(機能構成)
図17は、本実施形態における処理装置950の機能構成を表す図である。同図において、前記第1乃至第5実施形態と同様の機能を持つ構成要素には、各実施形態で付した番号と同一の番号を付しており、説明は省略する。統計モデル生成部1840は、正規化部1040および学習症例変形生成部1220およびスケール算出部1230が実行する処理の結果に基づき、学習症例の統計モデルを生成する。対象症例変形生成部1850は、統計モデル生成部1840が生成する統計モデルおよび、解剖学特徴抽出部1020が抽出する対象症例の解剖学特徴およびスケール算出部1230が算出するスケールに基づいて、対象症例の変形を生成する。
(処理フロー)
次に、処理装置950が行う全体の動作について説明する。本実施形態では、主メモリ212に格納されている各部の機能を実現するプログラムをCPU211が実行することにより実現される。また、以下に説明する処理装置200が行う各処理の結果は、主メモリ212に格納することにより記録される。本実施形態における処理装置950の処理は、第5実施形態と同様に、学習フェーズの処理と変形推定フェーズの処理から成り、まず学習フェーズの処理を実行し、その後に変形推定フェーズの処理を実行する。学習フェーズの処理では、多症例の伏臥位および仰臥位のMRI画像間の変形を学習し、統計モデルを生成する処理を実行する。変形推定フェーズの処理では学習フェーズで生成した統計モデルを用いて対象症例の伏臥位、仰臥位間の変形位置合わせを実行する。
(学習フェーズの処理)
図18は、本実施形態における処理装置950が行う、学習フェーズの処理の手順を説明するフローチャートである。以下、このフローチャートに示す処理手順に従って、本実施形態の学習フェーズの処理について詳しく説明する。
(S7000)から(S7040)
ステップS7000からステップS7040において処理装置950は、第5実施形態における処理装置900が実行するステップS6000からステップS6040と同様の処理を実行する。詳細な説明は省略する。
(S7050)
ステップS7050において統計モデル生成部1840は、統計モデルを生成する処理を実行する。この処理について詳しく説明する。まず統計モデル生成部1840は、学習症例の伏臥位および仰臥位の夫々に対して、第4実施形態におけるステップS740およびS760の一部と同様の処理を実行し、MRI画像座標系から正規化座標系への変換に関する離散ベクトルを算出する。本実施形態では、伏臥位のMRI画像座標系から伏臥位の正規化座標系への変換に関する離散ベクトルをqp_x,j、qp_y,j、qp_z,j、同様に仰臥位のMRI画像座標系から仰臥位の正規化座標系への変換に関する離散ベクトルをqs_x,j、qs_y,j、qs_z,j、と記す。
次に、統計モデル生成部1840は、学習症例の伏臥位および仰臥位の夫々の正規化座標系の間の変形に関する離散ベクトルpx,j、py,j、pz,jを算出する。この処理は、第3実施形態におけるステップS5800およびステップS5820と同様の処理により実行する。詳細な説明は省略する。
次に、統計モデル生成部1840は、各症例毎に上記の方法で算出したベクトルを結合したベクトルを生成する。そして、統計モデル生成部1840は、各症例に関する結合ベクトルについて主成分分析を行い、平均ベクトルと複数の固有ベクトルとを算出する。以後、この処理によって算出した平均ベクトルと複数の固有ベクトルを統計モデルと称する。この統計モデルは、後述する変形推定フェーズにおいて利用される。統計モデルが持つ平均ベクトルと複数の固有ベクトルの加重和により、MRI画像座標系から正規化座標系への変換と、各正規化座標系の間の変形の統計的な特性が記述される。
以上に説明したステップS7000からステップS7040の処理により、本実施形態の学習フェーズの処理が実行される。この処理の結果、統計モデルが生成される。
(変形推定フェーズの処理)
図19は、本実施形態による処理装置950が行う、変形推定フェーズの処理の手順を説明するフローチャートである。以下、このフローチャートに示す処理手順に従って、本実施形態の変形推定フェーズの処理について詳しく説明する。
(S7100)
ステップS6700において処理装置950は、学習フェーズの処理で生成した統計モデルを、主メモリ212に読みだす処理を実行する。
(S7120)から(S7140)
ステップS7120からステップS7140において処理装置950は、第5実施形態における処理装置900が実行するステップS6120からステップS6140と同様の処理を実行する。詳細な説明は省略する。
(S7160)
ステップS7160において、対象症例変形生成部1850は、ステップS7100で取得した統計モデル、ステップS7130で取得した対象症例の体表面および大胸筋面の形状、およびステップS7140で算出したスケールに基づいて、対象症例の変形を生成する処理を実行する。この処理について詳しく説明する。
本処理ステップにおいて、対象症例変形生成部1850は、ステップS7100で取得した統計モデルに関するパラメータを、対象症例に対して最適化することによって変形を推定する。ここでパラメータとは、生成された統計モデルの固有ベクトルに対する重み係数を表すベクトルである。対象症例変形生成部1850は、このパラメータに基づいて統計モデルの平均ベクトルと固有ベクトルとの加重線形和を演算することにより、伏臥位および仰臥位のMRI画像座標系から正規化座標系への変換と、各正規化座標系の間の変形を生成する。
統計モデルに関するパラメータの最適化は、以下の評価基準を組み合わせた評価に基づいて実行するようにできる。すなわち、第5実施形態のステップS6150の処理として説明した統計形状モデルの最適化に関する基準(正規化評価の基準と称する)、および第3実施形態のステップS640の処理として説明した評価関数G(変形評価の基準と称する)を組み合わせて実行できる。具体的には、対象症例変形生成部1850は、正規化評価の基準と変形評価の基準の両方の基準により算出される評価値を最小化、または構成によっては最大化、するように統計モデルに関するパラメータを最適化する。
以上の処理により、統計モデルのパラメータが最適化され、対象症例の伏臥位および仰臥位のMRI画像座標系から各正規化座標系への変換、および各正規化座標系の間の変形が推定される。
(S7170)から(S7180)
ステップS7170からステップS7180において処理装置950は、第5実施形態における処理装置900が実行するステップS6170からステップS6180と同様の処理を実行する。詳細な説明は省略する。
以上に説明したステップS7100からステップS7180の処理により、本実施形態の変形推定フェーズの処理が実行される。この処理の結果、対象症例に関する伏臥位MRI画像と仰臥位MRI画像の間の変形推定処理が実行される。そして、仰臥位MRIと対応するように伏臥位MRI画像を変形した変形MRI画像が生成され、仰臥位MRI画像と共に表示されることで、比較容易な形態で入力画像が提示される。
本実施形態の処理装置によれば、学習フェーズにおいては、学習症例に関する伏臥位および仰臥位のMRI画像座標系から正規化座標系への変換と、各正規化座標系の間の変形との両方の統計的な特性を持つ統計モデルが生成される。ここで、MRI画像座標系から正規化座標系への変換は、主に各症例の形状の違いを吸収する変換であり、この変換自体には各症例の形状に関する情報が含まれている。本実施形態では、このMRI画像座標系から正規化座標系への変換に関する情報と、各正規化座標系の間の変形とを対にした学習データを主成分分析することによって統計モデルを生成する。これにより、各症例の形状と変形との間の統計的な特性を反映したモデルが生成される。したがって、変形推定フェーズにおける前記統計モデルを用いた変形推定により、対象症例の形状と変形の両方を同時に考慮し、統計的に尤もらしい変形を推定できる。これにより、第5実施形態の処理装置と比べて、高い精度で変形推定を行える効果がある。
(第6実施形態の変形例1)階層モデルにする
本実施形態では、ステップS7050の処理では、統計モデル生成部1840は、学習症例の夫々についてMRI画像座標系から各正規化座標系への変換と、正規化座標系の間の変形の夫々に関する離散ベクトルを算出し、算出した離散ベクトルを結合したベクトルに基づいて統計モデルを生成する場合を例として説明した。しかし、変形例として、例えば、第5実施形態で説明した処理装置900と同様に、統計形状モデルと統計変形モデルを夫々算出し、それに加えて両モデルの上位モデルを構築するようにしても良い。
具体的には、まず、第5実施形態と同様にして、統計モデル生成部1840は、統計形状モデルと統計変形モデルとを生成する。そして、統計モデル生成部1840は、伏臥位および仰臥位のMRI画像座標系から各正規化座標系への変換を、統計形状モデルを用いて表現した際のパラメータを学習症例の夫々について算出する。算出したパラメータをベクトルbshape,j(1≦j≦Nsamples)とする。また、統計モデル生成部1840は、伏臥位および仰臥位の正規化座標系の間の変形を、統計変形モデルを用いて表現した際のパラメータを学習症例の夫々について生成する。算出したパラメータをベクトルbdeform,j(1≦j≦Nsamples)とする。そして、統計モデル生成部1840は、学習症例毎にbshape,jとbdeform,jを結合したベクトルを生成し、結合したベクトルに対して主成分分析を行う。これにより、算出される平均ベクトルおよび複数の固有ベクトルを以て上位モデルとすることができる。この場合、変形推定フェーズでは、上位モデルに関するパラメータを推定することにより、対象症例の変形を推定するようにできる。
以上の方法によれば、統計形状モデルと統計変形モデルとで異なる数の固有ベクトルを利用するようにできる。そのため、第6実施形態に記載した処理装置950が持つ効果に加え、より柔軟に統計モデルを生成するようにできる効果がある。
[第7実施形態]
第1実施形態から第6実施形態では、人体の乳房を処理の対象とする場合を、一例として説明した。しかし、処理の対象は人体の乳房に限定されない。例えば、処理の対象は、人体以外の動物の乳房であっても良い。また、処理の対象は、他の臓器であっても良い。例えば心臓を対象とする場合には、心臓の外壁と内壁によって囲まれた領域を心臓領域とし、心臓の外壁形状から尖点を基準点として抽出することで上記の実施形態を適用できる。この場合、該実施形態による処理装置は、例えば診断対象の心臓の形状と正常な心臓の形状とを位置合わせすることで形状間の比較を行い、形状に現れる心疾患に関する解析情報を生成することができる。
また、上記の実施形態による処理装置は、拍動する心臓の時系列撮像画像の間の位置合わせを行い、心臓形状を時系列にトラッキングすることで、心臓の形状変動に現れる心疾患に関する解析情報を生成することができる。また、該実施形態による処理装置は、同一症例の過去の撮像画像と現在の撮像画像との間、もしくは時期の異なる過去の複数の画像との間の位置合わせに適用し、心疾患等の進行度などを示す解析情報を生成することもできる。本発明の実施はこれ以外にも、例えば肝臓や肺などの他の臓器に適用することができる。
また、上記の実施形態は必ずしも人体の臓器を対象とした医療目的に用いられる場合に限らず、例えば工業用部品等の形状解析や精度分析などに適用することもできる。一例としてはプレス加工による部品成型を行う際に、成形された部品の形状と、金型の形状との間の形状比較に本発明を実施することもできる。これによれば、成形された部品の形状のばらつきが比較的に大きい場合であっても、形状間の比較を頑健に行える効果が期待できる。
[その他の実施形態]
また、本発明は、以下の処理を実行することによっても実現される。即ち、上述した実施形態の機能を実現するソフトウェア(プログラム)を、ネットワーク又は各種記憶媒体を介してシステム或いは装置に供給し、そのシステム或いは装置のコンピュータ(またはCPUやMPU等)がプログラムを読み出して実行する処理である。
110 MRI画像撮像装置; 120 データサーバ; 100、200、800、900、950 処理装置

Claims (19)

  1. 対象物の三次元のボリュームデータである三次元画像データを取得する画像取得手段と、
    前記対象物における少なくとも第一の解剖学特徴と第二の解剖学特徴を取得する特徴取得手段と、
    前記三次元画像データから、前記第一の解剖学特徴における第一の面と前記第二の解剖学特徴における第二の面との距離が所定の条件を満たす位置の集合で構成される断面に対応する断面画像を生成する生成手段と、
    を有することを特徴とする処理装置。
  2. 前記第一の解剖学特徴と第二の解剖学特徴それぞれ、前記対象物の体表面と前記対象物の大胸筋面であり
    前記生成手段は、前記三次元画像データから、前記対象物の体表面からの距離と前記対象物の大胸筋面からの距離が前記所定の条件を満たす曲面形状上の前記断面に対応する前記断面画像を生成する
    ことを特徴とする請求項1に記載の処理装置。
  3. 前記第一の解剖学特徴と第二の解剖学特徴それぞれ、前記対象物の体表面と前記対象物の大胸筋面であり
    前記生成手段は、前記三次元画像データから、前記対象物の体表面からの距離と前記対象物の大胸筋面からの距離の割合が一定である曲面形状上の断面に対応する前記断面画像を生成する
    ことを特徴とする請求項1に記載の処理装置。
  4. 前記第一の解剖学特徴と第二の解剖学特徴それぞれ、前記対象物の体表面と前記対象物の大胸筋面であり
    前記生成手段は、
    前記対象物の体表面を含む面と前記対象物の大胸筋面を含む面とが対向する矩形形状に前記三次元画像データを変形する座標変換に基づいて、前記三次元画像データから、前記対象物の体表面を含む面または前記対象物の大胸筋面を含む面と平行する断面に対応する前記断面画像を生成する
    ことを特徴とする請求項1に記載の処理装置。
  5. 前記生成手段は、
    前記三次元画像データ中の前記対象物の表面の形状を平面に変形する座標変換に基づいて、前記三次元画像データから、前記第一の解剖学特徴における第一の面と前記第二の解剖学特徴における第二の面との距離が所定の条件を満たす位置の集合で構成される断面に対応する前記断面画像を生成する
    ことを特徴とする請求項1に記載の処理装置。
  6. 前記特徴取得手段は、前記三次元画像データから前記対象物における少なくとも前記第一の解剖学特徴と前記第二の解剖学特徴を取得する
    ことを特徴とする請求項1から5のいずれか1項に記載の処理装置。
  7. 前記特徴取得手段は、前記三次元画像データに対するユーザーの指示に基づいて、前記対象物における少なくとも前記第一の解剖学特徴と前記第二の解剖学特徴を取得することを特徴とする請求項1から5のいずれか1項に記載の処理装置。
  8. 前記第一の解剖学特徴と第二の解剖学特徴の少なくとも一つは、前記対象物の乳頭位置を含む
    ことを特徴とする請求項1記載の処理装置。
  9. 前記生成手段により生成された前記断面画像を表示する表示手段を更に有し、
    前記画像取得手段は、複数の前記三次元画像データを取得し、
    前記特徴取得手段は、複数の前記三次元画像データのそれぞれから前記対象物における少なくとも第一の解剖学特徴と第二の解剖学特徴を取得し、
    前記生成手段は、複数の前記三次元画像データのそれぞれから、前記第一の解剖学特徴における第一の面と前記第二の解剖学特徴における第二の面との距離が所定の条件を満たす位置の集合で構成される断面に対応する断面画像を生成し、
    前記表示手段は、前記生成手段により生成された複数の前記断面画像を並べて表示する
    ことを特徴とする請求項1に記載の処理装置。
  10. 対象物の三次元のボリュームデータである三次元画像データから断面画像を生成する処理方法であって、
    前記対象物における少なくとも第一の解剖学特徴と第二の解剖学特徴を取得し、
    前記三次元画像データから、前記第一の解剖学特徴における第一の面と前記第二の解剖学特徴における第二の面との距離が所定の条件を満たす位置の集合で構成される断面に対応する断面画像を生成する
    ことを特徴とする処理方法。
  11. 前記第一の解剖学特徴と第二の解剖学特徴それぞれ、前記対象物の体表面と前記対象物の大胸筋面であり
    前記三次元画像データから、前記対象物の体表面からの距離と前記対象物の大胸筋面からの距離が前記所定の条件を満たす曲面形状上の前記断面に対応する前記断面画像を生成する
    ことを特徴とする請求項10に記載の処理方法。
  12. 前記第一の解剖学特徴と第二の解剖学特徴それぞれ、前記対象物の体表面と前記対象物の大胸筋面であり
    前記三次元画像データから、前記対象物の体表面からの距離と前記対象物の大胸筋面からの距離の割合が一定である曲面形状上の断面に対応する前記断面画像を生成する
    ことを特徴とする請求項10に記載の処理方法。
  13. 前記第一の解剖学特徴と第二の解剖学特徴それぞれ、前記対象物の体表面と前記対象物の大胸筋面であり
    前記対象物の体表面を含む面と前記対象物の大胸筋面を含む面とが対向する矩形形状に前記三次元画像データを変形する座標変換に基づいて、前記三次元画像データから、前記対象物の体表面を含む面または前記対象物の大胸筋面を含む面と平行する断面に対応する前記断面画像を生成する
    ことを特徴とする請求項10に記載の処理方法。
  14. 前記三次元画像データ中の前記対象物の表面の形状を平面に変形する座標変換に基づいて、前記三次元画像データから、前記第一の解剖学特徴における第一の面と前記第二の解剖学特徴における第二の面との距離が所定の条件を満たす位置の集合で構成される断面に対応する前記断面画像を生成する
    ことを特徴とする請求項10に記載の処理方法。
  15. 前記三次元画像データから前記対象物の前記少なくとも第一の解剖学特徴と第二の解剖学特徴とを取得する
    ことを特徴とする請求項10から14のいずれか1項に記載の処理方法。
  16. 前記三次元画像データに対するユーザーの指示に基づいて、前記少なくとも第一の解剖学特徴と第二の解剖学特徴を取得する
    ことを特徴とする請求項10から14のいずれか1項に記載の処理方法。
  17. 前記第一の解剖学特徴と第二の解剖学特徴の少なくとも一つは、前記対象物の乳頭位置を含む
    ことを特徴とする請求項10に記載の処理方法。
  18. 複数の前記三次元画像データを取得し、
    複数の前記三次元画像データのそれぞれから前記対象物における少なくとも第一の解剖学特徴と第二の解剖学特徴を取得し、
    複数の前記三次元画像データのそれぞれから、前記第一の解剖学特徴における第一の面と前記第二の解剖学特徴における第二の面との距離が所定の条件を満たす位置の集合で構成される断面に対応する断面画像を生成し、
    複数の前記断面画像を並べて表示する
    ことを特徴とする請求項10に記載の処理方法。
  19. 請求項1からのいずれか1項に記載の処理装置としてコンピュータを機能させることを特徴とするプログラム。
JP2018145407A 2018-08-01 2018-08-01 処理装置、処理方法、およびプログラム Active JP6660428B2 (ja)

Priority Applications (1)

Application Number Priority Date Filing Date Title
JP2018145407A JP6660428B2 (ja) 2018-08-01 2018-08-01 処理装置、処理方法、およびプログラム

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
JP2018145407A JP6660428B2 (ja) 2018-08-01 2018-08-01 処理装置、処理方法、およびプログラム

Related Parent Applications (1)

Application Number Title Priority Date Filing Date
JP2014003699A Division JP6383153B2 (ja) 2014-01-10 2014-01-10 処理装置、処理方法、およびプログラム

Publications (2)

Publication Number Publication Date
JP2018167077A JP2018167077A (ja) 2018-11-01
JP6660428B2 true JP6660428B2 (ja) 2020-03-11

Family

ID=64018120

Family Applications (1)

Application Number Title Priority Date Filing Date
JP2018145407A Active JP6660428B2 (ja) 2018-08-01 2018-08-01 処理装置、処理方法、およびプログラム

Country Status (1)

Country Link
JP (1) JP6660428B2 (ja)

Family Cites Families (13)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP3758894B2 (ja) * 1999-04-03 2006-03-22 コニカミノルタホールディングス株式会社 マンモグラム画像診断支援装置
JP5121399B2 (ja) * 2007-11-06 2013-01-16 株式会社東芝 画像表示装置
JP5159301B2 (ja) * 2007-12-28 2013-03-06 株式会社東芝 医用画像表示装置および画像表示方法
JP5145170B2 (ja) * 2008-08-27 2013-02-13 富士フイルム株式会社 医用画像の評価装置
JP5538862B2 (ja) * 2009-12-18 2014-07-02 キヤノン株式会社 画像処理装置、画像処理システム、画像処理方法、及びプログラム
JP5606832B2 (ja) * 2010-03-05 2014-10-15 富士フイルム株式会社 画像診断支援装置、方法およびプログラム
JP6012931B2 (ja) * 2010-03-30 2016-10-25 東芝メディカルシステムズ株式会社 画像処理装置及び画像処理装置の制御プログラム
JP5737858B2 (ja) * 2010-04-21 2015-06-17 キヤノン株式会社 画像処理装置、画像処理方法、及びプログラム
CN102947841A (zh) * 2010-04-30 2013-02-27 沃康普公司 放射摄影图像中的毛刺恶性肿块检测和分级
EP2646980B1 (en) * 2010-11-30 2017-01-11 Volpara Health Technologies Limited An imaging technique and imaging system
JP2012213558A (ja) * 2011-04-01 2012-11-08 Canon Inc 画像処理装置、画像処理方法およびプログラム
US8417005B1 (en) * 2011-09-15 2013-04-09 Sunnybrook Health Sciences Centre Method for automatic three-dimensional segmentation of magnetic resonance images
US9058647B2 (en) * 2012-01-16 2015-06-16 Canon Kabushiki Kaisha Information processing apparatus, information processing method, and storage medium

Also Published As

Publication number Publication date
JP2018167077A (ja) 2018-11-01

Similar Documents

Publication Publication Date Title
JP6346445B2 (ja) 処理装置、処理装置の制御方法、およびプログラム
JP5546230B2 (ja) 情報処理装置、情報処理方法、及びプログラム
US8908944B2 (en) Information processing apparatus, information processing method, and program
JP6000705B2 (ja) データ処理装置及びデータ処理方法
KR101267759B1 (ko) 화상 처리 장치, 화상 처리 방법 및 저장 매체
US8437521B2 (en) Systems and methods for automatic vertebra edge detection, segmentation and identification in 3D imaging
US9218542B2 (en) Localization of anatomical structures using learning-based regression and efficient searching or deformation strategy
EP3625768B1 (en) Determining a clinical target volume
KR102509659B1 (ko) 딥러닝을 이용한 3차원 전신 골격 모델 생성 장치 및 방법
JP5832938B2 (ja) 画像処理装置、方法及びプログラム
JP7214434B2 (ja) 医用画像処理装置及び医用画像処理プログラム
JP6383153B2 (ja) 処理装置、処理方法、およびプログラム
JP6541334B2 (ja) 画像処理装置、画像処理方法、およびプログラム
RU2669680C2 (ru) Инициализация модели на основе классификации видов
US9129392B2 (en) Automatic quantification of mitral valve dynamics with real-time 3D ultrasound
JP6905323B2 (ja) 画像処理装置、画像処理方法、およびプログラム
US9286688B2 (en) Automatic segmentation of articulated structures
JP2022111705A (ja) 学習装置、画像処理装置、医用画像撮像装置、学習方法およびプログラム
JP2022111704A (ja) 画像処理装置、医用画像撮像装置、画像処理方法、およびプログラム
JP6660428B2 (ja) 処理装置、処理方法、およびプログラム
US11341661B2 (en) Method and apparatus for registering live medical image with anatomical model
JP2019500114A (ja) 位置合わせ精度の決定
Mao Three-dimensional Ultrasound Fusion for Transesophageal Echocardiography
WO2020138136A1 (ja) 画像処理装置、画像処理方法及びプログラム
JP2022111706A (ja) 画像処理装置および画像処理方法、医用撮像装置、プログラム

Legal Events

Date Code Title Description
A521 Request for written amendment filed

Free format text: JAPANESE INTERMEDIATE CODE: A523

Effective date: 20180827

A621 Written request for application examination

Free format text: JAPANESE INTERMEDIATE CODE: A621

Effective date: 20180827

A977 Report on retrieval

Free format text: JAPANESE INTERMEDIATE CODE: A971007

Effective date: 20190520

A131 Notification of reasons for refusal

Free format text: JAPANESE INTERMEDIATE CODE: A131

Effective date: 20190701

A521 Request for written amendment filed

Free format text: JAPANESE INTERMEDIATE CODE: A523

Effective date: 20190805

TRDD Decision of grant or rejection written
A01 Written decision to grant a patent or to grant a registration (utility model)

Free format text: JAPANESE INTERMEDIATE CODE: A01

Effective date: 20200110

A61 First payment of annual fees (during grant procedure)

Free format text: JAPANESE INTERMEDIATE CODE: A61

Effective date: 20200207

R151 Written notification of patent or utility model registration

Ref document number: 6660428

Country of ref document: JP

Free format text: JAPANESE INTERMEDIATE CODE: R151