JP6661196B2 - Image processing apparatus, method, and program - Google Patents
Image processing apparatus, method, and program Download PDFInfo
- Publication number
- JP6661196B2 JP6661196B2 JP2016544238A JP2016544238A JP6661196B2 JP 6661196 B2 JP6661196 B2 JP 6661196B2 JP 2016544238 A JP2016544238 A JP 2016544238A JP 2016544238 A JP2016544238 A JP 2016544238A JP 6661196 B2 JP6661196 B2 JP 6661196B2
- Authority
- JP
- Japan
- Prior art keywords
- shape
- subject
- input image
- eigenspace
- specific object
- 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
Links
- 238000012545 processing Methods 0.000 title claims description 88
- 238000000034 method Methods 0.000 title claims description 67
- 230000006870 function Effects 0.000 claims description 109
- 230000011218 segmentation Effects 0.000 claims description 94
- 238000005070 sampling Methods 0.000 claims description 30
- 238000010845 search algorithm Methods 0.000 claims description 17
- 239000000284 extract Substances 0.000 claims description 13
- 238000003672 processing method Methods 0.000 claims description 12
- 230000014509 gene expression Effects 0.000 claims description 5
- 210000002487 multivesicular body Anatomy 0.000 claims description 4
- 210000004027 cell Anatomy 0.000 claims 2
- 238000004422 calculation algorithm Methods 0.000 description 24
- 210000000496 pancreas Anatomy 0.000 description 19
- 230000008569 process Effects 0.000 description 17
- 238000010586 diagram Methods 0.000 description 16
- 238000004364 calculation method Methods 0.000 description 15
- 230000003187 abdominal effect Effects 0.000 description 11
- 238000004891 communication Methods 0.000 description 9
- 210000004185 liver Anatomy 0.000 description 8
- 238000013507 mapping Methods 0.000 description 8
- 238000005457 optimization Methods 0.000 description 8
- 238000007796 conventional method Methods 0.000 description 7
- 238000007781 pre-processing Methods 0.000 description 6
- 230000005540 biological transmission Effects 0.000 description 4
- 238000000513 principal component analysis Methods 0.000 description 4
- 238000013179 statistical model Methods 0.000 description 4
- 238000012886 linear function Methods 0.000 description 3
- 238000003860 storage Methods 0.000 description 3
- 238000004458 analytical method Methods 0.000 description 2
- 230000008859 change Effects 0.000 description 2
- 238000002591 computed tomography Methods 0.000 description 2
- 238000010191 image analysis Methods 0.000 description 2
- 210000000056 organ Anatomy 0.000 description 2
- 210000003240 portal vein Anatomy 0.000 description 2
- 238000001356 surgical procedure Methods 0.000 description 2
- 230000009471 action Effects 0.000 description 1
- 238000013459 approach Methods 0.000 description 1
- 230000001174 ascending effect Effects 0.000 description 1
- 230000008901 benefit Effects 0.000 description 1
- 230000000747 cardiac effect Effects 0.000 description 1
- 238000013461 design Methods 0.000 description 1
- 238000009826 distribution Methods 0.000 description 1
- 230000000694 effects Effects 0.000 description 1
- 238000002474 experimental method Methods 0.000 description 1
- 238000003384 imaging method Methods 0.000 description 1
- 230000002452 interceptive effect Effects 0.000 description 1
- 238000002955 isolation Methods 0.000 description 1
- 210000004072 lung Anatomy 0.000 description 1
- 238000004519 manufacturing process Methods 0.000 description 1
- 238000012986 modification Methods 0.000 description 1
- 230000004048 modification Effects 0.000 description 1
- 238000003909 pattern recognition Methods 0.000 description 1
- 230000009467 reduction Effects 0.000 description 1
- 238000012552 review Methods 0.000 description 1
- 210000000952 spleen Anatomy 0.000 description 1
- 238000007619 statistical method Methods 0.000 description 1
- 238000012800 visualization Methods 0.000 description 1
Classifications
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T7/00—Image analysis
Landscapes
- Engineering & Computer Science (AREA)
- Computer Vision & Pattern Recognition (AREA)
- Physics & Mathematics (AREA)
- General Physics & Mathematics (AREA)
- Theoretical Computer Science (AREA)
- Image Analysis (AREA)
Description
本発明は、画像処理装置、方法、及びプログラムに係り、特に、被写体の領域を抽出する画像処理装置、方法、及びプログラムに関するものである。 The present invention relates to an image processing device, a method, and a program, and more particularly, to an image processing device, a method, and a program for extracting an area of a subject.
デジタル画像の中からコンピュータを用いて図形を認識する処理はセグメンテーションと呼ばれるが、対象図形のSNが低く、かつ、形状が統計的に変動する場合、正確なセグメンテーションは非常に難しくなる。 The process of recognizing a figure from a digital image using a computer is called segmentation. If the SN of the target figure is low and the shape fluctuates statistically, accurate segmentation becomes extremely difficult.
グラフカット(例えば、Boykov, Y., Veksler, O., Zabih, R.,“Fast approximate energy minimization via graph cuts.”, Pattern Analysis and Machine Intelligence, IEEE Transactions on 23 (11), 2001, p.1222‐1239.を参照。)などの最適化理論に基づくセグメンテーションアルゴリズムは、目的関数を真に最適化できるという意味で従来のアルゴリズムより優れており、SNが低い場合にもしばしば非常に高い性能を示した。しかし、それでも形状が変動する場合には正確に認識できないことがしばしばあった。 Graph cuts (eg, Boykov, Y., Veksler, O., Zabih, R., “Fast approximate energy minimization via graph cuts.”, Pattern Analysis and Machine Intelligence, IEEE Transactions on 23 (11), 2001, p. 1222 Segmentation algorithms based on optimization theory such as -1239.) Are superior to conventional algorithms in that they can truly optimize the objective function, and often exhibit very high performance even at low SN. Was. However, even when the shape fluctuates, it is often impossible to accurately recognize the shape.
そこで、最近は、最適化理論に基づくアルゴリズムにおいて図形の形状情報を利用する方法が主流になってきたが、それらは、大まかに次の3通りに分けられる。 Therefore, recently, methods using shape information of a figure in an algorithm based on an optimization theory have become mainstream, but they are roughly divided into the following three types.
1つ目の方法としては、単一の形状テンプレートを利用する方法である。例えば、一般的形状テンプレートとして、楕円形状のテンプレート(例えば、Slabaugh, G., Unal, G., Sept, “Graph cuts segmentation using an elliptical shape prior.”, In: IEEE International Conference on Image Processing, 2005, Vol. 2. p.II‐1222‐5.を参照。)を用いる手法や、塊状図形のテンプレート(例えば、Funka-Lea, G., Boykov, Y., Florin, C., Jolly, M.‐P., Moreau-Gobard, R., Ramaraj, R., Rinck, D., 2006. “Automatic heart isolation for CT coronary visualization using graph-cuts.”, In: Biomedical Imaging: Nano to Macro, 3rd IEEE International Symposium on. IEEE, 2006, p.614‐617.を参照。)を利用する方法が知られている。
また、特定の形状テンプレートとして、ユーザ定義の任意形状(例えば、Freedman, D., Zhang, T., “Interactive graph cut based segmentation with shape priors.”, In: IEEE Computer Society Conference on Computer Vision and Pattern Recognition. Vol. 1. IEEE, 2005, p.755‐762.を参照。)を用いる手法や、統計モデルから選択した形状(例えば、Grosgeorge, D., Petitjean, C., Dacher, J.-N., Ruan, S., “Graph cut segmentation with a statistical shape model in cardiac mri.”, Computer Vision and Image Understanding 117 (9), 2013, p.1027‐1035.,Akinobu Shimizu, Keita Nakagomi, Takuya Narihira, Hidefumi Kobatake, Shigeru Nawano, Kenji Shinozaki, Koich Ishizu, and Kaori Togashi, “Automated Segmentation of 3D CT Images based on Statistical Atlas and Graph Cuts”, Proc. of MICCAI workshop MCV, 2010, p.129-138.,Malcolm, J., Rathi, Y., Tannenbaum, A., “Graph cut segmentation with nonlinear shape priors.”, In: Image Processing, 2007. ICIP 2007. IEEE International Conference on. Vol. 4. IEEE, 2007, p.IV‐365.を参照。)を用いる手法が知られている。The first method is to use a single shape template. For example, as a general shape template, an elliptical template (for example, Slabaugh, G., Unal, G., Sept, “Graph cuts segmentation using an elliptical shape prior.”, In: IEEE International Conference on Image Processing, 2005, Vol. 2. p.II-1222-5.), Or a template of a massive figure (for example, Funka-Lea, G., Boykov, Y., Florin, C., Jolly, M.- P., Moreau-Gobard, R., Ramaraj, R., Rinck, D., 2006. “Automatic heart isolation for CT coronary visualization using graph-cuts.”, In: Biomedical Imaging: Nano to Macro, 3rd IEEE International Symposium on. IEEE, 2006, pp. 614-617.).
In addition, as a specific shape template, a user-defined arbitrary shape (for example, Freedman, D., Zhang, T., “Interactive graph cut based segmentation with shape priors.”, In: IEEE Computer Society Conference on Computer Vision and Pattern Recognition) Vol. 1. IEEE, 2005, pp. 755-762.) And shapes selected from statistical models (eg, Grosgeorge, D., Petitjean, C., Dacher, J.-N. , Ruan, S., “Graph cut segmentation with a statistical shape model in cardiac mri.”, Computer Vision and Image Understanding 117 (9), 2013, pp. 1027-1035., Akinobu Shimizu, Keita Nakagomi, Takuya Narihira, Hidefumi Kobatake, Shigeru Nawano, Kenji Shinozaki, Koich Ishizu, and Kaori Togashi, “Automated Segmentation of 3D CT Images based on Statistical Atlas and Graph Cuts”, Proc. Of MICCAI workshop MCV, 2010, p.129-138., Malcolm, J ., Rathi, Y., Tannenbaum, A., “Graph cut segmentation with nonlinear shape priors.”, In: Image Processing, 2007. A method using ICIP 2007. IEEE International Conference on Vol. 4. IEEE, 2007, p.IV-365.) Is known.
2つ目の方法としては、複数(数個)の形状テンプレートや形状の確率的表現を利用する方法である。例えば、統計モデルから選択した複数形状(例えば、Nakagomi, K., Shimizu, A., Kobatake, H., Yakami, M., Fujimoto, K., Togashi, K., “Multi-shape graph cuts with neighbor prior constraints and its application to lung segmentation from a chest CT volume.”, Medical image analysis 17 (1), 2013, p.62‐77.を参照。)を用いる手法や、形状の確率的表現(例えば、Linguraru, M. G., Pura, J. A., Pamulapati, V., Summers, R. M., “Statistical 4d graphs for multi-organ abdominal segmentation from multiphase CT.”, Medical image analysis 16 (4), 2012, p.904‐914.を参照。)を用いる手法が知られている。 A second method is to use a plurality of (several) shape templates or stochastic expressions of shapes. For example, a plurality of shapes selected from a statistical model (for example, Nakagomi, K., Shimizu, A., Kobatake, H., Yakami, M., Fujimoto, K., Togashi, K., “Multi-shape graph cuts with neighbor prior constraints and its application to lung segmentation from a chest CT volume. ”, Medical image analysis 17 (1), 2013, p.62-77.), and stochastic expression of shapes (for example, Linguraru , MG, Pura, JA, Pamulapati, V., Summers, RM, “Statistical 4d graphs for multi-organ abdominal segmentation from multiphase CT.”, Medical image analysis 16 (4), 2012, p.904-914. ) Is known.
3つ目の方法としては、大量(数十以上)の形状テンプレートを利用する方法である。例えば、特定の問題用(例えば、Kohli, P., Rihan, J., Bray, M., Torr, P. H., “Simultaneous segmentation and pose estimation of humans using dynamic graph cuts.”, International Journal of Computer Vision 79 (3), 2008, p.285‐298.を参照。)のテンプレートを用いる手法や、任意の大量の図形形状(例えば、米国特許第8249349号明細書を参照。)を用いる手法が知られている。
A third method is to use a large number (several tens or more) of shape templates. For example, for specific problems (eg, Kohli, P., Rihan, J., Bray, M., Torr, PH, “Simultaneous segmentation and pose estimation of humans using dynamic graph cuts.”, International Journal of Computer Vision 79 ( 3), 2008, pp. 285-298.) And a method using an arbitrary large number of graphic shapes (for example, see US Pat. No. 8,249,349). .
上述した従来の技術は、1つ目の手法から3つ目の手法の順に統計的に変動をする図形をうまく扱えるようになる。最も先端的な方法は、上記3つ目の手法の、任意の大量の図形形状(上記米国特許第8249349号明細書)である。しかし、この方法も含めてすべての方法は事前に形状を選択するなどして、あらかじめ形状集合を準備する必要がある。そのため、事前に選択した形状集合の中にそもそも最適な形状が含まれていなければ、セグメンテーションの最適性は保証されず、性能は低下する。 The above-described conventional technique can handle a graphic that fluctuates statistically in order from the first technique to the third technique. The most advanced method is the arbitrarily large number of graphic shapes of the third method (U.S. Pat. No. 8,249,349). However, in all methods including this method, it is necessary to prepare a shape set in advance by selecting a shape in advance. Therefore, unless the optimal shape is originally included in the shape set selected in advance, the optimality of the segmentation is not guaranteed, and the performance is reduced.
一方、デジタル画像は有限個の画素の集合であることから、可能な全ての形状パターンを準備する方法も原理的には考えられ、その場合には最適性は保証される。しかし、そのために用意しなければならない形状数は膨大になり、前処理も含めた処理全体の必要なメモリと計算時間の観点から現実的ではない。例えば、従来最も優れていた上記米国特許第8249349号明細書のアルゴリズムでも、2次元画像で107個程度の数の形状テンプレートを扱うのが限界であった。しかし、統計モデルが生成しうる可能な形状パターン数は109個以上ある場合も多く、従来のアルゴリズムでは計算コストの観点から非常に困難である。さらに医用画像では3次元画像が主流であるが、その場合には必要なメモリは二桁以上増加し、ほぼ従来のアルゴリズムによって現実的な時間やメモリサイズで処理を実行することは不可能であった。On the other hand, since a digital image is a set of a finite number of pixels, a method of preparing all possible shape patterns can be considered in principle, and in that case, optimality is guaranteed. However, the number of shapes that must be prepared for that purpose is enormous, and it is not realistic from the viewpoint of the memory required for the entire processing including preprocessing and the calculation time. For example, even in a conventional best it has been the U.S. Patent No. 8249349 Pat algorithms that handle 107 or so the number of shape template with the two-dimensional image was limited. However, the number of possible shapes pattern statistical model may be generated in many cases there 10 9 or more, it is very difficult in terms of computational cost in the conventional algorithm. Furthermore, three-dimensional images are mainly used for medical images, but in that case, the required memory increases by more than two orders of magnitude, and it is impossible to execute the processing in a practically time and memory size using almost conventional algorithms. Was.
本発明の一実施形態は、上記事情に鑑みてなされたものである。 One embodiment of the present invention has been made in view of the above circumstances.
上記目的を達成するために、第1の態様の画像処理装置は、特定の物体である被写体を表す入力画像から、前記被写体の領域を抽出する画像処理装置であって、前記入力画像を受け付ける受付手段と、(A)前記被写体を表す学習用の複数の画像であって、前記被写体の領域が予め求められた前記複数の画像に基づいて予め計算された固有ベクトルを基底とする固有空間であって、かつ、(B)固有空間上の点が、前記特定の物体の形状の統計的変動を表す統計的形状モデルの形状パラメータを示す固有空間において、前記入力画像に基づいて、前記固有空間上の点が示す前記形状パラメータが表す前記特定の物体の形状の尤もらしさと前記入力画像中の隣接画素間の画素値の差とに応じた値を表す予め定められた目的関数を最適化するように、前記入力画像が表す前記被写体の形状を表す前記形状パラメータを推定し、前記推定された前記形状パラメータが表す前記特定の物体の形状を事前知識として、前記被写体の領域を、前記入力画像から抽出するセグメンテーション手段と、を備えている。 In order to achieve the above object, an image processing apparatus according to a first aspect is an image processing apparatus that extracts an area of a subject from an input image representing a subject that is a specific object, and that receives an input image. Means, and (A) a plurality of learning images representing the subject, wherein the subject region is an eigenspace based on eigenvectors calculated in advance based on the plurality of images obtained in advance. And (B) a point on the eigenspace indicates a shape parameter of a statistical shape model representing a statistical variation of the shape of the specific object in the eigenspace based on the input image. A predetermined objective function representing a value corresponding to the likelihood of the shape of the specific object represented by the shape parameter indicated by a point and a pixel value difference between adjacent pixels in the input image is optimized. Estimating the shape parameter representing the shape of the subject represented by the input image, and extracting the area of the subject from the input image using the shape of the specific object represented by the estimated shape parameter as prior knowledge. And a segmentation means for performing the segmentation.
また、第2の態様の画像処理方法は、受付手段、及びセグメンテーション手段を含み、特定の物体である被写体を表す入力画像から、前記被写体の領域を抽出する画像処理装置における画像処理方法であって、前記受付手段が、前記入力画像を受け付けるステップと、前記セグメンテーション手段が、(A)前記被写体を表す学習用の複数の画像であって、前記被写体の領域が予め求められた前記複数の画像に基づいて予め計算された固有ベクトルを基底とする固有空間であって、かつ、(B)固有空間上の点が、前記特定の物体の形状の統計的変動を表す統計的形状モデルの形状パラメータを示す固有空間において、前記入力画像に基づいて、前記固有空間上の点が示す前記形状パラメータが表す前記特定の物体の形状の尤もらしさと前記入力画像中の隣接画素間の画素値の差とに応じた値を表す予め定められた目的関数を最適化するように、前記入力画像が表す前記被写体の形状を表す前記形状パラメータを推定し、前記推定された前記形状パラメータが表す前記特定の物体の形状を事前知識として、前記被写体の領域を、前記入力画像から抽出するステップと、を含む。 An image processing method according to a second aspect is an image processing method in an image processing apparatus that includes a receiving unit and a segmentation unit, and extracts an area of the subject from an input image representing the subject as a specific object. Wherein the accepting means accepts the input image; and wherein the segmentation means includes: (A) a plurality of learning images representing the subject, wherein the region of the subject is obtained in advance; (B) a point on the eigenspace indicates a shape parameter of a statistical shape model representing a statistical variation of the shape of the specific object, based on an eigenvector calculated in advance based on the eigenvector. In the eigenspace, based on the input image, the likelihood of the shape of the specific object represented by the shape parameter indicated by a point on the eigenspace and the Estimating the shape parameter representing the shape of the subject represented by the input image, so as to optimize a predetermined objective function representing a value corresponding to a pixel value difference between adjacent pixels in the force image, Extracting the region of the subject from the input image, using the shape of the specific object represented by the estimated shape parameter as prior knowledge.
また、第3の態様の前記セグメンテーション手段は、前記固有空間において、前記入力画像に基づいて、前記目的関数を最適化するように、前記入力画像が表す前記被写体の形状を表す前記形状パラメータを推定すると同時に、前記推定された前記形状パラメータが表す前記特定の物体の形状を事前知識として、前記被写体の領域を、前記入力画像から抽出するようにすることができる。 In the third aspect, the segmentation unit estimates the shape parameter representing the shape of the subject represented by the input image in the eigenspace based on the input image so as to optimize the objective function. At the same time, the region of the subject can be extracted from the input image using the shape of the specific object represented by the estimated shape parameter as prior knowledge.
また、第4の態様の前記セグメンテーション手段は、最適な形状パラメータを示す点が表す前記特定の物体の形状を含む形状集合を表す前記固有空間上の凸多胞体を繰り返し分割して最適な形状パラメータを示す点を含む前記固有空間上の凸多胞体を探索する探索アルゴリズムに従って、前記目的関数を最適化するように、前記入力画像が表す前記被写体の形状を表す前記形状パラメータを推定し、前記推定された前記形状パラメータが表す前記特定の物体の形状を事前知識として、前記被写体の領域を、前記入力画像から抽出するようにすることができる。 Further, the segmentation means of the fourth aspect is configured to repeatedly divide the convex multivesicular body on the eigenspace representing a shape set including a shape of the specific object represented by a point indicating an optimal shape parameter to obtain an optimal shape parameter. Estimating the shape parameters representing the shape of the subject represented by the input image so as to optimize the objective function, according to a search algorithm for searching for convex polyvesicles on the eigenspace including points indicating The region of the subject may be extracted from the input image using the shape of the specific object represented by the obtained shape parameter as prior knowledge.
また、第5の態様の前記セグメンテーション手段は、前記探索アルゴリズムにおいて、前記凸多胞体の各頂点を調べることで、前記凸多胞体に含まれる形状集合に対する前記目的関数の下界を計算することにより、最適な形状パラメータを示す点を含む前記固有空間上の凸多胞体を探索して前記入力画像が表す前記被写体の形状を表す前記形状パラメータを推定し、前記推定された前記形状パラメータが表す前記特定の物体の形状を事前知識として、前記被写体の領域を、前記入力画像から抽出するようにすることができる。 Further, the segmentation means of the fifth aspect, in the search algorithm, by examining each vertex of the convex polycyst, by calculating the lower bound of the objective function for the shape set included in the convex polycyst, A search is made for convex polyvesicles on the eigenspace including points indicating the optimal shape parameters to estimate the shape parameters representing the shape of the subject represented by the input image, and the identification represented by the estimated shape parameters. The region of the subject may be extracted from the input image using the shape of the object as prior knowledge.
また、第6の態様の前記セグメンテーション手段は、分割後の2つの凸多胞体の体積が対応するように、前記固有空間上の凸多胞体を分割し、かつ、最適な形状パラメータを示す点を含む前記固有空間上の凸多胞体を探索することを繰り返して、前記入力画像が表す前記被写体の形状を表す前記形状パラメータを推定し、前記推定された前記形状パラメータが表す前記特定の物体の形状を事前知識として、前記被写体の領域を、前記入力画像から抽出するようにすることができる。 In the sixth aspect, the segmentation means divides the convex polyvesicles in the eigenspace so that the volumes of the two convex polyvesicles after division correspond to each other, and sets a point indicating an optimal shape parameter. Iterating the search for the convex multivesicular body on the eigenspace including, estimating the shape parameter representing the shape of the subject represented by the input image, the shape of the specific object represented by the estimated shape parameter May be extracted as the prior knowledge from the input image.
また、第7の態様の前記セグメンテーション手段は、前記固有空間上にサンプリング点を設定し、前記固有空間上に設定されたサンプリング点から決定される超平面を用いて前記凸多胞体を分割して最適な形状パラメータを示す点を含む前記固有空間上の凸多胞体を探索することを繰り返して、前記入力画像が表す前記被写体の形状を表す前記形状パラメータを推定し、前記推定された前記形状パラメータが表す前記特定の物体の形状を事前知識として、前記被写体の領域を、前記入力画像から抽出するようにすることができる。 In the seventh aspect, the segmentation means sets a sampling point on the eigenspace, and divides the convex polymorph using a hyperplane determined from the sampling points set on the eigenspace. Iteratively searching for convex polyvesicles on the eigenspace including points indicating the optimal shape parameters, estimates the shape parameters representing the shape of the subject represented by the input image, and estimates the shape parameters The area of the subject may be extracted from the input image using the shape of the specific object represented by as the prior knowledge.
また、第8の態様の前記セグメンテーション手段は、任意に、前記固有空間上にサンプリング点を設定するようにすることができる。 Further, the segmentation means of the eighth aspect may arbitrarily set a sampling point on the eigenspace.
また、第9の態様の前記目的関数は、前記特定の物体の形状の尤もらしさとして、前記形状パラメータが表す前記特定の物体の形状に対する画素における形状ラベルの値に関して単調に変化する単調関数を含むようにすることができる。 Further, the objective function of the ninth aspect includes a monotonic function that monotonically changes as a likelihood of the shape of the specific object with respect to a value of a shape label in a pixel corresponding to the shape of the specific object represented by the shape parameter. You can do so.
また、第10の態様の前記探索アルゴリズムは、Branch and bound法及びグラフカット法を用いることができる。 Further, the search algorithm according to the tenth aspect can use a Branch and bound method and a graph cut method.
また、第11の態様のプログラムは、コンピュータを、上記の画像処理装置の各手段として機能させるためのプログラムである。 Further, a program according to an eleventh aspect is a program for causing a computer to function as each unit of the image processing apparatus.
実施の形態に係る画像処理装置、方法、及びプログラムによれば、(A)被写体の領域が予め求められた複数の画像に基づいて予め計算された固有ベクトルを基底とする固有空間であって、(B)固有空間上の点が、特定の物体の形状を表す統計的形状モデルの形状パラメータを示す固有空間において、固有空間上の点が示す形状パラメータが表す特定の物体の形状の尤もらしさと入力画像中の隣接画素間の画素値の差とに応じた値を表す予め定められた目的関数を最適化するように、入力画像が表す被写体の形状を表す形状パラメータを推定する。推定された形状パラメータが表す特定の物体の形状を事前知識として、被写体の領域を、入力画像から抽出する。これにより、計算量の増大を抑制して、被写体の領域を精度よく抽出することができる、という効果が得られる。 According to the image processing apparatus, method, and program according to the embodiment, (A) the eigenspace based on the eigenvector calculated in advance based on a plurality of images obtained in advance, B) In the eigenspace in which the point on the eigenspace indicates the shape parameter of the statistical shape model representing the shape of the specific object, the likelihood and input of the shape of the specific object represented by the shape parameter indicated by the point on the eigenspace A shape parameter representing a shape of a subject represented by the input image is estimated so as to optimize a predetermined objective function representing a value corresponding to a pixel value difference between adjacent pixels in the image. The region of the subject is extracted from the input image using the shape of the specific object represented by the estimated shape parameter as prior knowledge. As a result, it is possible to obtain an effect of suppressing an increase in the amount of calculation and extracting a region of a subject with high accuracy.
以下、図面を参照して、実施の形態を詳細に説明する。なお、実施の形態では、特定の物体の統計的変動を表す図形の形状の統計モデル(以下、統計的形状モデルと称する。)を生成する統計的形状モデル生成装置10と、特定の物体である被写体を表す入力画像から、被写体の領域を抽出する画像処理装置100とについて説明する。また、実施の形態では、3次元腹部CT画像を入力画像とし、当該入力画像から、膵臓の領域を抽出する場合を例に説明する。 Hereinafter, embodiments will be described in detail with reference to the drawings. In the embodiment, the statistical shape model generation device 10 that generates a statistical model (hereinafter, referred to as a statistical shape model) of the shape of a figure representing the statistical variation of a specific object, and the specific object. An image processing apparatus 100 that extracts a subject area from an input image representing the subject will be described. In the embodiment, a case will be described as an example where a three-dimensional abdominal CT image is used as an input image and a pancreas region is extracted from the input image.
<概要>
本実施の形態では、統計的形状モデルを用いるセグメンテーションアルゴリズムを提案する。このアルゴリズムの特徴は以下の5点である。<Overview>
In the present embodiment, a segmentation algorithm using a statistical shape model is proposed. This algorithm has the following five features.
1.セグメンテーションの過程で、統計的形状モデルを用いて生成可能な全ての3次元形状(109個以上)を考慮可能であり、セグメンテーションの目的関数の観点から最適な形状を選択可能である。1. In the course of the segmentation are possible considering all three-dimensional shape can be generated using a statistical shape model (10 9 or more), it is possible to select the optimal shape in terms of the objective function segmentation.
2.大量の形状テンプレートを用いる従来手法では必須で、特に計算コストが高い前処理(形状の生成・選択とクラスタリング)が不要である。 2. The conventional method using a large number of shape templates is indispensable, and preprocessing (generation and selection of shapes and clustering), which is particularly expensive, is unnecessary.
3.Branch and bound法とグラフカット法との組み合わせに限らず、最適化理論に基づくアルゴリズムであれば適用可能である。 3. The algorithm is not limited to the combination of the branch and bound method and the graph cut method, but can be applied to any algorithm based on optimization theory.
4.異なる統計的形状モデルや固有空間の分割法に対しても適用可能である。 4. It can be applied to different statistical shape models and eigenspace division methods.
5.3次元腹部CT画像内の膵臓のセグメンテーションに関して世界最高の精度が得られる。 5. World's highest accuracy for pancreas segmentation in 3D abdominal CT images.
図19に従来技術の処理の概要を示す。図19は、上記米国特許第8249349号明細書の手法を、統計的形状モデルを利用した膵臓のセグメンテーションに適用させた場合の処理である。図19に示すように、従来手法では、形状テンプレート集合T(⊂S)の生成処理、及び集合Tに対するクラスタリング処理が必要であるが、本実施の形態では、上述したように、大量の形状テンプレートを用いる従来手法で必要であった前処理が不要である。なお、従来手法(上記米国特許第8249349号明細書)で扱われていた画像は2次元であり、形状の最大数は107であるのに対し、本実施の形態で用いるアルゴリズムは3次元形状(109個以上)を考慮可能である。FIG. 19 shows an outline of the processing of the related art. FIG. 19 shows processing in a case where the technique of the above-mentioned US Pat. No. 8,249,349 is applied to pancreas segmentation using a statistical shape model. As shown in FIG. 19, the conventional method requires a process of generating a shape template set T (⊂S) and a clustering process for the set T. In the present embodiment, as described above, a large number of shape templates are used. This eliminates the need for the pre-processing required by the conventional method using. Incidentally, the conventional method image which has been treated with (aforementioned U.S. Pat. No. 8,249,349) a two-dimensional, while the maximum number of shapes is 107, the algorithm used in this embodiment is a three-dimensional shape (10 9 or more) can be considered.
また、実施の形態では、統計的形状モデルとしてレベルセット分布モデル(Level set distribution model; LSDM。符号付距離モデル(Signed distance model)と呼ばれることもある)を利用する。以下では、まず、統計的形状モデル生成装置10において統計的形状モデル(LSDM)の固有ベクトルと固有値とを算出した後で、画像処理装置100において統計的形状モデルの形状パラメータを推定し、推定された形状パラメータを用いた統計的形状モデルに基づいた最適セグメンテーションアルゴリズムを示す。 In the embodiment, a level set distribution model (LSDM; sometimes referred to as a signed distance model) is used as a statistical shape model. In the following, first, the eigenvectors and eigenvalues of the statistical shape model (LSDM) are calculated in the statistical shape model generation device 10, and then the shape parameters of the statistical shape model are estimated in the image processing device 100, and the estimation is performed. 4 shows an optimal segmentation algorithm based on a statistical shape model using shape parameters.
[第1の実施の形態]
<第1の実施の形態に係る統計的形状モデル生成装置10の構成>
第1の実施の形態に係る統計的形状モデル生成装置10は、図1に示される機能ブロックで表すことができる。また、これらの機能ブロックは、図2に示されるコンピュータのハードウェア構成により実現することができる。図2を参照してコンピュータの構成を説明する。[First Embodiment]
<Configuration of Statistical Shape Model Generator 10 According to First Embodiment>
The statistical shape model generation device 10 according to the first embodiment can be represented by the functional blocks shown in FIG. These functional blocks can be realized by the hardware configuration of the computer shown in FIG. The configuration of the computer will be described with reference to FIG.
図2に示す第1の実施の形態に係る統計的形状モデル生成装置10は、プログラムに基づき統計的形状モデル生成装置10の本実施の形態に係る処理を行うCPU(Central Processing Unit;中央処理装置)21と、CPU21による各種プログラムの実行時のワークエリア等として用いられるRAM(Random Access Memory)22と、各種制御プログラムや各種パラメータ等が予め記憶された記録媒体であるROM(Read Only Memory)23と、各種情報を記憶するために用いられるハードディスク24(図中「HDD」と記載)と、キーボードやマウス等からなる入力装置25と、ディスプレイ等からなる表示装置26と、LAN(Local Area Network)等を用いて通信を行う通信装置27と、外部に接続された画像情報提供装置30との間の各種情報の授受を司る入出力インタフェース部(図中、「外部IF」と記載)28と、を備えており、これらがシステムバスBUS29により相互に接続されて構成されている。 A statistical shape model generation device 10 according to the first embodiment shown in FIG. 2 is a CPU (Central Processing Unit; Central Processing Unit) that performs processing according to the present embodiment of the statistical shape model generation device 10 based on a program. ) 21, a RAM (Random Access Memory) 22 used as a work area or the like when the CPU 21 executes various programs, and a ROM (Read Only Memory) 23 which is a recording medium in which various control programs, various parameters, and the like are stored in advance. And a hard disk 24 (described as “HDD” in the figure) used for storing various information, an input device 25 including a keyboard and a mouse, a display device 26 including a display, and a LAN (Local Area Network). And the like, and exchanges various kinds of information between the communication device 27 that performs communication using the communication device 27 and the image information providing device 30 that is connected to the outside. (In the figure, "external IF" and described) output interface unit 28 includes a, it is configured are connected to each other by a system bus BUS29.
CPU21は、RAM22、ROM23、及びハードディスク24に対するアクセス、入力装置25を介した各種情報の取得、表示装置26に対する各種情報の表示、通信装置27を用いた各種情報の通信処理、及び入出力インタフェース部28に接続された画像情報提供装置30を含む外部装置からの情報の入力等を、各々行うことができる。 The CPU 21 accesses the RAM 22, the ROM 23, and the hard disk 24, acquires various information via the input device 25, displays various information on the display device 26, performs communication processing of various information using the communication device 27, and an input / output interface unit. Input of information from an external device including the image information providing device 30 connected to the MFP 28 can be performed.
CPU21が、ハードディスク24に記憶された本実施形態に係る処理を制御するプログラムを、RAM22に読み込み実行することにより、図1に示す本実施の形態に係る統計的形状モデル生成装置10における、図1に示す各処理部の機能が実行される。 The CPU 21 reads the program for controlling the processing according to the present embodiment stored in the hard disk 24 into the RAM 22 and executes the program, thereby causing the statistical shape model generating apparatus 10 according to the present embodiment shown in FIG. Are executed.
このようなコンピュータ構成により、図1に示す本実施の形態に係る統計的形状モデル生成装置10が構成されている。なお、図1は機能ブロックとなる構成を表し、一方、図2はデバイス等の接続状態を表すものである。前記したように、機能ブロックとデバイス等とは有機的、かつ相互に関連して統計的形状モデル生成装置10を構成するため、図2と共に図1についても詳細に説明する。 With such a computer configuration, the statistical shape model generation device 10 according to the present embodiment shown in FIG. 1 is configured. FIG. 1 shows a configuration serving as a functional block, while FIG. 2 shows a connection state of devices and the like. As described above, since the functional blocks and the devices and the like constitute the statistical shape model generating apparatus 10 in an organic and mutually related manner, FIG. 1 will be described in detail together with FIG. 2.
統計的形状モデル生成装置10は、図2に示したコンピュータのハードウェア及び制御プログラムを含むソフトウェアを利用して構成される機能として、学習用受付部12と、学習部14と、統計的形状モデルデータベース16とを備えている。 The statistical shape model generation device 10 includes a learning reception unit 12, a learning unit 14, a statistical shape model, and a function configured using software including a computer hardware and a control program illustrated in FIG. A database 16 is provided.
学習用受付部12は、特定の物体である被写体を表す複数の画像であって、被写体の領域が予め求められた複数の画像を学習データとして受け付ける。本実施の形態で被写体となる特定の物体は膵臓であるため、学習用受付部12は、膵臓の領域が予め求められた複数の画像を受け付ける。 The learning accepting unit 12 accepts, as learning data, a plurality of images representing a subject which is a specific object, in which a subject area is obtained in advance. Since the specific object to be the subject in the present embodiment is the pancreas, the learning reception unit 12 receives a plurality of images in which the region of the pancreas is obtained in advance.
学習部14は、学習用受付部12によって受け付けた膵臓の領域が予め求められた複数の画像に基づいて、固有ベクトルと固有値とを算出する。例えば、学習部14は、膵臓の領域が予め求められた複数の画像に基づいて、主成分分析によって、固有ベクトルと固有値とを算出する。固有ベクトルが算出されることにより、固有ベクトルを基底とする固有空間が生成される。また、学習部14は、算出された固有ベクトルと固有値とを統計的形状モデルデータベース16に格納する。 The learning unit 14 calculates an eigenvector and an eigenvalue based on a plurality of images in which the pancreas region received by the learning reception unit 12 is obtained in advance. For example, the learning unit 14 calculates an eigenvector and an eigenvalue by principal component analysis based on a plurality of images in which the region of the pancreas is obtained in advance. By calculating the eigenvector, an eigenspace based on the eigenvector is generated. The learning unit 14 stores the calculated eigenvectors and eigenvalues in the statistical shape model database 16.
ここで、図3に、肝臓の形状の一例を示す(なお、図3に示す例は肝臓であるが、問題の本質は膵臓の場合と同様である。)。図3に示すように、肝臓の形状は様々である。そのため、実施の形態では、特定の物体の形状を表す統計的形状モデルを用いて、形状のばらつきを少数の形状パラメータにより表現する。また、図4に、主成分分析によって得られた第1主成分に対応する肝臓の形状の固有ベクトルを示す。図4に示すように、固有ベクトルに対応する形状パラメータαの値に応じて、表わされる肝臓の形状が変化することがわかる。なお、σは、統計的ばらつきを表すパラメータである。 Here, an example of the shape of the liver is shown in FIG. 3 (note that the example shown in FIG. 3 is a liver, but the nature of the problem is the same as that of the pancreas). As shown in FIG. 3, the shape of the liver varies. Therefore, in the embodiment, a variation in shape is represented by a small number of shape parameters using a statistical shape model representing the shape of a specific object. FIG. 4 shows an eigenvector of the shape of the liver corresponding to the first principal component obtained by the principal component analysis. As shown in FIG. 4, it can be seen that the shape of the represented liver changes according to the value of the shape parameter α corresponding to the eigenvector. Here, σ is a parameter representing statistical variation.
図5に、統計的形状モデルの一例を示す。図5に示すように、統計的形状モデルとしては様々なものが存在するが、本実施の形態では、統計的形状モデルとして、LSDM(例えば、参考文献(Cremers, D., Rousson, M., Deriche, R., “A review of statistical approaches to level set segmentation: integrating color, texture, motion and shape.”, International journal of computer vision 72 (2), 2007, p.195‐215.)を参照。)を用いる。 FIG. 5 shows an example of a statistical shape model. As shown in FIG. 5, there are various statistical shape models, but in the present embodiment, LSDM (for example, a reference (Cremers, D., Rousson, M., See Deriche, R., “A review of statistical approaches to level set segmentation: integrating color, texture, motion and shape.”, International journal of computer vision 72 (2), 2007, p.195-215.) Is used.
LSDMは、代表的な統計的形状モデルであり、学習データの図形ラベルに対する符号付の距離画像を作成し、符号付の距離画像に対して線形の統計解析(例えば、PCAやICAなど)を行うことで任意形状を以下の式(1)で表現する方法である。 The LSDM is a representative statistical shape model, creates a signed distance image for a graphic label of learning data, and performs linear statistical analysis (for example, PCA or ICA) on the signed distance image. This is a method of expressing an arbitrary shape by the following equation (1).
ここで、上記式(1)は画素p(∈P)に関する式である。なお、Pは画素集合を表し、φp(α)は、画素pに関する任意形状に対応するレベルセット関数を表し、μpは、画素pに関する平均のレベルセット関数を表し、λiはi番目の固有値を表し、{u1 p,...,ud p}はd個の固有ベクトルにおける画素pに対する成分を表し、α=[α1,... ,αd]Tは、固有空間上のRα={r∈Rd|||r||∞≦w}の領域(通常は±3σの範囲)で定義される形状パラメータを表す。また、wは、正の定数を示す。デジタル画像における形状は、パラメータαを以下の式(2)に示す関数gで写像することによって得られる。なお、実施の形態では、固有空間の次元d=2である場合を例に説明する。Here, the above equation (1) is an equation relating to the pixel p (∈P). Note that P represents a pixel set, φ p (α) represents a level set function corresponding to an arbitrary shape relating to the pixel p, μ p represents an average level set function relating to the pixel p, and λ i represents the i-th represents the eigenvalues, {u 1 p, ..., u d p} represents the component for the pixel p at the d eigenvectors, α = [α 1, ... , α d] T , the eigenspace on of R α = {r∈R d ||| r || ∞ ≦ w} of area (usually in the range of ± 3 [sigma]) representing the shape parameters defined in. W indicates a positive constant. The shape in the digital image is obtained by mapping the parameter α with a function g shown in the following equation (2). In the embodiment, a case where the dimension of the eigenspace is d = 2 will be described as an example.
ここで、L={0,1}はラベル集合で、0は背景であることを示し、1は図形であることを示す。yは画像における1つの形状を表す。また、LSDMの場合、以下の式(3)に示すように、パラメータαを形状に写像する関数gは、Heaviside function H(・)である。 Here, L = {0, 1} is a label set, 0 indicates a background, and 1 indicates a figure. y represents one shape in the image. In the case of LSDM, a function g that maps the parameter α to a shape is a Heaviside function H (•) as shown in the following equation (3).
図6に、d=2である場合のLSDMの固有空間(固有ベクトルが張る空間)とデジタル画像における形状集合Sとの関係を示す。関数gは固有空間からデジタル画像の形状集合S上への写像である。また、形状集合S⊂L|P|はgの像である。図6に示すように、固有空間では形状ラベルの集合は多角形を構成する。各形状y∈Sは固有空間内に原像を持つ。FIG. 6 shows the relationship between the LSDM eigenspace (space with eigenvectors) and the shape set S in the digital image when d = 2. The function g is a mapping from the eigenspace onto the shape set S of the digital image. The shape set S 集合 L | P | is an image of g. As shown in FIG. 6, a set of shape labels forms a polygon in the eigenspace. Each shape y∈S has an original image in the eigenspace.
統計的形状モデルデータベース16には、学習部14によって算出された固有ベクトルと固有値とが格納される。 The statistical shape model database 16 stores the eigenvectors and eigenvalues calculated by the learning unit 14.
<画像処理装置100の構成>
本実施の形態に係る画像処理装置100は、図7に示される機能ブロックで表すことができる。また、これらの機能ブロックは、上記図2に示されるコンピュータのハードウェア構成により実現することができる。<Configuration of Image Processing Apparatus 100>
Image processing apparatus 100 according to the present embodiment can be represented by the functional blocks shown in FIG. Further, these functional blocks can be realized by the hardware configuration of the computer shown in FIG.
画像処理装置100は、上記図2に示したコンピュータのハードウェア及び制御プログラムを含むソフトウェアを利用して構成される機能として、図7に示すように、受付部102と、演算部104と、出力部120とを備えている。 As shown in FIG. 7, the image processing apparatus 100 includes, as functions shown in FIG. 7, a receiving unit 102, an arithmetic unit 104, and an output unit. A part 120.
受付部102は、入力画像として、3時相の3次元腹部CT画像を受け付ける。なお、受付部102は、早期相画像、門脈相画像、及び晩期相画像を、3時相の3次元腹部CT画像として受け付ける。 The receiving unit 102 receives a three-dimensional three-dimensional abdominal CT image as an input image. The receiving unit 102 receives the early phase image, the portal phase image, and the late phase image as a three-dimensional three-dimensional abdominal CT image.
演算部104は、受付部102によって受け付けた3次元腹部CT画像から、膵臓の領域を抽出する。演算部104は、統計的形状モデルデータベース106と、画像処理部108とを備えている。 The calculation unit 104 extracts a pancreas region from the three-dimensional abdominal CT image received by the reception unit 102. The calculation unit 104 includes a statistical shape model database 106 and an image processing unit 108.
統計的形状モデルデータベース106には、統計的形状モデルデータベース16と同じ固有ベクトルと固有値とが格納されている。 The statistical shape model database 106 stores the same eigenvectors and eigenvalues as the statistical shape model database 16.
画像処理部108は、画像間位置合わせ部110と、空間的標準化部112と、セグメンテーション部114とを備えている。 The image processing unit 108 includes an inter-image registration unit 110, a spatial standardization unit 112, and a segmentation unit 114.
画像間位置合わせ部110は、受付部102によって受け付けた、早期相画像と、門脈相画像と、晩期相画像との間の位置合わせを行い(例えば、Shimizu, A., Kimoto, T., Kobatake, H., Nawano, S., Shinozaki, K., “Automated pancreas segmentation from three-dimensional contrast-enhanced computed tomography.”, International journal of computer assisted radiology and surgery 5 (1), 2010, p.85‐98.を参照。)、位置合わせ画像を生成する。 The inter-image registration unit 110 performs registration between the early phase image, the portal vein phase image, and the late phase image received by the reception unit 102 (for example, Shimizu, A., Kimoto, T., Kobatake, H., Nawano, S., Shinozaki, K., “Automated pancreas segmentation from three-dimensional contrast-enhanced computed tomography.”, International journal of computer assisted radiology and surgery 5 (1), 2010, p.85- 98.), and generate an alignment image.
空間的標準化部112は、画像間位置合わせ部110によって生成された位置合わせ画像に基づいて、例えば予め定められた非線形関数を用いて、位置合わせ画像を標準となる画像に合わせるように空間的標準化を行い(例えば、Shimizu, A., Kimoto, T., Kobatake, H., Nawano, S., Shinozaki, K., “Automated pancreas segmentation from three-dimensional contrast-enhanced computed tomography.”, International journal of computer assisted radiology and surgery 5 (1), 2010, p.85‐98.を参照。)、空間的標準化画像を生成する。 The spatial standardization unit 112 performs spatial standardization based on the alignment image generated by the inter-image alignment unit 110, for example, using a predetermined nonlinear function so as to match the alignment image with a standard image. (For example, Shimizu, A., Kimoto, T., Kobatake, H., Nawano, S., Shinozaki, K., “Automated pancreas segmentation from three-dimensional contrast-enhanced computed tomography.”, International journal of computer assisted radiology and surgery 5 (1), 2010, p.85-98.), and generate a spatially standardized image.
セグメンテーション部114は、統計的形状モデルデータベース106に格納された予め計算された固有ベクトルを基底とする固有空間において、受付部102によって受け付けた入力画像に基づいて、目的関数を最適化するように、空間的標準化部112によって生成された空間的標準化画像が表す被写体の形状を表す形状パラメータを推定し、空間的標準化画像が表す被写体の領域を求める。具体的には、セグメンテーション部114は、形状パラメータを推定すると同時に、被写体の領域を抽出する。 In the eigenspace based on the pre-computed eigenvectors stored in the statistical shape model database 106, the segmentation unit 114 optimizes the objective function based on the input image received by the reception unit 102. The shape parameter indicating the shape of the subject represented by the spatially standardized image generated by the spatial standardization unit 112 is estimated, and the area of the subject represented by the spatially standardized image is obtained. Specifically, the segmentation unit 114 extracts a region of the subject while estimating the shape parameters.
なお、上記図6に示すように、固有ベクトルを基底とする固有空間上の点は、統計的形状モデルの形状パラメータを示す。また、目的関数は、固有空間上の点が示す形状パラメータが表す特定の物体の形状の尤もらしさと入力画像中の隣接画素間の画素値の差とに応じた値を表すように、予め定められている。 As shown in FIG. 6, a point on the eigenspace based on the eigenvector indicates a shape parameter of the statistical shape model. Further, the objective function is determined in advance so as to represent a value corresponding to the likelihood of the shape of the specific object represented by the shape parameter indicated by the point on the eigenspace and the difference in pixel value between adjacent pixels in the input image. Have been.
一般に、上述した形状集合L|P|のサイズ(=デジタル画像上の形状数)は膨大であり、従来のアルゴリズムでは扱うことができなかった。本実施の形態におけるアルゴリズムでは、目的関数の最適化を固有空間において実行する点が従来とは異なり、これによって膨大な数の形状を効率的に扱うことが可能である。In general, the size (= the number of shapes on a digital image) of the shape set L | P | is enormous and cannot be handled by the conventional algorithm. The algorithm in the present embodiment is different from the conventional one in that the optimization of the objective function is performed in the eigenspace, whereby an enormous number of shapes can be efficiently handled.
セグメンテーション部114は、固有空間上の凸多角形を繰り返し分割して最適な形状パラメータを示す点を含む凸多角形を探索する探索アルゴリズムに従って、以下の式(4)の目的関数を最適化するように、入力画像が表す被写体の形状を表す形状パラメータを推定する。なお、探索される固有空間上の凸多角形は、最適な形状パラメータを示す点が表す特定の物体の形状を表している。
また、セグメンテーション部114は、探索アルゴリズムにおいて、凸多角形の各頂点を調べることで、その凸多角形に含まれる形状集合に対する目的関数の下界を計算することにより、最適な形状パラメータを示す点を含む固有空間上の凸多角形を探索する。
なお、本実施の形態では、固有空間が2次元である場合を例に説明するため、セグメンテーション部114は固有空間上の凸多角形を探索するが、固有空間が3次元以上の場合には、セグメンテーション部114は固有空間上の凸多胞体を探索する。また、探索時に固有空間を分割する直線は、固有空間が3次元以上の場合には超平面となる。The segmentation unit 114 optimizes the objective function of the following equation (4) according to a search algorithm for repeatedly dividing a convex polygon on the eigenspace and searching for a convex polygon including a point indicating an optimal shape parameter. Next, a shape parameter representing the shape of the subject represented by the input image is estimated. The convex polygon on the eigenspace to be searched represents the shape of a specific object represented by a point indicating the optimal shape parameter.
Further, the segmentation unit 114 searches the vertices of the convex polygon in the search algorithm, calculates the lower bound of the objective function for the set of shapes included in the convex polygon, and finds a point indicating the optimal shape parameter. Search for a convex polygon on the included eigenspace.
Note that, in the present embodiment, the segmentation unit 114 searches for a convex polygon on the eigenspace in order to describe a case where the eigenspace is two-dimensional, as an example. The segmentation unit 114 searches for a convex polyhedron in the eigenspace. In addition, a straight line that divides the eigenspace at the time of search becomes a hyperplane when the eigenspace is three-dimensional or more.
本実施の形態における重要なアイディアは、上記図6に示したように、デジタル画像上の形状集合が、固有空間上の凸多角形に対応する事実を利用することにある。本発明の実施の形態では、当該事実を利用することにより、多数の形状を一度に効率的に扱うことが可能である。 An important idea in the present embodiment is to use the fact that a shape set on a digital image corresponds to a convex polygon on an eigenspace as shown in FIG. In the embodiment of the present invention, by utilizing this fact, it is possible to efficiently handle many shapes at once.
セグメンテーション部114は、デジタル画像上の複数の形状集合が、固有空間上の凸多角形に対応する事実に基づいて、予め設定した目的関数を最小化する形状に対応する、固有空間上の凸多角形を探索する。そして、セグメンテーション部114は、探索された固有空間上の凸多角形に基づいて、目的関数を最小化する形状を形状集合Sから見つけ出し、見つけ出された形状を事前知識として当該形状を被写体の領域として抽出する。本実施の形態では、グラフカット(上記Boykov, Y., Veksler, O., Zabih, R.,“Fast approximate energy minimization via graph cuts.”, Pattern Analysis and Machine Intelligence, IEEE Transactions on 23 (11), 2001, p.1222‐1239.を参照。)で良く用いられるエネルギー関数を、予め設定した目的関数として用いる場合を例に説明する。本実施の形態で用いる目的関数を、以下の式(4)に示す。 Based on the fact that a plurality of shape sets on the digital image correspond to convex polygons on the eigenspace, the segmentation unit 114 corresponds to a convex shape on the eigenspace corresponding to a shape that minimizes a preset objective function. Search for a polygon. Then, the segmentation unit 114 finds a shape that minimizes the objective function from the shape set S based on the searched convex polygon on the eigenspace, and uses the found shape as prior knowledge to define the shape as the area of the subject. Extract as In this embodiment, the graph cut (the above-mentioned Boykov, Y., Veksler, O., Zabih, R., “Fast approximate energy minimization via graph cuts.”, Pattern Analysis and Machine Intelligence, IEEE Transactions on 23 (11), 2001, p. 1222-1239.). An example in which an energy function often used as an objective function set in advance will be described. The objective function used in the present embodiment is shown in the following equation (4).
ここで、Nは、隣接する画素のペアの集合である。xpは画素pのラベルを表し、背景であれば0の値をとり、図形であれば1の値をとる。また、Ipは、画素pの画素値を表し、ypは、画素pにおける推定形状の値を表し、推定形状が背景であれば0の値をとり、図形であれば1の値をとる。Here, N is a set of pairs of adjacent pixels. x p represents the label of the pixel p, a value of 0 if the background, taking a value of 1 if figure. Also, I p denotes a pixel value of the pixel p, y p represents the value of the estimated shape of the pixel p, the estimated shape takes a value of 0 if the background, taking a value of 1 if the figure .
また、上記式(4)におけるFp(Ip,yp)と、Bp(Ip,yp)とは、以下の式(5)、(6)で定義される。ここで、以下の式(5)は、画素pが図形として割り当てられるコストを表し、以下の式(6)は、画素pが背景として割り当てられるコストを表す。 Further, F p (I p, y p) in the above formula (4) and, B p (I p, y p) and the following equation (5) is defined by (6). Here, the following equation (5) represents the cost of assigning the pixel p as a figure, and the following equation (6) represents the cost of assigning the pixel p as a background.
Ipは画素pの画素値を表し、Iqは画素qの画素値を表す。また、λ1、λ2は予め定められた正の定数である。また、上記式(4)におけるPpq(Ip,Iq)は、以下の式(7)で定義され、互いに隣接する画素pと画素qの間の画素値の差を評価している。 Ip represents the pixel value of pixel p, and Iq represents the pixel value of pixel q. Λ 1 and λ 2 are predetermined positive constants. In addition, P pq (I p , I q ) in the above equation (4) is defined by the following equation (7), and evaluates a difference in pixel value between the adjacent pixels p and q.
本実施の形態では、上記式(4)に示す目的関数の最適化手法として、一般的な最適化アルゴリズムの一つであるbranch and bound探索アルゴリズムを用いる。本実施の形態では、形状集合に対する目的関数の下界を効率的に計算する方法を提案する。形状集合に対する目的関数の下界を効率的に計算するための重要なアイディアは、凸多角形の頂点のみを調べれば、その凸多角形に含まれるすべての形状集合に対する目的関数の下界が分かる点である。 In the present embodiment, a branch and bound search algorithm, which is one of general optimization algorithms, is used as an objective function optimization method represented by the above equation (4). This embodiment proposes a method for efficiently calculating the lower bound of an objective function for a shape set. An important idea for efficiently calculating the lower bound of the objective function for a set of shapes is that examining only the vertices of a convex polygon gives the lower bound of the objective function for all the sets of shapes included in the convex polygon. is there.
上記式(4)を、形状yの代わりに、関数g(α)を用いて表すと、以下の式(8)で表すことができる。 If the above equation (4) is expressed using a function g (α) instead of the shape y, it can be expressed by the following equation (8).
[branching処理]
branch and bound探索アルゴリズムのbranching処理では、親ノードH0(∈Rα)が与えられると、親ノードH0は子ノードH1及びH2へ分解され、上記式(8)に示した目的関数は、以下の式(9)で表わされる。[Branching process]
In the branching process of the branch and bound search algorithm, when a parent node H 0 (∈R α ) is given, the parent node H 0 is decomposed into child nodes H 1 and H 2 , and the objective function shown in the above equation (8) Is represented by the following equation (9).
ここで、branching処理における分割処理は、様々な分割方法が可能である。本実施の形態では、上記式(1)に示したφk(α)の符号によって、以下の式(10)及び(11)に示すように、親ノードH0は子ノードH1及びH2へ分割される。なお、画素kは、集合Qからサンプリングによって選択された画素を表す。Here, for the dividing process in the branching process, various dividing methods are possible. In the present embodiment, as shown in the following equations (10) and (11), the parent node H 0 becomes the child nodes H 1 and H 2 according to the sign of φ k (α) shown in the above equation (1). Divided into Note that the pixel k represents a pixel selected by sampling from the set Q.
ここで、集合Qは、以下の式(12)に示すように、固有空間内のノードH0を切る直線φk(α)=0であるような画素kの集合である。なお、固有空間が3次元以上である場合、集合Qは、固有空間内のノードH0を切る超平面φk(α)=0であるような画素kの集合となる。Here, the set Q is a set of pixels k such that a straight line φ k (α) = 0 that cuts the node H 0 in the eigenspace as shown in the following equation (12). When the eigenspace is three-dimensional or more, the set Q is a set of pixels k such that a hyperplane φ k (α) = 0 that cuts the node H 0 in the eigenspace.
また、V0は、ノードH0が表す凸多角形の頂点集合を表す。V0は、H0のエッジに対応するφe(α)の連立方程式を解析的に解くことで得ることができる。V 0 represents a vertex set of the convex polygon represented by the node H 0 . V 0 can be obtained by analytically solving simultaneous equations of φ e (α) corresponding to the edge of H 0 .
[bounding処理]
branch and bound探索アルゴリズムのbranching処理においては、下界L(Hi)(i={1,2})は、以下の式(13)〜(14)で与えられる。[Bounding process]
In the branching process of the branch and bound search algorithm, the lower bound L (H i ) (i = {1, 2}) is given by the following equations (13) to (14).
上記式(13)は、最小化問題のためのJensenの不等式である。上記式(14)における「max」と「min」との変化は、上記式(5)のマイナス符号に基づいている。また、α∈Hiにおけるφ(α)の最大値と最小値とは、線形計画法の基礎理論に基づき、Hiの頂点集合Viから取得され、以下の式(15)及び(16)で表すことができる。Equation (13) above is Jensen's inequality for the minimization problem. The change between “max” and “min” in equation (14) is based on the minus sign in equation (5). Further, the maximum value and the minimum value of phi (alpha) in Arufa∈H i, based on the basic theory of linear programming, is obtained from the vertex set V i of H i, the following equation (15) and (16) Can be represented by
また、セグメンテーション部114は、分割後の2つの凸多角形の面積が対応するように(面積がほぼ等しくなるように)、固有空間上の凸多角形を分割することを繰り返すことにより、最適な形状パラメータを示す点を含む固有空間上の凸多角形を探索する。 In addition, the segmentation unit 114 repeatedly divides the convex polygon on the eigenspace so that the areas of the two convex polygons after division correspond to each other (so that the areas are substantially equal), so that the optimal polygon is repeated. Search for a convex polygon on the eigenspace containing points indicating shape parameters.
セグメンテーション部114によるセグメンテーションアルゴリズムの疑似コードを表1に示す。まず、対象画像Iが与えられると、セグメンテーション部114は、固有空間全体Rα={r∈Rd|||r||∞≦w}を親ノードH0として処理を開始する。すなわち、セグメンテーション部114は、Rαをルートとして処理を開始し、デジタル画像から関数Select(Q)を用いて選択した画素kによって親ノードH0を分割する。分割後の子ノードH1とH2とは、Queueに入れられる。ただし、Queueに入れられたノードは、下界の値が昇順になるように並べられているとする。次に、セグメンテーション部114は、Queueから、最も低い下界のノードを選択して親ノードH0とし、上記の処理をQが空集合になるまで反復する。反復終了時の結果から、最適な形状パラメータが得られ、最適な形状が求まる。そして、セグメンテーション部114は、最適な形状を用いて、最適なセグメンテーション結果を得る。Table 1 shows the pseudo code of the segmentation algorithm by the segmentation unit 114. First, when the target image I is given, the segmentation unit 114, the entire eigenspace R α = {r∈R d ||| r || ∞ ≦ w} to start processing as a parent node H 0. That is, the segmentation unit 114 starts processing the R alpha as root, divides the parent node H 0 by the pixel k selected by using the function the Select (Q) from the digital image. Child nodes H 1 and H 2 after the division are placed in Queue. However, it is assumed that the nodes put in Queue are arranged such that the lower bound values are in ascending order. Next, the segmentation unit 114, a Queue, and the parent node H 0 Select the lowest lower bound of the node is repeated until the processing Q is an empty set. From the result at the end of the repetition, an optimal shape parameter is obtained and an optimal shape is obtained. Then, the segmentation unit 114 obtains an optimal segmentation result using the optimal shape.
なお、H0は、固有空間における凸多角形を表す。なお、H0の初期値はRαに対応する。また、H1、H2は、H0を分割してできた二つの凸多角形を表す。Qは、H0を分割可能なデジタル画像上の画素集合である。Queueは、凸多角形ノードが入っているキュー(Branch and bound探索アルゴリズムで利用)を表す。Select(Q)は、分割のための画素kを選択する関数を表す。Note that H 0 represents a convex polygon in the eigenspace. The initial value of H 0 corresponds to R alpha. H 1 and H 2 represent two convex polygons formed by dividing H 0 . Q is a set of pixels on the dividable digital image H 0. Queue represents a queue containing convex polygon nodes (used in the Branch and bound search algorithm). Select (Q) represents a function for selecting a pixel k for division.
また、効率的にアルゴリズムを動作させるためのSelect(Q)に関する新しい提案も、本実施の形態における重要な特徴である。ここで、Select(Q)は画像から凸多角形を分割するための画素を選択する関数であるが、探索を効率的に行なうために、分割後の子ノードH1,H2に対応する凸多角形の面積がほぼ等しくなるような分割法を本実施の形態で提案し、その効率を実際の画像で検証した。Further, a new proposal on Select (Q) for efficiently operating the algorithm is also an important feature in the present embodiment. Here, Select (Q) is a function for selecting a pixel for dividing a convex polygon from an image. In order to perform a search efficiently, convex (Q) corresponding to child nodes H 1 and H 2 after division is selected. In the present embodiment, a division method in which the areas of polygons are almost equal was proposed, and the efficiency was verified with an actual image.
本実施の形態では、以下の式(17)に示す関数を、Select(Q)として用いる。 In the present embodiment, the function shown in the following equation (17) is used as Select (Q).
ここで、 here,
である。 It is.
図8に、上記式(17)の関数の作用の概念図を示す。図8に示すように、上記式(17)の作用は、直感的には、H0を面積が等しい2つの領域(H1,H2)に分割するような画素pを選択する。FIG. 8 shows a conceptual diagram of the operation of the function of the above equation (17). As shown in FIG. 8, the action of the above formula (17), Intuitively, selects pixels p that divides the H 0 into two regions having the same area as that of (H 1, H 2).
セグメンテーション部114は、探索された特定の物体の最適な形状パラメータが表す形状を事前知識として、目的関数を最適化するように、特定の物体の形状である被写体の領域を、空間的標準化部112によって生成された空間的標準化画像から抽出する。 The segmentation unit 114 uses the shape represented by the optimum shape parameter of the searched specific object as a priori knowledge, and converts the region of the subject having the shape of the specific object into the spatial standardization unit 112 so as to optimize the objective function. From the spatially standardized image generated by.
出力部120は、セグメンテーション部114によって抽出された被写体の領域を結果として出力する。 The output unit 120 outputs the region of the subject extracted by the segmentation unit 114 as a result.
<統計的形状モデル生成装置10の作用>
次に、統計的形状モデル生成装置10の作用を説明する。<Operation of Statistical Shape Model Generation Apparatus 10>
Next, the operation of the statistical shape model generation device 10 will be described.
図9に示す統計的形状モデル生成装置10の処理は、上記図2におけるCPU21によりHDD24等に記憶されたプログラムに基づき行なわれるものである。統計的形状モデル生成装置10に、膵臓の領域が予め求められた複数の画像が入力されると、図9に示す学習処理ルーチンが実行される。 The processing of the statistical shape model generation device 10 shown in FIG. 9 is performed by the CPU 21 shown in FIG. 2 based on a program stored in the HDD 24 or the like. When a plurality of images in which pancreatic regions are obtained in advance are input to the statistical shape model generation device 10, a learning processing routine shown in FIG. 9 is executed.
まず、ステップS100において、学習用受付部12は、膵臓の領域が予め求められた複数の画像を受け付ける。 First, in step S100, the learning accepting unit 12 accepts a plurality of images in which pancreatic regions have been obtained in advance.
次に、ステップS102において、学習部14は、上記ステップS100で受け付けた膵臓の領域が予め求められた複数の画像に基づいて、主成分分析によって、固有ベクトルと固有値とを算出する。 Next, in step S102, the learning unit 14 calculates an eigenvector and an eigenvalue by principal component analysis based on a plurality of images in which the pancreas region received in step S100 is obtained in advance.
ステップS104において、学習部14は、上記ステップS102で算出された固有ベクトルと固有値とを統計的形状モデルデータベース16に格納し、学習処理ルーチンを終了する。 In step S104, the learning unit 14 stores the eigenvectors and eigenvalues calculated in step S102 in the statistical shape model database 16, and ends the learning processing routine.
<画像処理装置100の作用>
次に、画像処理装置100の作用を説明する。<Operation of Image Processing Apparatus 100>
Next, the operation of the image processing apparatus 100 will be described.
図10に示す画像処理装置100の処理は、上記図2におけるCPU21によりHDD24等に記憶されたプログラムに基づき行なわれるものである。まず、統計的形状モデル生成装置10の統計的形状モデルデータベース16に格納されている固有ベクトルと固有値とが、画像処理装置100に入力されると、統計的形状モデルデータベース106に格納される。そして、画像処理装置100に、3時相の3次元腹部CT画像が入力画像として入力されると、上記図10に示す画像処理ルーチンが実行される。 The processing of the image processing apparatus 100 shown in FIG. 10 is performed by the CPU 21 shown in FIG. 2 based on a program stored in the HDD 24 or the like. First, when the eigenvectors and eigenvalues stored in the statistical shape model database 16 of the statistical shape model generation device 10 are input to the image processing device 100, they are stored in the statistical shape model database 106. Then, when a three-dimensional three-dimensional abdominal CT image is input to the image processing apparatus 100 as an input image, the image processing routine shown in FIG. 10 is executed.
ステップS200において、受付部102は、入力画像として、3時相の3次元腹部CT画像を受け付ける。 In step S200, the receiving unit 102 receives a three-dimensional three-dimensional abdominal CT image as an input image.
ステップS202において、画像間位置合わせ部110は、上記ステップS200で3時相の3次元腹部CT画像として受け付けた、早期相画像と、門脈相画像と、晩期相画像との間の位置合わせを行い、位置合わせ画像を生成する。 In step S202, the inter-image positioning unit 110 performs positioning between the early phase image, the portal vein phase image, and the late phase image received as the three-phase three-dimensional abdominal CT image in step S200. Then, an alignment image is generated.
ステップS204において、空間的標準化部112は、上記ステップS202で生成された位置合わせ画像に基づいて、例えば非線形関数を用いて、空間的標準化画像を生成する。 In step S204, the spatial standardization unit 112 generates a spatial standardized image based on the alignment image generated in step S202, for example, using a nonlinear function.
ステップS206において、セグメンテーション部114は、統計的形状モデルデータベース106に格納された予め計算された固有ベクトルを基底とする固有空間において、上記ステップS200で受け付けた入力画像に基づいて、目的関数を最適化するように、上記ステップS204で生成された空間的標準化画像が表す被写体の形状を表す形状パラメータを推定し、空間的標準化画像が表す被写体の形状を求める。ステップS206は、図11に示すセグメンテーション処理ルーチンによって実現される。 In step S206, the segmentation unit 114 optimizes the objective function based on the input image received in step S200 in the eigenspace based on the pre-calculated eigenvector stored in the statistical shape model database 106. As described above, the shape parameter representing the shape of the subject represented by the spatial standardized image generated in step S204 is estimated, and the shape of the subject represented by the spatial standardized image is obtained. Step S206 is realized by the segmentation processing routine shown in FIG.
<セグメンテーション処理ルーチン>
ステップS300において、セグメンテーション部114は、固有空間全体Rαを親ノードH0として設定する。<Segmentation processing routine>
In step S300, the segmentation unit 114 sets the entire eigenspace R alpha as a parent node H 0.
ステップS302において、セグメンテーション部114は、上記式(12)に示すように、親ノードH0を分割可能な画素kの集合を、集合Qとして設定する。In step S302, the segmentation unit 114, as shown in the equation (12), the set of splittable pixel k parent node H 0, is set as the set Q.
ステップS304において、セグメンテーション部114は、Queueを初期化する。 In step S304, the segmentation unit 114 initializes Queue.
ステップS306において、セグメンテーション部114は、上記ステップS300で設定した親ノードH0と、上記式(14)に従って算出される当該親ノードH0に対応する下界L(H0;I)との組み合わせを、上記ステップS304で初期化されたQueueに格納する。In step S306, the segmentation unit 114 combines the combination of the parent node H 0 set in step S300 with the lower bound L (H 0 ; I) corresponding to the parent node H 0 calculated according to the above equation (14). , Stored in the Queue initialized in step S304.
ステップS308において、セグメンテーション部114は、上記ステップS302で設定された集合Q、又は後述するステップS316で更新された集合Qから、上記式(17)に従って、画素kを選択する。 In step S308, the segmentation unit 114 selects a pixel k from the set Q set in step S302 or the set Q updated in step S316 described below according to the above equation (17).
ステップS310において、セグメンテーション部114は、上記ステップS300で設定された親ノードH0、又は後述するステップS314で更新された親ノードH0を、上記ステップS308で選択された画素kを用いて上記式(10)及び(11)に従って分割し、子ノードH1及びH2を設定する。In step S310, the segmentation unit 114, the above formula parent node H 0 is set in the step S300, the or parent node H 0 updated in step S314 to be described later, by using a pixel k selected in step S308 (10) and divided according to (11), sets the child nodes H 1 and H 2.
ステップS312において、セグメンテーション部114は、上記ステップS310で設定された子ノードH1及び上記式(14)に従って算出される当該子ノードH1に対応する下界L(H1;I)との組み合わせを、Queueに格納する。また、上記ステップS310で設定された子ノードH2及び上記式(14)に従って算出される当該子ノードH2に対応する下界L(H2;I)との組み合わせを、Queueに格納する。In step S312, the segmentation unit 114, the lower bound L corresponding to the child node H 1 calculated in accordance with the set child nodes H 1 and the formula in the above step S310 (14); a combination of (H 1 I) , Queue. Also, the lower bound L corresponding to the child node H 2 calculated in accordance with the set child nodes H 2 and the equation in step S310 (14); a combination of (H 2 I), and stored in Queue.
ステップS314において、セグメンテーション部114は、Queueに格納されたノードのうち、最も低い下界のノードを選択し、親ノードH0を、選択されたノードに更新する。In step S314, the segmentation unit 114, among the stored node in Queue, select the lowest lower bound of the node, the parent node H 0, updates the selected node.
ステップS316において、セグメンテーション部114は、上記ステップS314で更新された親ノードH0に基づいて、上記式(12)に示すように、集合Qを、親ノードH0を分割可能な画素kの集合に更新する。なお、親ノードH0を分割可能な画素kが存在しない場合には、集合Qが空集合に更新される。In step S316, the segmentation unit 114, based on the parent node H 0 updated in the step S314, the as shown in the equation (12), the set Q, the set of splittable pixel k of the parent node H 0 Update to In the case where there is no pixel k dividable parent node H 0 is the set Q is updated to an empty set.
ステップS318において、セグメンテーション部114は、集合Qが空集合であるか否かを判定する。集合Qが空集合である場合には、ステップS308へ戻る。一方、集合Qが空集合でない場合には、ステップS320へ進む。 In step S318, the segmentation unit 114 determines whether the set Q is an empty set. If the set Q is an empty set, the process returns to step S308. On the other hand, if the set Q is not an empty set, the process proceeds to step S320.
ステップS320において、セグメンテーション部114は、上記ステップS314で最終的に更新された親ノードH0に含まれる形状パラメータα*を選択し、選択された形状パラメータα*を関数gに代入し、特定の物体の形状y*を決定する。In step S320, the segmentation unit 114 is finally updated by selecting the shape parameter alpha * in the parent node H 0 above step S314, the substituted a selected shape parameter alpha * the function g, the specific Determine the shape y * of the object.
ステップS322において、セグメンテーション部114は、上記ステップS320で決定された特定の物体の形状y*及び目的関数に基づいて、目的関数を最適化するように、特定の物体の形状である被写体の領域x*を、上記ステップS204で生成された空間的標準化画像から抽出する。In step S322, based on the specific object shape y * and the objective function determined in step S320, the segmentation unit 114 subjects the object region x having the specific object shape to optimize the objective function. * Is extracted from the spatially standardized image generated in step S204.
ステップS324において、セグメンテーション部114は、上記ステップS322で抽出された被写体の領域x*を結果として出力する。In step S324, the segmentation unit 114 outputs the region x * of the subject extracted in step S322 as a result.
次に画像処理ルーチンに戻り、出力部120は、ステップS208において、上記ステップS206で抽出された被写体の領域x*を結果として出力して、画像処理ルーチンを終了する。Next, returning to the image processing routine, in step S208, the output unit 120 outputs as a result the region x * of the subject extracted in step S206, and ends the image processing routine.
以上説明したように、本実施形態の画像処理装置100では、プログラムに基づくコンピュータ処理による機能として、(A)被写体の領域が予め求められた複数の画像に基づいて予め計算された固有ベクトルを基底とする固有空間であって、(B)固有空間上の点が、特定の物体の形状を表す統計的形状モデルの形状パラメータを示す固有空間において、固有空間上の点が示す形状パラメータが表す特定の物体の形状の尤もらしさと入力画像中の隣接画素間の画素値の差とに応じた値を表す予め定められた目的関数を最適化するように、入力画像が表す被写体の形状を表す形状パラメータを推定する。推定された形状パラメータが表す特定の物体の形状を事前知識として、被写体の領域を、入力画像から抽出する。これにより、計算量の増大を抑制して、被写体の領域を精度よく抽出することができる。 As described above, in the image processing apparatus 100 of the present embodiment, as a function of the computer processing based on the program, (A) the eigenvector calculated in advance based on a plurality of images in which the region of the subject is obtained in advance is used as the basis. (B) In the eigenspace indicating the shape parameter of the statistical shape model representing the shape of the specific object, the point on the eigenspace indicates the specific parameter represented by the shape parameter indicated by the point on the eigenspace. Shape parameters representing the shape of the subject represented by the input image so as to optimize a predetermined objective function representing a value corresponding to the likelihood of the shape of the object and the difference between pixel values between adjacent pixels in the input image Is estimated. The region of the subject is extracted from the input image using the shape of the specific object represented by the estimated shape parameter as prior knowledge. As a result, it is possible to suppress an increase in the amount of calculation and accurately extract the region of the subject.
また、本実施形態の画像処理装置100は、被写体の領域を高速に抽出することができる。 Further, the image processing apparatus 100 according to the present embodiment can extract a region of a subject at high speed.
また、提案するアルゴリズムと従来のアルゴリズムとの重要な違いの一つに、前処理が不要である点が挙げられる。従来は、あらかじめ形状を用意し、クラスタリングをする前処理が必要であり、前処理にかなりの時間とメモリを必要としていたが、本実施の形態における提案アルゴリズムでは前処理は一切不要である。これにより、コストの削減を実現した。 One of the important differences between the proposed algorithm and the conventional algorithm is that no preprocessing is required. Conventionally, a pre-process for preparing a shape in advance and performing clustering is required, and a considerable amount of time and memory are required for the pre-process. However, the pre-process is not required at all in the proposed algorithm in the present embodiment. As a result, cost reduction was realized.
また、本実施形態の形態では、140例の3次元腹部CT画像から膵臓をセグメンテーションする実験を行って処理精度や処理時間を検証し、以下の結果が得られた。 Further, in the embodiment of the present invention, an experiment for segmenting the pancreas from 140 three-dimensional abdominal CT images was performed to verify processing accuracy and processing time, and the following results were obtained.
1.3次元腹部CT画像からの膵臓のセグメンテーションに関して世界最高の精度が得られた。
2.提案アルゴリズムによる目的関数の最小化後の値は、従来手法と比較すると必ず小さくなることが理論的に保証されているが、そのことを実際の画像でも確認した。
3.従来手法(上記米国特許第8249349号明細書)で扱われていた形状の最大数が107であるのに比べ、提案法では109以上の数を扱うことができるようになった。また、画像も従来は2次元の画像に限られていたのに対して、画素数が二桁近く多い3次元画像も扱えるようになった。1. The world's highest accuracy was obtained for pancreas segmentation from 3D abdominal CT images.
2. It is theoretically guaranteed that the value after minimization of the objective function by the proposed algorithm is always smaller than that of the conventional method, but this was confirmed in actual images.
3. Compared to the maximum number of conventional techniques shapes covered in (the U.S. Pat. No. 8,249,349) is 10 7, the proposed method it has become possible to handle the number of more than 10 9. Further, the image is conventionally limited to a two-dimensional image, but a three-dimensional image having a number of pixels nearly two digits larger can be handled.
[第2の実施の形態]
次に、第2の実施の形態について説明する。なお、第2の実施の形態に係る統計的形状モデル生成装置及び画像処理装置の構成は、第1の実施の形態と同様の構成となるため、同一符号を付して説明を省略する。[Second embodiment]
Next, a second embodiment will be described. Note that the configurations of the statistical shape model generation device and the image processing device according to the second embodiment are the same as those of the first embodiment, and thus the same reference numerals are given and the description is omitted.
第2の実施の形態では、近似解法を用いて、被写体の形状を表す形状パラメータを推定する点が、第1の実施の形態と異なる。 The second embodiment is different from the first embodiment in that a shape parameter representing the shape of a subject is estimated using an approximate solution.
図12に、第2の実施の形態における近似解法を説明するための図を示す。図12に示すように、上記第1の実施の形態で用いた厳密解法とは異なり、第2の実施の形態では近似解法によって形状パラメータを推定する。第2の実施の形態では、図12に示すように、固有空間上に格子状にサンプリング点を設定し、サンプリング点から決定される直線を用いて凸多角形を分割することを繰り返す。なお、サンプリング点の設定は、予め定められたサンプリングパラメータkに応じて行われる。例えば、サンプリングパラメータk=4の場合、固有空間の次元ごとに24個のサンプリング点が設定される。FIG. 12 is a diagram for explaining an approximate solution method according to the second embodiment. As shown in FIG. 12, unlike the strict solution method used in the first embodiment, the shape parameter is estimated by an approximate solution in the second embodiment. In the second embodiment, as shown in FIG. 12, the sampling points are set in a grid on the eigenspace and the division of the convex polygon using a straight line determined from the sampling points is repeated. The setting of the sampling point is performed according to a predetermined sampling parameter k. For example, if the sampling parameters k = 4, 2 4 sampling points for each dimension of the eigenspace is set.
第2の実施の形態に係る画像処理装置のセグメンテーション部114は、まず、固有空間上に格子状にサンプリング点を設定する。 First, the segmentation unit 114 of the image processing apparatus according to the second embodiment sets sampling points in a lattice on the eigenspace.
そして、セグメンテーション部114は、固有空間上に設定されたサンプリング点から、分割後の2つの凸多角形の面積が対応するように、直線を設定し、設定された直線で凸多角形を分割して、最適な形状パラメータを示す点を含む凸多角形を探索することを繰り返し、入力画像が表す被写体の形状を表す形状パラメータを推定する。
具体的には、上記表1に示したセグメンテーションアルゴリズムの疑似コードにおいて、関数Select(Q)を用いて選択した画素kによって親ノードH0を分割する処理の代わりに、サンプリング点から決定される直線を用いて、親ノードH0を分割する処理を行う。Then, the segmentation unit 114 sets a straight line from the sampling points set on the eigenspace so that the areas of the two convex polygons after division correspond to each other, and divides the convex polygon by the set straight line. Then, the search for the convex polygon including the point indicating the optimal shape parameter is repeated to estimate the shape parameter representing the shape of the subject represented by the input image.
Linear Specifically, in pseudo-code segmentation algorithm shown in Table 1, instead of the process of dividing the parent node H 0 by the pixel k selected by using the function the Select (Q), which is determined from the sampling point using, performs processing for dividing a parent node H 0.
第2の実施の形態では、近似解法を用いることにより、固有空間の次元数を増加させることができ、形状パラメータを精度良く推定することができる。 In the second embodiment, by using the approximate solution method, the number of dimensions of the eigenspace can be increased, and the shape parameter can be accurately estimated.
図13に、第2の実施の形態に係る画像処理装置による140例の膵臓の認識結果を示す。図13において、破線の円は上記第1の実施の形態に係る画像処理装置による膵臓の認識結果を表し、実線の円は第2の実施の形態に係る画像処理装置による膵臓の認識結果を表す。
図13に示すように、近似解法を用いることにより、固有空間の次元数を上げることで厳密解よりも高いセグメンテーション精度が達成されていることがわかる。FIG. 13 shows the results of 140 cases of pancreas recognition performed by the image processing apparatus according to the second embodiment. In FIG. 13, a broken-line circle represents the result of pancreatic recognition by the image processing device according to the first embodiment, and a solid-line circle represents the result of pancreatic recognition by the image processing device according to the second embodiment. .
As shown in FIG. 13, it can be seen that by using the approximate solution method, a higher segmentation accuracy than the exact solution is achieved by increasing the number of dimensions of the eigenspace.
また、図14に、第2の実施の形態の近似解法、上記米国特許第8249349号明細書に記載の手法を膵臓のセグメンテーションに適用させた場合の解法、及び上記第1の実施の形態の厳密解法の計算時間の比較結果を示す。図14は、固有空間の次元dを増加させた場合の計算時間の変化を示している。なお、厳密解法については計算時間の上限を100hとした。また、近似解法についてはサンプリングパラメータの値をk=4とした。また、CPUは、3.1GHzのIntel(R)Xeon(R)を用い、前処理に1CPU、最適化に2CPUを用いた。 FIG. 14 shows an approximate solution of the second embodiment, a solution in a case where the method described in the above-mentioned US Pat. No. 8,249,349 is applied to pancreas segmentation, and an exact solution of the first embodiment. The result of comparing the calculation time of the solution is shown. FIG. 14 shows a change in calculation time when the dimension d of the eigenspace is increased. In the strict solution method, the upper limit of the calculation time was set to 100 hours. For the approximate solution, the value of the sampling parameter was set to k = 4. As the CPU, 3.1 GHz Intel (R) Xeon (R) was used, and 1 CPU was used for preprocessing and 2 CPU was used for optimization.
図14に示すように、近似解法を用いた場合には、上記第1の実施の形態の厳密解法を用いた場合及び上記米国特許第8249349号明細書に記載の手法を膵臓のセグメンテーションに適用させた場合に比べ、固有空間の次元数を上げた場合の計算コストを大幅に削減可能であることがわかる。 As shown in FIG. 14, when the approximate solution is used, the case where the exact solution of the first embodiment is used and the method described in US Pat. No. 8,249,349 are applied to the pancreas segmentation. It can be seen that the calculation cost when the number of dimensions of the eigenspace is increased can be greatly reduced as compared with the case where
なお、第2の実施の形態に係る画像処理装置の他の構成及び作用については、第1の実施の形態と同様であるため、説明を省略する。 Note that other configurations and operations of the image processing apparatus according to the second embodiment are the same as those of the first embodiment, and a description thereof will not be repeated.
以上説明したように、第2の実施の形態に係る画像処理装置によれば、固有空間上の凸多角形上に、格子状にサンプリング点を設定し、固有空間上の凸多角形に設定された直線を用いて凸多角形を分割して探索することを繰り返して、入力画像が表す被写体の形状を表す形状パラメータを推定する。推定された形状パラメータが表す特定の物体の形状を事前知識として、被写体の領域を、入力画像から抽出する。これにより、計算量の増大を抑制して、被写体の領域を精度よく抽出することができる。 As described above, according to the image processing apparatus according to the second embodiment, sampling points are set in a grid on convex polygons in the eigenspace, and are set as convex polygons in the eigenspace. By repeatedly performing the search by dividing the convex polygon using the straight line, the shape parameter representing the shape of the subject represented by the input image is estimated. The region of the subject is extracted from the input image using the shape of the specific object represented by the estimated shape parameter as prior knowledge. As a result, it is possible to suppress an increase in the amount of calculation and accurately extract the region of the subject.
また、近似解法により、次元数を増加させたより高度な統計的形状モデルが利用可能となる。 Further, a more advanced statistical shape model with an increased number of dimensions can be used by the approximate solution method.
[第3の実施の形態]
次に、第3の実施の形態について説明する。なお、第3の実施の形態に係る統計的形状モデル生成装置及び画像処理装置の構成は、第1の実施の形態と同様の構成となるため、同一符号を付して説明を省略する。[Third Embodiment]
Next, a third embodiment will be described. Since the configurations of the statistical shape model generation device and the image processing device according to the third embodiment are the same as those of the first embodiment, the same reference numerals are given and the description is omitted.
第3の実施の形態では、ランダムにサンプリング点を設定する近似解法を用いて、被写体の形状を表す形状パラメータを推定する点が、第2の実施の形態と異なる。 The third embodiment differs from the second embodiment in that a shape parameter representing the shape of a subject is estimated using an approximate solution in which sampling points are randomly set.
図15に、第3の実施の形態における近似解法を説明するための図を示す。図15に示すように、第3の実施の形態では、固有空間上にランダムにサンプリング点を設定し、サンプリング点から決定される直線を用いて凸多角形を分割することを繰り返す。 FIG. 15 is a diagram for explaining an approximate solution method according to the third embodiment. As shown in FIG. 15, in the third embodiment, the sampling points are randomly set on the eigenspace, and the process of dividing the convex polygon using a straight line determined from the sampling points is repeated.
第3の実施の形態に係る画像処理装置のセグメンテーション部114は、まず、固有空間上にランダムにサンプリング点を設定する。 The segmentation unit 114 of the image processing device according to the third embodiment first sets sampling points randomly on the eigenspace.
そして、セグメンテーション部114は、固有空間上に設定されたサンプリング点から、分割後の2つの凸多角形の面積が対応するように、直線を設定し、設定された直線で凸多角形を分割して、最適な形状パラメータを示す点を含む凸多角形を探索することを繰り返し、入力画像が表す被写体の形状を表す形状パラメータを推定する。 Then, the segmentation unit 114 sets a straight line from the sampling points set on the eigenspace so that the areas of the two convex polygons after division correspond to each other, and divides the convex polygon by the set straight line. Then, the search for the convex polygon including the point indicating the optimal shape parameter is repeated to estimate the shape parameter representing the shape of the subject represented by the input image.
具体的には、上記表1に示したセグメンテーションアルゴリズムの疑似コードにおいて、関数Select(Q)を用いて選択した画素kによって親ノードH0を分割する処理の代わりに、サンプリング点から決定される直線を用いて親ノードH0を分割する処理を行う。Linear Specifically, in pseudo-code segmentation algorithm shown in Table 1, instead of the process of dividing the parent node H 0 by the pixel k selected by using the function the Select (Q), which is determined from the sampling point It performs a process of dividing the parent node H 0 using.
なお、第3の実施の形態に係る画像処理装置の他の構成及び作用については、第1又は第2の実施の形態と同様であるため、説明を省略する。 Note that other configurations and operations of the image processing apparatus according to the third embodiment are the same as those of the first or second embodiment, and thus description thereof is omitted.
以上説明したように、第3の実施の形態に係る画像処理装置によれば、固有空間上の凸多角形上に、ランダムにサンプリング点を設定し、固有空間上の凸多角形に設定された直線を用いて凸多角形を分割して探索することを繰り返して、入力画像が表す被写体の形状を表す形状パラメータを推定する。推定された形状パラメータが表す特定の物体の形状を事前知識として、被写体の領域を、入力画像から抽出する。これにより、計算量の増大を抑制して、被写体の領域を精度よく抽出することができる。 As described above, according to the image processing apparatus according to the third embodiment, the sampling points are randomly set on the convex polygon in the eigenspace, and are set to the convex polygon in the eigenspace. The shape parameter representing the shape of the subject represented by the input image is estimated by repeating the search by dividing the convex polygon using the straight line. The region of the subject is extracted from the input image using the shape of the specific object represented by the estimated shape parameter as prior knowledge. As a result, it is possible to suppress an increase in the amount of calculation and accurately extract the region of the subject.
[第4の実施の形態]
次に、第4の実施の形態について説明する。なお、第4の実施の形態に係る統計的形状モデル生成装置及び画像処理装置の構成は、第1の実施の形態と同様の構成となるため、同一符号を付して説明を省略する。[Fourth Embodiment]
Next, a fourth embodiment will be described. Since the configurations of the statistical shape model generation device and the image processing device according to the fourth embodiment are the same as those of the first embodiment, the same reference numerals are given and the description is omitted.
第4の実施の形態では、目的関数を拡張させる点が、第1〜第3の実施の形態と異なる。 The fourth embodiment differs from the first to third embodiments in that the objective function is extended.
上記第1の実施の形態では、上記式(4)を目的関数として用い、上記式(4)におけるFp(Ip,yp)及びBp(Ip,yp)は、上記式(5)、(6)で定義される場合について説明した。In the first embodiment, using the above formula (4) as the objective function, F p (I p, y p) in the above formula (4) and B p (I p, y p ) is the formula ( The cases defined in 5) and (6) have been described.
第4の実施の形態では、上記式(4)に示す目的関数は、特定の物体の形状の尤もらしさとして、形状パラメータが表す特定の物体の形状に対する画素における形状ラベルの値に関して単調に変化する単調関数を含むように定義し、目的関数を拡張する。 In the fourth embodiment, the objective function represented by the above equation (4) monotonically changes as the likelihood of the shape of a specific object with respect to the value of the shape label in a pixel corresponding to the shape of the specific object represented by the shape parameter. Define to include a monotonic function and extend the objective function.
具体的には、第4の実施の形態では、上記式(4)におけるFp(Ip,yp)とBp(Ip,yp)とを、以下の式(18)、(19)に示すように定義し、目的関数を拡張する。Specifically, in the fourth embodiment, F p (I p, y p) in the above formula (4) and B p (I p, y p ) and the following formula (18), (19 ) And extend the objective function.
図16に、目的関数の下界が単調性を満たすための十分条件を表す図を示す。この十分条件は,hF p(y),hB p(y)(∀p∈P)がそれぞれyq(∀q∈P)に関して単調減少,単調増加する単調関数であることである。図16に示すように、ラベル対{y,y’}∈L|P|がFIG. 16 is a diagram showing a sufficient condition for the lower bound of the objective function to satisfy monotonicity. The sufficient condition is that h F p (y) and h B p (y) (∀p∈P) are monotonic functions that decrease and increase monotonically with respect to y q (∀q∈P), respectively. As shown in FIG. 16, the label pair {y, y '} L | P |
であるとき、以下の式(20)に示す関係が成り立つ。 , The relationship shown in the following equation (20) holds.
また、図17に、上記式(20)の関係を満たすhF p(y)とhB p(y)の一例を示す。図17に示す例は、以下の式(21)に示す距離関数によってhF p(y)とhB p(y)とを定義したものである。FIG. 17 shows an example of h F p (y) and h B p (y) satisfying the relationship of the above equation (20). In the example shown in FIG. 17, h F p (y) and h B p (y) are defined by the distance function shown in the following equation (21).
上記式(21)に示す符号付き距離関数によって単調関数hF p(y)とhB p(y)とを定義し、形状の輪郭付近でのペナルティを小さくすることで、事前形状の誤りによるセグメンテーションの誤差を低減することができる。The monotone functions h F p (y) and h B p (y) are defined by the signed distance function shown in the above equation (21), and the penalty near the contour of the shape is reduced. Segmentation errors can be reduced.
以上説明したように、第4の実施の形態に係る画像処理装置によれば、目的関数における特定の物体の形状の尤もらしさとして、形状パラメータが表す特定の物体の形状に対する画素における形状ラベルの値に関して単調に変化する単調関数を含むように定義する。これにより、被写体の領域を精度よく抽出することができる。 As described above, according to the image processing apparatus according to the fourth embodiment, as the likelihood of the shape of a specific object in the objective function, the value of the shape label in the pixel with respect to the shape of the specific object represented by the shape parameter Is defined to include a monotone function that changes monotonically with respect to. As a result, the region of the subject can be accurately extracted.
また、第4の実施の形態では、特定の物体の形状の正解の輪郭でエネルギーが最小になるような単調関数hF p(y)とhB p(y)とを目的関数に導入することで、より優れた目的関数が利用可能となる。In the fourth embodiment, a monotone function h F p (y) and h B p (y) that minimizes energy at the correct contour of the shape of a specific object is introduced into the objective function. Thus, a better objective function can be used.
[第5の実施の形態]
次に、第5の実施の形態について説明する。なお、第5の実施の形態に係る統計的形状モデル生成装置及び画像処理装置の構成は、第1の実施の形態と同様の構成となるため、同一符号を付して説明を省略する。[Fifth Embodiment]
Next, a fifth embodiment will be described. Note that the configurations of the statistical shape model generation device and the image processing device according to the fifth embodiment are the same as those of the first embodiment, and thus the same reference numerals are given and the description is omitted.
第5の実施の形態では、固有空間上のパラメータαを形状に写像する関数gを拡張させる点が、第1〜第4の実施の形態と異なる。 The fifth embodiment is different from the first to fourth embodiments in that a function g for mapping the parameter α on the eigenspace to a shape is extended.
上記第1の実施の形態において、上記式(2)に写像の関数gは、上記式(3)のHeaviside function H(・)によって定義される場合について説明した。 In the first embodiment, the case where the function g of the mapping in the above equation (2) is defined by the Heaviside function H (·) in the above equation (3) has been described.
第5の実施の形態では、形状パラメータを特定の物体の形状に写像する関数gは、予め定められた関数fを含むように定義し、写像する関数gを拡張する。 In the fifth embodiment, a function g for mapping a shape parameter to the shape of a specific object is defined to include a predetermined function f, and the function g for mapping is extended.
また、写像する関数gpを以下の式(23)のように拡張した結果、レベルセット関数以外の形状表現も利用可能となる(例えば、Kilian M. Pohl et al.,“Logarithm Odds Maps for Shape Representation”, Proc. of Medical Image Computing and Computer-Assisted Intervention, Vol.4191, pp.955-963, 2006 に記載されているLogOddsなど)。Also, extended result as equation (23) below the functions g p mapping, also made available shape representations other than the level set function (e.g., Kilian M. Pohl et al., "Logarithm Odds Maps for Shape Representation ”, Proc. Of Medical Image Computing and Computer-Assisted Intervention, Vol.4191, pp.955-963, 2006.
ここでfは予め定められた単調関数を表し、tは予め定められた閾値を表す。 Here, f represents a predetermined monotone function, and t represents a predetermined threshold.
例えば、統計的形状モデルが表す可能な形状表現の一例として、LogOddsを用いる場合について説明する。形状を表すクラスが2である場合、統計的形状モデルとしてLogOddsを用いるときは、以下の式で表される。 For example, a case where LogOdds is used as an example of a possible shape expression represented by a statistical shape model will be described. When the class representing the shape is 2, and LogOdds is used as a statistical shape model, it is expressed by the following equation.
また、図18に、統計的形状モデルとしてLogOddsを用いたとき例を示す。図18に示すように、LogOddsを用いることにより、形状表現が拡張される。 FIG. 18 shows an example in which LogOdds is used as a statistical shape model. As shown in FIG. 18, the shape representation is extended by using LogOdds.
以上説明したように、第5の実施の形態に係る統計的形状モデル生成装置によれば、写像する関数gの拡張(予め定められた関数fと閾値tとの導入)により、さらに高度な統計的形状モデルが利用可能となる。 As described above, according to the statistical shape model generation device according to the fifth embodiment, a more advanced statistical data is obtained by expanding the function g to be mapped (introducing the predetermined function f and the threshold value t). The dynamic shape model becomes available.
なお、本発明は、上述した例に限定されるものではなく、この発明の要旨を逸脱しない範囲内で様々な変形や応用が可能である。 Note that the present invention is not limited to the above-described example, and various modifications and applications can be made without departing from the gist of the present invention.
本実施の形態は、基礎的であることから、その応用範囲は広いと考えられる。具体的には、形状などが統計的に変動する図形を対象とした画像認識処理に幅広く利用可能である。例えば、他の臓器や他の医用画像(例えばMRやPETなど)はもちろん、それ以外に、顔画像からの顔の認識や、一般の風景画像から特定の図形を認識する処理などに利用可能である。本実施の形態では、膵臓を対象としたが、他の臓器として例えば脾臓を対象とすることもできる。 Since this embodiment is fundamental, its application range is considered to be wide. Specifically, the present invention can be widely used for image recognition processing for a graphic whose shape or the like fluctuates statistically. For example, it can be used not only for other organs and other medical images (for example, MR and PET), but also for recognition of a face from a face image and processing for recognizing a specific figure from a general landscape image. is there. In the present embodiment, the pancreas is targeted, but the spleen may be targeted as another organ, for example.
また、第1の実施の形態では、統計的形状モデルとしてLSDMを用いる場合を例に説明したが、これに限定されるものではない。LSDMでなく、他の統計的形状モデルの場合にも利用可能である。統計的形状モデルは、上記式(1)のような線形の関数(厳密には画素についての関数φが線形の関数)で表現できていれば、他の統計的形状モデルに本実施の形態を適用することができる。目的関数や最適化手法によっては、実施の形態で示した上記式(1)のような線形の関数で表現される統計的形状モデル以外の統計的形状モデルに対しても、本実施の形態を適用することができる。
また、上記第5の実施の形態に示すように、写像の関数gを拡張させることで、様々な統計的形状モデルを用いることができる。Further, in the first embodiment, the case where LSDM is used as a statistical shape model has been described as an example, but the present invention is not limited to this. It can be used in the case of other statistical shape models instead of LSDM. As long as the statistical shape model can be expressed by a linear function (strictly, the function φ for a pixel is a linear function) as in the above equation (1), the present embodiment can be applied to other statistical shape models. Can be applied. Depending on the objective function and the optimization method, the present embodiment can be applied to a statistical shape model other than the statistical shape model represented by the linear function such as the equation (1) shown in the embodiment. Can be applied.
Further, as shown in the fifth embodiment, various statistical shape models can be used by expanding the mapping function g.
また、第1の実施の形態では、固有空間の次元は2(d=2)である場合を例に説明したが、これに限定されるものではない。例えば、より高次元の固有空間を対象としてもよい。例えば、第2及び第3の実施の形態で示したように、近似解法を用いることによって、固有空間の次元を増加させることができる。 Further, in the first embodiment, the case where the dimension of the eigenspace is 2 (d = 2) has been described as an example, but the present invention is not limited to this. For example, a higher-dimensional eigenspace may be targeted. For example, as shown in the second and third embodiments, the dimension of the eigenspace can be increased by using the approximate solution.
上記実施の形態において、固有空間の次元を3以上に増加させた場合、凸多角形は凸多胞体に対応し、凸多角形の面積は、凸多胞体の体積に対応する。
従って固有空間の次元を3以上に増加させた場合、セグメンテーション部114は、固有空間上の凸多胞体を繰り返し分割して最適な形状パラメータを示す点を含む凸多胞体を探索する探索アルゴリズムに従って、上記式(4)に示す目的関数を最適化するように、入力画像が表す被写体の形状を表す形状パラメータを推定する。
また、セグメンテーション部114は、探索アルゴリズムにおいて、凸多胞体に含まれる形状集合に対する目的関数の下界を計算することにより、最適な形状パラメータを示す点を含む固有空間上の凸多胞体を探索して、形状パラメータを推定する。
また、セグメンテーション部114は、分割後の2つの凸多胞体の体積が対応するように、固有空間上の凸多胞体を分割し、かつ、最適な形状パラメータを示す点を含む固有空間上の凸多胞体を探索することを繰り返して、形状パラメータを推定する。In the above embodiment, when the dimension of the eigenspace is increased to 3 or more, the convex polygon corresponds to the convex polycell, and the area of the convex polygon corresponds to the volume of the convex polycell.
Therefore, when the dimension of the eigenspace is increased to 3 or more, the segmentation unit 114 repeatedly divides the convex polyvesicle on the eigenspace and searches for a convex polyvesicle including a point indicating an optimal shape parameter according to a search algorithm. The shape parameter representing the shape of the subject represented by the input image is estimated so as to optimize the objective function represented by the above equation (4).
In addition, the segmentation unit 114 searches for a convex multivessel on an eigenspace including a point indicating an optimal shape parameter by calculating a lower bound of an objective function for a shape set included in the convex polyvessel in a search algorithm. , And estimate the shape parameters.
In addition, the segmentation unit 114 divides the convex polyvesicle on the eigenspace so that the volumes of the two convex polyvesicles after the division correspond to each other, and generates the convex polyvesicle on the eigenspace including points indicating the optimal shape parameters. The search for multivesicles is repeated to estimate shape parameters.
また、実施の形態では、上記式(17)に従って選択された集合Qに含まれる画素kによって固有空間を分割する場合を例に説明したが、これに限定されるものではなく、他の任意の分割法で作成した凸多角形に対しても適用可能である。例えば、第2の実施の形態で示したように、格子状など、規則的に固有空間を分割した凸多角形にも適用可能である。
この場合、分割の細かさによってはTightness(最適性保証の3つの条件の内の一つ。)が成立しないことがある。Further, in the embodiment, the case where the eigenspace is divided by the pixels k included in the set Q selected according to the above equation (17) has been described as an example. However, the present invention is not limited to this, and the present invention is not limited to this. The present invention is also applicable to convex polygons created by the division method. For example, as shown in the second embodiment, the present invention can be applied to a convex polygon in which an eigenspace is regularly divided, such as a lattice.
In this case, Tightness (one of the three conditions for guaranteeing optimality) may not be satisfied depending on the degree of division.
ここで、Tightnessについて簡単に説明する。 Here, Tightness will be briefly described.
上記式(14)を、以下の式(24)のように示す。 The above equation (14) is shown as the following equation (24).
この場合、すべての葉ノードHiについて以下の式(25)が成立する場合には、Tightnessが成立している。In this case, when the following equation for all the leaf nodes H i (25) is satisfied, tightness is established.
分割の細かさが粗い場合には、上記式(25)が成立せずに、Tightnessが成立しないこともあり、近似解しか得られないことがあるが、固有空間の次元dを大きくした場合の計算コストのさらなる削減が可能であるというメリットがある。注意すべき点は、最適性をやや犠牲にすることで、計算コストと計算精度のバランスをユーザが選ぶことが可能になる点であり、本実施の形態の提案法の拡張に関する重要な点である。なお、このとき、最適性が保証されない点では上記米国特許第8249349号明細書などと同様であるが、上記米国特許第8249349号明細書と比べて、前処理が一切不要である点に注意されたい。 When the fineness of the division is coarse, Expression (25) does not hold, and Tightness may not hold, and only an approximate solution may be obtained. However, when the dimension d of the eigenspace is increased, There is an advantage that the calculation cost can be further reduced. It should be noted that by slightly sacrificing the optimality, the user can select the balance between the calculation cost and the calculation accuracy, which is an important point regarding the extension of the proposed method of the present embodiment. is there. At this time, the point that the optimality is not guaranteed is the same as in the above-mentioned US Pat. No. 8,249,349, but it should be noted that no pre-processing is required as compared with the above-mentioned US Pat. No. 8,249,349. I want to.
また、最小化の目的関数やアルゴリズムについても、今回は上記実施の形態で示した目的関数、および、branch and bound法とグラフカット法とを用いたが、最小化が可能であれば、今回示したものに限らず適用可能である。 In addition, the objective function and algorithm for minimization used the objective function shown in the above embodiment, and the branch and bound method and the graph cut method. The present invention is not limited to this and can be applied.
また、上記第2又は第3の実施の形態に係る近似解法を実行した後に、上記第1の実施の形態に係る厳密解法を実行してもよい。 After the approximate solution according to the second or third embodiment is executed, the exact solution according to the first embodiment may be executed.
また、図2に示したコンピュータ構成において、実施の形態に係る各処理部の機能を実現するためのプログラムをコンピュータ読み取り可能な記録媒体に記録して、この記録媒体に記録されたプログラムをコンピュータシステムに読み込ませ、実行することにより、各構成による処理が実行されてもよいし、図示されていない通信機能を用いて、当該プログラムを読み込ませることでもよい。 In the computer configuration shown in FIG. 2, a program for realizing the function of each processing unit according to the embodiment is recorded on a computer-readable recording medium, and the program recorded on the recording medium is stored in a computer system. The processing by each component may be executed by reading and executing the program, or the program may be read by using a communication function (not shown).
なお、コンピュータ読み取り可能な記録媒体とは、フレキシブルディスク、光磁気ディスク、ROM、CD−ROM等の可搬媒体、コンピュータシステムに内蔵されるハードディスク等の記憶装置のことをいう。 Note that the computer-readable recording medium refers to a portable medium such as a flexible disk, a magneto-optical disk, a ROM, and a CD-ROM, and a storage device such as a hard disk built in a computer system.
また、上記プログラムは、このプログラムを記憶装置等に格納したコンピュータシステムから、伝送媒体を介して、あるいは、伝送媒体中の伝送波により他のコンピュータシステムに伝送されてもよい。ここで、プログラムを伝送する「伝送媒体」は、インターネット等のネットワーク(通信網)や電話回線等の通信回線(通信線)のように情報を伝送する機能を有する媒体のことをいう。 Further, the above program may be transmitted from a computer system storing the program in a storage device or the like to another computer system via a transmission medium or by a transmission wave in the transmission medium. Here, the "transmission medium" for transmitting a program refers to a medium having a function of transmitting information, such as a network (communication network) such as the Internet or a communication line (communication line) such as a telephone line.
また、上記プログラムは、前述した機能の一部を実現するためのものであってもよい。さらに、前述した機能を、コンピュータシステムにすでに記録されているプログラムとの組み合わせで実現できるもの、いわゆる差分ファイル(差分プログラム)であってもよい。 Further, the program may be for realizing a part of the functions described above. Furthermore, what can implement | achieve the function mentioned above in combination with the program already recorded on the computer system, and what is called a difference file (difference program) may be sufficient.
また、本実施の形態では、上記図1に示す統計的形状モデル生成装置10及び上記図7に示す画像処理装置100の各処理部を、プログラムにより各機能の実行が可能なコンピュータで構成するものとしているが、論理素子回路からなるハードウェア構成とすることでも良い。 In the present embodiment, each processing unit of the statistical shape model generation device 10 shown in FIG. 1 and the image processing device 100 shown in FIG. 7 is configured by a computer capable of executing each function by a program. However, a hardware configuration including a logic element circuit may be adopted.
また、上記図2に示す統計的形状モデル生成装置10及び画像処理装置100のコンピュータ構成に関しても適宜にその構成内容を変更しても良い。 Further, the computer configuration of the statistical shape model generation device 10 and the image processing device 100 shown in FIG. 2 may be changed as appropriate.
このように、実施する形態例を、図面を参照して詳述してきたが、具体的な構成はこの実施の形態例に限られるものではなく、要旨を逸脱しない範囲の設計等も含まれる。 As described above, the embodiment has been described in detail with reference to the drawings. However, the specific configuration is not limited to the embodiment and includes a design and the like within a range not departing from the gist.
また、実施の形態のプログラムは、記憶媒体に格納して提供するようにしてもよい。 Further, the program of the embodiment may be provided by being stored in a storage medium.
実施の形態のコンピュータ可読媒体は、特定の物体である被写体を表す入力画像から、前記被写体の領域を抽出するプログラムであって、コンピュータを、前記入力画像を受け付ける受付手段、及び(A)前記被写体を表す学習用の複数の画像であって、前記被写体の領域が予め求められた前記複数の画像に基づいて予め計算された固有ベクトルを基底とする固有空間であって、かつ、(B)固有空間上の点が、前記特定の物体の形状の統計的変動を表す統計的形状モデルの形状パラメータを示す固有空間において、前記入力画像に基づいて、前記固有空間上の点が示す前記形状パラメータが表す前記特定の物体の形状の尤もらしさと前記入力画像中の隣接画素間の画素値の差とに応じた値を表す予め定められた目的関数を最適化するように、前記入力画像が表す前記被写体の形状を表す前記形状パラメータを推定し、前記推定された前記形状パラメータが表す前記特定の物体の形状を事前知識として、前記被写体の領域を、前記入力画像から抽出するセグメンテーション手段として機能させるためのプログラムを記憶したコンピュータ可読媒体である。 A computer-readable medium according to an embodiment is a program for extracting a region of a subject from an input image representing a subject, which is a specific object. And (B) an eigenspace in which the region of the subject is an eigenspace based on eigenvectors calculated in advance based on the plurality of images obtained in advance. In the eigenspace where the upper point indicates the shape parameter of the statistical shape model representing the statistical variation of the shape of the specific object, the shape parameter indicated by the point on the eigenspace represents based on the input image. To optimize a predetermined objective function representing a value corresponding to the likelihood of the shape of the specific object and a difference between pixel values of adjacent pixels in the input image. Estimating the shape parameter representing the shape of the subject represented by the input image, and extracting the area of the subject from the input image, using the shape of the specific object represented by the estimated shape parameter as prior knowledge. It is a computer-readable medium storing a program for causing it to function as a segmentation unit.
日本出願2014−169911の開示はその全体が参照により本明細書に取り込まれる。 The disclosure of Japanese Application 2014-169911 is incorporated herein by reference in its entirety.
本明細書に記載された全ての文献、特許出願、及び技術規格は、個々の文献、特許出願、及び技術規格が参照により取り込まれることが具体的かつ個々に記載された場合と同程度に、本明細書中に参照により取り込まれる。 All documents, patent applications, and technical standards mentioned herein are to the same extent as if each individual document, patent application, and technical standard were specifically and individually stated to be incorporated by reference. Incorporated herein by reference.
Claims (19)
前記入力画像を受け付ける受付手段と、
(A)前記被写体を表す学習用の複数の画像であって、前記被写体の領域が予め求められた前記複数の画像に基づいて予め計算された固有ベクトルを基底とする固有空間であって、かつ、(B)固有空間上の点が、前記特定の物体の形状の統計的変動を表す統計的形状モデルの形状パラメータを示す固有空間において、
前記入力画像に基づいて、前記固有空間上の点が示す前記形状パラメータが表す前記特定の物体の形状の尤もらしさと前記入力画像中の隣接画素間の画素値の差とに応じた値を表す予め定められた目的関数を最適化するように、前記入力画像が表す前記被写体の形状を表す前記形状パラメータを推定し、
前記推定された前記形状パラメータが表す前記特定の物体の形状を事前知識として、前記被写体の領域を、前記入力画像から抽出するセグメンテーション手段と、
を含む画像処理装置。 An image processing apparatus for extracting an area of the subject from an input image representing a subject which is a specific object,
Receiving means for receiving the input image;
(A) a plurality of learning images representing the subject, wherein the subject area is an eigenspace based on eigenvectors calculated in advance based on the plurality of images obtained in advance, and (B) In the eigenspace in which a point on the eigenspace indicates a shape parameter of a statistical shape model representing a statistical variation of the shape of the specific object,
Based on the input image, a value corresponding to a likelihood of the shape of the specific object represented by the shape parameter indicated by the point on the eigenspace and a difference in pixel value between adjacent pixels in the input image. Estimating the shape parameters representing the shape of the subject represented by the input image so as to optimize a predetermined objective function,
As a priori knowledge of the shape of the specific object represented by the estimated shape parameters, a region of the subject, a segmentation unit that extracts from the input image,
An image processing apparatus including:
前記入力画像に基づいて、前記目的関数を最適化するように、前記入力画像が表す前記被写体の形状を表す前記形状パラメータを推定すると同時に、
前記推定された前記形状パラメータが表す前記特定の物体の形状を事前知識として、前記被写体の領域を、前記入力画像から抽出する請求項1記載の画像処理装置。 The segmentation means includes:
Based on the input image, so as to optimize the objective function, while estimating the shape parameters representing the shape of the subject represented by the input image,
The image processing apparatus according to claim 1, wherein a region of the subject is extracted from the input image, using the shape of the specific object represented by the estimated shape parameter as prior knowledge.
最適な形状パラメータを示す点が表す前記特定の物体の形状を含む形状集合を表す前記固有空間上の凸多胞体を繰り返し分割して最適な形状パラメータを示す点を含む前記固有空間上の凸多胞体を探索する探索アルゴリズムに従って、前記目的関数を最適化するように、前記入力画像が表す前記被写体の形状を表す前記形状パラメータを推定し、
前記推定された前記形状パラメータが表す前記特定の物体の形状を事前知識として、前記被写体の領域を、前記入力画像から抽出する請求項1又は請求項2記載の画像処理装置。 The segmentation means comprises:
The convex polymorphism on the eigenspace including points indicating the optimal shape parameters is obtained by repeatedly dividing the convex multicellular body on the eigenspace representing the shape set including the shape of the specific object represented by the point indicating the optimal shape parameter. According to a search algorithm for searching for a cell, to optimize the objective function, estimate the shape parameter representing the shape of the subject represented by the input image,
3. The image processing apparatus according to claim 1, wherein the region of the subject is extracted from the input image, using the shape of the specific object represented by the estimated shape parameter as prior knowledge. 4.
前記推定された前記形状パラメータが表す前記特定の物体の形状を事前知識として、前記被写体の領域を、前記入力画像から抽出する請求項4又は請求項5記載の画像処理装置。 The segmentation means sets a sampling point on the eigenspace, and divides the convex polyhedron using a hyperplane determined from the sampling points set on the eigenspace to indicate an optimal shape parameter. Iterating the search for convex multivesicular bodies on the eigenspace including, estimating the shape parameters representing the shape of the subject represented by the input image,
The image processing device according to claim 4, wherein the region of the subject is extracted from the input image, using the shape of the specific object represented by the estimated shape parameter as prior knowledge.
前記受付手段が、前記入力画像を受け付けるステップと、
前記セグメンテーション手段が、(A)前記被写体を表す学習用の複数の画像であって、前記被写体の領域が予め求められた前記複数の画像に基づいて予め計算された固有ベクトルを基底とする固有空間であって、かつ、(B)固有空間上の点が、前記特定の物体の形状の統計的変動を表す統計的形状モデルの形状パラメータを示す固有空間において、前記入力画像に基づいて、前記固有空間上の点が示す前記形状パラメータが表す前記特定の物体の形状の尤もらしさと前記入力画像中の隣接画素間の画素値の差とに応じた値を表す予め定められた目的関数を最適化するように、前記入力画像が表す前記被写体の形状を表す前記形状パラメータを推定し、前記推定された前記形状パラメータが表す前記特定の物体の形状を事前知識として、前記被写体の領域を、前記入力画像から抽出するステップと、
を含む画像処理方法。 An image processing method in an image processing device that includes an accepting unit and a segmentation unit, and extracts an area of the subject from an input image representing a subject that is a specific object,
The receiving means receiving the input image;
The segmentation means includes: (A) a plurality of learning images representing the subject, wherein a region of the subject is an eigenspace based on eigenvectors calculated in advance based on the plurality of images obtained in advance. And (B) in the eigenspace indicating a shape parameter of a statistical shape model representing a statistical variation of the shape of the specific object, based on the input image, Optimize a predetermined objective function representing a value corresponding to the likelihood of the shape of the specific object represented by the shape parameter indicated by the upper point and a difference in pixel value between adjacent pixels in the input image. As described above, the shape parameter representing the shape of the object represented by the input image is estimated, and the shape of the specific object represented by the estimated shape parameter is determined in advance as prior knowledge. The area of the object, extracting from the input image,
An image processing method including:
前記固有空間において、
前記入力画像に基づいて、前記目的関数を最適化するように、前記入力画像が表す前記被写体の形状を表す前記形状パラメータを推定すると同時に、
前記推定された前記形状パラメータが表す前記特定の物体の形状を事前知識として、前記被写体の領域を、前記入力画像から抽出する請求項10記載の画像処理方法。 The step of extracting the area of the subject from the input image,
In the eigenspace,
Based on the input image, so as to optimize the objective function, while estimating the shape parameters representing the shape of the subject represented by the input image,
The image processing method according to claim 10, wherein a region of the subject is extracted from the input image, using the shape of the specific object represented by the estimated shape parameter as prior knowledge.
最適な形状パラメータを示す点が表す前記特定の物体の形状を含む形状集合を表す前記固有空間上の凸多胞体を繰り返し分割して最適な形状パラメータを示す点を含む前記固有空間上の凸多胞体を探索する探索アルゴリズムに従って、前記目的関数を最適化するように、前記入力画像が表す前記被写体の形状を表す前記形状パラメータを推定し、前記推定された前記形状パラメータが表す前記特定の物体の形状を事前知識として、前記被写体の領域を、前記入力画像から抽出する請求項11記載の画像処理方法。 The step of extracting the area of the subject from the input image,
The convex polymorphism on the eigenspace including points indicating the optimal shape parameters is obtained by repeatedly dividing the convex multicellular body on the eigenspace representing the shape set including the shape of the specific object represented by the point indicating the optimal shape parameter. According to a search algorithm for searching for a cell, the shape parameter representing the shape of the subject represented by the input image is estimated so as to optimize the objective function, and the specific object represented by the estimated shape parameter is The image processing method according to claim 11, wherein a region of the subject is extracted from the input image using a shape as prior knowledge.
前記推定された前記形状パラメータが表す前記特定の物体の形状を事前知識として、前記被写体の領域を、前記入力画像から抽出する請求項13又は請求項14記載の画像処理方法。 The step of extracting the region of the subject from the input image by the segmentation means includes setting a sampling point on the eigenspace, and using a hyperplane determined from the sampling point set on the eigenspace. The shape parameters representing the shape of the subject represented by the input image are estimated by repeating the search for the convex polyvesicles on the eigenspace including the points indicating the optimal shape parameters by dividing the convex polyvesicles. ,
The image processing method according to claim 13, wherein the region of the subject is extracted from the input image, using the shape of the specific object represented by the estimated shape parameter as prior knowledge.
Applications Claiming Priority (3)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
JP2014169911 | 2014-08-22 | ||
JP2014169911 | 2014-08-22 | ||
PCT/JP2015/073277 WO2016027840A1 (en) | 2014-08-22 | 2015-08-19 | Image processing device, method, and program |
Publications (2)
Publication Number | Publication Date |
---|---|
JPWO2016027840A1 JPWO2016027840A1 (en) | 2017-06-01 |
JP6661196B2 true JP6661196B2 (en) | 2020-03-11 |
Family
ID=55350782
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
JP2016544238A Active JP6661196B2 (en) | 2014-08-22 | 2015-08-19 | Image processing apparatus, method, and program |
Country Status (2)
Country | Link |
---|---|
JP (1) | JP6661196B2 (en) |
WO (1) | WO2016027840A1 (en) |
Family Cites Families (3)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
JP2004188201A (en) * | 2002-12-10 | 2004-07-08 | Eastman Kodak Co | Method to automatically construct two-dimensional statistical form model for lung area |
US8131038B2 (en) * | 2007-08-21 | 2012-03-06 | Siemens Aktiengesellschaft | System and method for global-to-local shape matching for automatic liver segmentation in medical imaging |
US9495752B2 (en) * | 2012-09-27 | 2016-11-15 | Siemens Product Lifecycle Management Software Inc. | Multi-bone segmentation for 3D computed tomography |
-
2015
- 2015-08-19 JP JP2016544238A patent/JP6661196B2/en active Active
- 2015-08-19 WO PCT/JP2015/073277 patent/WO2016027840A1/en active Application Filing
Also Published As
Publication number | Publication date |
---|---|
WO2016027840A1 (en) | 2016-02-25 |
JPWO2016027840A1 (en) | 2017-06-01 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
Bitter et al. | Penalized-distance volumetric skeleton algorithm | |
CN109741343B (en) | T1WI-fMRI image tumor collaborative segmentation method based on 3D-Unet and graph theory segmentation | |
US8135189B2 (en) | System and method for organ segmentation using surface patch classification in 2D and 3D images | |
EP2082370B1 (en) | Generalized rigid alignment of 3d ear impression models | |
Chan et al. | Volumetric parametrization from a level set boundary representation with PHT-splines | |
US20130127847A1 (en) | System and Method for Interactive Image-based Modeling of Curved Surfaces Using Single-view and Multi-view Feature Curves | |
JP2008511366A (en) | Feature-weighted medical object contour detection using distance coordinates | |
JP2003180654A (en) | Method for forming three-dimensional statistical shape model for left ventricle from non-dense two-dimensional contour input value, and program storing device which stores program for performing the method | |
US9836834B2 (en) | Topology-preserving ROI remapping method between medical images | |
US9984311B2 (en) | Method and system for image segmentation using a directed graph | |
Hu et al. | Surface segmentation for polycube construction based on generalized centroidal Voronoi tessellation | |
CN107516314B (en) | Medical image hyper-voxel segmentation method and device | |
CN112419343A (en) | System and method for image segmentation | |
Heitz et al. | Statistical shape model generation using nonrigid deformation of a template mesh | |
US20160019435A1 (en) | Image processing apparatus, non-transitory computer-readable recording medium having stored therein image processing program, and operation method of image processing apparatus | |
CA2643865A1 (en) | Method and system for locating landmarks on 3d models | |
Hosseini-Asl et al. | Lung segmentation based on nonnegative matrix factorization | |
Pereañez et al. | A framework for the merging of pre-existing and correspondenceless 3D statistical shape models | |
JP6661196B2 (en) | Image processing apparatus, method, and program | |
Korez et al. | Intervertebral disc segmentation in MR images with 3D convolutional networks | |
EP2534613B1 (en) | Image analysis | |
JP2006031114A (en) | Image division processing system | |
CN112581513B (en) | Cone beam computed tomography image feature extraction and corresponding method | |
Wu et al. | Real‐Time Symmetry‐Preserving Deformation | |
Zhang et al. | Automatic liver segmentation from CT scans based on a statistical shape model |
Legal Events
Date | Code | Title | Description |
---|---|---|---|
A621 | Written request for application examination |
Free format text: JAPANESE INTERMEDIATE CODE: A621 Effective date: 20180518 |
|
A131 | Notification of reasons for refusal |
Free format text: JAPANESE INTERMEDIATE CODE: A131 Effective date: 20190709 |
|
A521 | Request for written amendment filed |
Free format text: JAPANESE INTERMEDIATE CODE: A523 Effective date: 20190904 |
|
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: 20200107 |
|
A61 | First payment of annual fees (during grant procedure) |
Free format text: JAPANESE INTERMEDIATE CODE: A61 Effective date: 20200204 |
|
R150 | Certificate of patent or registration of utility model |
Ref document number: 6661196 Country of ref document: JP Free format text: JAPANESE INTERMEDIATE CODE: R150 |
|
R250 | Receipt of annual fees |
Free format text: JAPANESE INTERMEDIATE CODE: R250 |
|
R250 | Receipt of annual fees |
Free format text: JAPANESE INTERMEDIATE CODE: R250 |