JP5037705B2 - 画像処理装置および方法並びにプログラム - Google Patents

Description

本発明は、画像中の線状構造物または面状構造物を判別するための画像処理装置および方法並びに画像処理方法をコンピュータに実行させるためのプログラムに関するものである。

近年、医療機器（例えば多検出器型ＣＴ等）の進歩により質の高い３次元画像が画像診断に用いられる。３次元画像は多数の２次元画像から構成され情報量が多いため、医師が所望の観察部位を見つけ診断することに時間を要する場合がある。そこで、注目する臓器を抽出しＭＩＰ、ＶＲ、ＣＰＲ等の表示を行うことにより、臓器全体や病変の視認性を高め診断の効率化を図ることが行われている。

一方、医用画像中の血管および骨を抽出する手法として、ヘッセ行列を用いたヘッセ解析が提案されている（非特許文献１参照）。ヘッセ解析は、ガウシアンカーネル等の所定のフィルタについての二次微分カーネルを用いることにより算出した、二階の偏微分係数を要素とするヘッセ行列の固有値を解析することにより、画像中の局所構造が点、線および面のいずれであるかを判別するものである。ヘッセ解析を用いることにより血管は線状構造物として、骨は面状構造物として判別することが可能である。

A. F. Frangi, et. al., "Multiscale vessel enhancement filtering", MICCAI 1998.

しかしながら、血管および骨には様々な太さおよび大きさを有するものがあるため、ヘッセ解析を用いるのみでは、血管および骨を判別できない場合がある。例えば、面状構造物である脊椎等の皮質骨を抽出する場合において、皮質骨の面部分を検出可能なサイズのフィルタカーネルを使用した場合、図８の斜線部分に示すように骨１００の面状構造となる部分のみならず、太い血管１０２の表面を面状構造と判別してしまう場合がある。これは、二次偏微分が太い血管の表面に対して垂直な方向には高く、表面に平行な２方向には低いことから、面状構造と同じ特徴を示すために生じる問題である。また逆に、線状構造物である血管を抽出する場合において、血管の線状部分を検出可能なサイズのフィルタカーネルを使用した場合、図９の斜線部分に示すように血管１０２の線状構造となる部分のみならず、骨１００の角部を線状構造と判別してしまう場合がある。

本発明は、上記事情に鑑みなされたものであり、画像に含まれる線状構造物と面状構造物との誤判別を防止することを目的とする。

本発明による第１の画像処理装置は、画像中の各画素位置における画素値の二階偏微分行列および一階偏微分値を算出する微分値算出手段と、
前記二階偏微分行列値に基づいて、前記画素位置における線状構造らしさおよび／または面状構造らしさの評価値を算出する評価手段とを備え、
前記評価手段は、前記一階偏微分値が大きいほど、前記評価値を小さくする手段であることを特徴とするものである。

なお、本発明による第１の画像処理装置においては、前記微分値算出手段を、複数サイズのフィルタを用いて前記二階偏微分行列および前記一階偏微分値を算出するに際し、前記二階偏微分行列を算出するフィルタのサイズよりも大きいサイズのフィルタにより、前記一階偏微分値を算出する手段としてもよい。

また、本発明による第１の画像処理装置においては、前記微分値算出手段を、前記画像を多重解像度変換して解像度が異なる複数の解像度画像を取得し、該複数の解像度画像の対応する画素位置において所定サイズのフィルタを用いて前記二階偏微分行列および前記一階偏微分値を算出するに際し、前記二階偏微分行列を算出する解像度画像よりも低い解像度の解像度画像により、前記一階偏微分値を算出する手段としてもよい。

また、本発明による第１の画像処理装置においては、前記微分値算出手段を、１次元の基本ガウシアンカーネル、該基本ガウシアンカーネルを一階微分した一次微分カーネル、および該基本ガウシアンカーネルを二階微分した二次微分カーネルを用いて、前記二階偏微分行列および前記一階偏微分値を算出する手段としてもよい。

また、本発明による第１の画像処理装置においては、前記評価値に基づいて、前記画像の前記各画素に対して、対象領域らしさ、背景領域らしさ、および隣接する画素が同一の領域であることの確からしさを設定して、該対象領域および該背景領域を分割する分割手段をさらに備えるものとしてもよい。

この場合、前記分割手段を、前記対象領域らしさを前記線状構造らしさの評価値に基づいて設定し、前記背景領域らしさを前記面状構造らしさの評価値に基づいて設定する手段としてもよい。

また、前記分割手段を、前記同一の領域であることの確からしさを、前記二階偏微分行列値にのみ基づいて算出された、前記画素位置における線状構造らしさおよび面状構造らしさの評価値、または前記線状構造らしさの評価値および前記面状構造らしさの評価値における、前記一階偏微分値の影響を小さくした評価値に基づいて設定する手段としてもよい。

本発明による第１の画像処理方法は、画像中の各画素位置における画素値の二階偏微分行列および一階偏微分値を算出し、前記二階偏微分行列値に基づいて、前記画素位置における線状構造らしさおよび／または面状構造らしさの評価値を算出するに際し、前記一階偏微分値が大きいほど、前記評価値を小さくすることを特徴とするものである。

本発明による第２の画像処理装置は、画像中の各画素位置における画素値の二階偏微分行列および一階偏微分値を算出する微分値算出手段と、
前記二階偏微分行列値に基づいて、前記画素位置における線状構造らしさの評価値を算出する評価手段とを備え、
前記評価手段は、前記一階偏微分値の大きさに基づいて、前記評価値を変更する手段であることを特徴とするものである。

なお、本発明による第２の画像処理装置においては、前記評価手段を、前記二階偏微分行列値について、該二階偏微分方向に一致する前記一階偏微分値の大きさに応じて前記二階偏微分行列値を補正し、該補正された二階偏微分行列値に基づいて、前記評価値を算出する手段としてもよい。

また、本発明による第２の画像処理装置においては、前記微分値算出手段を、複数サイズのフィルタを用いて前記二階偏微分行列および前記一階偏微分値を算出するに際し、前記二階偏微分行列を算出するフィルタのサイズよりも大きいサイズのフィルタにより、前記一階偏微分値を算出する手段としてもよい。

また、本発明による第２の画像処理装置においては、前記微分値算出手段を、前記画像を多重解像度変換して解像度が異なる複数の解像度画像を取得し、該複数の解像度画像の対応する画素位置において所定サイズのフィルタを用いて前記二階偏微分行列および前記一階偏微分値を算出するに際し、前記二階偏微分行列を算出する解像度画像よりも低い解像度の解像度画像により、前記一階偏微分値を算出する手段としてもよい。

また、本発明による第２の画像処理装置においては、前記微分値算出手段を、１次元の基本ガウシアンカーネル、該基本ガウシアンカーネルを一階微分した一次微分カーネル、および該基本ガウシアンカーネルを二階微分した二次微分カーネルを用いて、前記二階偏微分行列および前記一階偏微分値を算出する手段としてもよい。

また、本発明による第２の画像処理装置においては、前記評価値に基づいて、前記画像の前記各画素に対して、対象領域らしさ、背景領域らしさ、および隣接する画素が同一の領域であることの確からしさを設定して、該対象領域および該背景領域を分割する分割手段をさらに備えるものとしてもよい。

この場合、前記分割手段を、前記対象領域らしさを前記線状構造らしさの評価値に基づいて設定する手段としてもよい。

また前記分割手段を、前記同一の領域であることの確からしさを、前記二階偏微分行列値にのみ基づいて算出された、前記画素位置における線状構造らしさの評価値、または前記線状構造らしさの評価値における、前記一階偏微分値の影響を小さくした評価値に基づいて設定する手段としてもよい。

本発明による第２の画像処理方法は、画像中の各画素位置における画素値の二階偏微分行列および一階偏微分値を算出し、
前記二階偏微分行列値に基づいて、前記画素位置における線状構造らしさの評価値を算出するに際し、前記一階偏微分値の大きさに基づいて、前記評価値を変更することを特徴とするものである。

なお、本発明による第１および第２の画像処理方法をコンピュータに実行させるためのプログラムとして提供してもよい。

医用画像に含まれる血管等の線状構造物をヘッセ行列に基づいて判別した際の骨構造の誤抽出、および皮質骨等の面状構造物を判別した際の血管の誤抽出においては、誤抽出した部分の輝度値のパターンに１次元的な偏りがある。例えば、図８の誤抽出は、血管の内側の輝度値が高く外側の輝度値が低いという１次元的な偏りがあることがわかる。このように輝度値のパターンに１次元的な偏りがあると、一階偏微分値が大きくなる。一方、理想的な線状構造および面状構造は、その中心を基準として対称な構造をなしているため一階偏微分値は０となる。本発明の第１の画像処理装置および方法によれば、画像中の任意の画素位置における画素値の二階偏微分行列および一階偏微分値を算出し、二階偏微分行列値に基づいて、画素位置における線状構造らしさおよび／または面状構造らしさの評価値を算出するに際して、一階偏微分値が大きいほど評価値を小さくするようにしたものである。このため、画像に含まれる線状構造物と面状構造物との誤判別を防止して、線状構造物と面状構造物とを精度良く判別できる。

一方、理想的な線状構造においては、二階偏微分行列値を固有値分解することにより得られた３方向の固有値のうち、線状構造からなる組織の主軸に直交する方向の２つの固有値の大きさは略等しくなる。しかしながら、例えば心臓の辺縁を走行する冠動脈等のように、面状構造物の近傍に線状構造物が存在する場合、線状構造物から面状構造物の方向（すなわち面状構造物の法線方向）における固有値が大きくなるため、線状構造らしさの評価値が低くなり、その結果、線状構造物の判別に失敗するおそれがある。本発明の第２の画像処理装置および方法によれば、画像中の任意の画素位置における画素値の二階偏微分行列および一階偏微分値を算出し、二階偏微分行列値に基づいて、画素位置における線状構造らしさの評価値を算出するに際して、一階偏微分値の大きさに基づいて評価値を変更するようにしたものである。とくに二階偏微分行列値について、二階偏微分方向に一致する一階偏微分値の大きさに応じて二階偏微分行列値を補正し、補正された二階偏微分行列値に基づいて評価値を算出するようにしたものである。このため、面状構造物の近傍に線状構造物が存在する場合であっても、面状構造物が存在する方向およびこれに直交する方向における二階偏微分行列値を略等しくすることができ、これにより、線状構造らしさの評価値が低くなることを防止でき、その結果、線状構造物を精度よく判別することができる。

また、複数サイズのフィルタを用いて二階偏微分行列および一階偏微分値を算出するに際し、二階偏微分行列を算出するフィルタのサイズよりも大きいサイズのフィルタにより、一階偏微分値を算出する、あるいは画像を多重解像度変換して解像度が異なる複数の解像度画像を取得し、複数の解像度画像の対応する画素位置に対して所定サイズのフィルタを用いて二階偏微分行列および一階偏微分値を算出するに際し、二階偏微分行列を算出する解像度画像の解像度よりも低い解像度の解像度画像により一階偏微分値を算出することにより、輝度値の１次元的な偏りを捉え易くなるため、より精度良く線状構造らしさおよび面状構造らしさの評価値を算出することができる。

また、１次元の基本ガウシアンカーネル、基本ガウシアンカーネルを一階微分した一次微分カーネル、および基本ガウシアンカーネルを二階微分した二次微分カーネルを用いて、二階偏微分行列および一階偏微分値を算出することにより、比較的簡易な演算により二階偏微分行列および一階偏微分値を算出することができるため、演算速度を向上させることができる。

また、評価値に基づいて、画像の各画素に対して、対象領域らしさ、背景領域らしさ、および隣接する画素が同一の領域であることの確からしさを設定して、対象領域および背景領域を分割することにより、対象領域および背景領域を精度良く分割できる。

とくに、対象領域らしさを線状構造らしさの評価値に基づいて設定し、背景領域らしさを面状構造らしさの評価値に基づいて設定することにより、画像中の血管等の線状構造物と、骨等の面状構造物とを精度良く分割することができる。

本発明の実施形態による画像処理装置の構成を示す概略ブロック図 多重解像度変換を説明するための図 ガウシアンカーネルを示す図 線状構造の固有値を説明するための図 面状構造の固有値を説明するための図 Graph Cut領域分割方法を説明するための図 本実施形態において行われる処理を示すフローチャート 面状構造物の誤検出を説明するための図 線状構造物の誤検出を説明するための図 面状構造物の近傍に線状構造物が存在する場合における線状構造物の誤検出を説明するための図

以下、図面を参照して本発明の実施形態について説明する。図１は本発明の実施形態による画像処理装置の構成を示す概略ブロック図である。なお、図１のような画像処理装置１の構成は、補助記憶装置（不図示）に読み込まれたプログラムをコンピュータ（例えばパーソナルコンピュータ等）上で実行することにより実現される。また、このプログラムは、ＣＤ−ＲＯＭ等の情報記憶媒体に記憶され、もしくはインターネット等のネットワークを介して配布され、コンピュータにインストールされることになる。

画像処理装置１は、例えばＸ線ＣＴ装置２により撮像された複数の２次元画像を用いて、３次元画像Ｍ０を生成し、この３次元画像Ｍ０に含まれる線状構造物と面状構造物とを自動的に分割するものであって、画像取得部１０、検出領域設定部２０、判別部３０、分割部４０、表示部５０および入力部６０を備える。

画像取得部１０は、例えばＸ線ＣＴ装置２により撮像された複数のＣＴ画像（２次元画像）を取得し、複数の２次元画像から３次元画像Ｍ０を生成する。なお、画像取得部１０は、ＣＴ画像のみならず、いわゆるＭＲＩ画像、ＲＩ画像、ＰＥＴ画像、Ｘ線画像等の２次元画像を取得するものであってもよい。

検出領域設定部２０は、まず３次元画像Ｍ０のボクセルサイズを等方化する。例えば、３次元画像Ｍ０のボクセルサイズが、３次元画像Ｍ０におけるＸ，Ｙ，Ｚ方向のそれぞれにおいて、０．３ｍｍ、０．３ｍｍ、０．６ｍｍであった場合、（Ｘ，Ｙ，Ｚ）＝（０．５，０．５，０．５）（ｍｍ）に等方化する。このように、３次元画像Ｍ０のボクセルサイズを等方化することにより、後述する判別部３０において、Ｘ，Ｙ，Ｚ方向のそれぞれに同一サイズのカーネルを適用することができるため、演算を簡易に行うことが可能となる。

検出領域設定部２０は、等方化後に３次元画像Ｍ０を多重解像度変換して図２に示すように、解像度が異なる複数の３次元多重解像度画像Ｍｓｉ（ｉ＝０〜ｎ）（ガウシアンピラミッド）を生成する。なお、ｉ＝０は３次元画像Ｍ０と同一解像度、ｉ＝ｎは最低解像度を表す。３次元多重解像度画像Ｍｓｉのボクセルサイズは、解像度が高い順に、（Ｘ，Ｙ，Ｚ）＝（０．５，０．５，０．５）、（１．０，１．０，１．０）、（２．０，２．０，２．０）…となる。

判別部３０は、フィルタリング部３２および評価部３４を備える。フィルタリング部３２は、ヘッセ行列を用いたヘッセ解析を行うために、３次元多重解像度画像Ｍｓｉのそれぞれに、ガウシアンカーネルを用いたフィルタリングを行う。すなわち、解像度が異なる３次元多重解像度画像Ｍｓｉのそれぞれに対して、同一サイズのフィルタカーネル（σ＝１．０とする）を畳み込む。これにより、３次元画像Ｍ０に対して、実質的に異なるサイズのフィルタカーネルを適用することとなるため、異なるサイズの線状構造物（例えば、血管）および面状構造物（例えば皮質骨等の骨）を検出することができる。

以下、ヘッセ解析について説明する。ヘッセ解析に使用されるヘッセ行列は３次元画像に対しては下記の式（１）に示すように３×３の行列となる。

また、ガウシアンカーネル関数ｆを用いた場合、ヘッセ行列を得るためのフィルタ係数は、式（２）に示す１次元の基本カーネル、基本カーネルを一階微分した一次微分カーネル、および基本カーネルを二階微分した二次微分カーネルよって求められる。なお、式（２）はＸ方向におけるフィルタ係数のみを示すが、Ｙ方向およびＺ方向についても同様に基本カーネル、一次微分カーネルおよび二次微分カーネルを求めることができる。ここで、基本カーネルが図３の実線に示すものである場合、一次微分カーネルは図３の破線に、二次微分カーネルは図３の一点鎖線に示すものとなる。

ここで、例えばヘッセ行列における二階偏微分値である要素Ｉｘｘは、３次元多重解像度画像Ｍｓｉの処理の対象となる画素におけるＸ方向に二次微分カーネルを、Ｙ方向およびＺ方向のそれぞれに基本カーネルを畳み込む、すなわちフィルタリングすることにより算出できる。また、ヘッセ行列における要素Ｉｘｙは、処理の対象となる画素におけるＸ方向およびＹ方向のそれぞれに一次微分カーネルを、Ｚ方向に基本カーネルを畳み込むことにより算出できる。

このようにして算出したヘッセ行列を固有値分解して固有値を得たとき、線状構造は図４に示すように、３つのうち２つの固有値の値が大きく、１つが０に近い特徴を持つことが知られている。例えば、式（１）の固有値は、線状構造からなる対象組織に対して、式（３）に示す関係を有する。

また、面状構造は図５に示すように、３つのうち１つの固有値の値が大きく、２つが０に近い特徴を持つことが知られている。例えば、式（１）の固有値は、面状構造からなる対象組織に対して、式（４）のような関係を有する。なお、図４，５において、ｅ１，ｅ２，ｅ３は、固有値λ１，λ２，λ３の固有ベクトルの方向を示す。

したがって、固有値から線状構造らしさおよび面状構造らしさを判別することができ、判別結果を用いて、３次元画像Ｍ０において、線状構造物である血管領域、および面状構造物である骨領域を分割することができる。

フィルタリング部３２は、まず、３次元多重解像度画像Ｍｓｉの各画素のＸ，Ｙ，Ｚ方向のそれぞれについての一階偏微分値を算出する。ここで、Ｘ方向の一階偏微分値は、３次元多重解像度画像Ｍｓｉの各画素のＸ方向に一次微分カーネルを、Ｙ方向およびＺ方向に基本カーネルを畳み込むことにより算出する。また、Ｙ方向の一階偏微分値は、３次元多重解像度画像Ｍｓｉの各画素のＹ方向に一次微分カーネルを、Ｘ方向およびＺ方向に基本カーネルを畳み込むことにより算出する。また、Ｚ方向の一階偏微分値は、３次元多重解像度画像Ｍｓｉの各画素のＺ方向に一次微分カーネルを、Ｘ方向およびＹ方向に基本カーネルを畳み込むことにより算出する。算出したＸ方向、Ｙ方向およびＺ方向の一階偏微分値をそれぞれρｘ，ρｙ，ρｚとする。

また、フィルタリング部３２は、３次元多重解像度画像Ｍｓｉの各画素のＸ，Ｙ，Ｚ方向のそれぞれについての二階偏微分値を算出することにより、ヘッセ行列の要素を算出する。ここで上述したようにヘッセ行列における要素Ｉｘｘは、３次元多重解像度画像Ｍｓｉの各画素のＸ方向に二次微分カーネルを、Ｙ方向およびＺ方向のそれぞれに基本カーネルを畳み込むことにより算出する。また、ヘッセ行列における要素Ｉｘｙは、３次元多重解像度画像Ｍｓｉの各画素のＸ方向およびＹ方向のそれぞれに一次微分カーネルを、Ｚ方向に基本カーネルを畳み込むことにより算出する。

評価部３４は、フィルタリング部３２が算出したヘッセ行列を固有値分解し、３つの固有値λ１，λ２，λ３を算出する。なお、固有値λ１，λ２，λ３は、｜λ１｜≦｜λ２｜≦｜λ３｜であるものとする。そして、下記の式（５）、（６）に示すように、３次元多重解像度画像Ｍｓｉの各画素における線状構造らしさの評価値Ｌ０（Lineness）および面状構造らしさの評価値Ｐ０（Planeness）を算出する。

なお、式（５）、（６）におけるａ〜ｈは定数である。また、Ｒ_A，Ｒ_B，Ｒ_Cは下記の式（７）〜（９）により算出する。また、Ｓ_2ndおよびＳ_1stは二階偏微分値および一階偏微分値のパワーであり、下記の式（１０）、（１１）により算出する。また、式（５）、（６）において使用する一階偏微分値ρｘ，ρｙ，ρｚは、二階偏微分値よりもフィルタサイズが大きい、すなわち低解像度の３次元多重解像度画像Ｍｓｉを用いて算出した値を用いる。具体的には、ガウシアンピラミッドにおいて一段階程度解像度が低い３次元多重解像度画像Ｍｓｉを用いて算出した一階偏微分値ρｘ，ρｙ，ρｚを用いるものとする。これにより、３次元多重解像度画像Ｍｓｉにおける輝度値の１次元的な偏りを捉え易くなるため、より精度良く線状構造らしさおよび面状構造らしさの評価値Ｌ０，Ｐ０を算出することができる。

ここで、医用画像に含まれる血管等の線状構造物をヘッセ行列を用いて判別する際の骨構造の誤抽出、および皮質骨等の面状構造物を判別する際の血管の誤抽出においては、誤抽出した部分の輝度値のパターンに１次元的な偏りがある。このように輝度値のパターンに１次元的な偏りがあると、一階偏微分値が大きくなる。一方、理想的な線状構造および面状構造は、その中心を基準として対称な構造をなしているため一階偏微分値は０となる。このため、一階偏微分値の値、すなわちパワーＳ_1stが大きい場合には、線状構造らしさおよび面状構造らしさは小さくなる。したがって、式（５）、（６）に示す最終項exp（−Ｓ_1st／２ｄ²）、exp（−Ｓ_1st／２ｈ²）の値が小さいほど、線状構造らしさおよび面状構造らしさの評価値Ｌ０，Ｐ０は小さくなる。

なお、本実施形態においては、解像度が異なる多重解像度３次元画像Ｍｓｉにおいて、それぞれ線状構造らしさの評価値Ｌ０および面状構造らしさの評価値Ｐ０が算出される。算出された評価値Ｌ０，Ｐ０は、元の３次元画像Ｍ０の対応する画素位置の評価値となるが、すべての多重解像度３次元画像Ｍｓｉの対応する画素位置において評価値が算出されることとなる。このため、本実施形態においては、対応する画素位置において算出された評価値のうち、最大値となる評価値を、元の３次元画像Ｍ０における各画素位置の評価値とする。

分割部４０は、判別部３０が算出した線状構造らしさの評価値Ｌ０および面状構造らしさの評価値Ｐ０に基づいて、３次元画像Ｍ０の血管領域および骨を含む血管以外の領域を領域分割する。具体的には、血管領域を対象領域、血管領域以外の領域を背景領域に設定し、３次元画像Ｍ０内の全画素位置において所定画素サイズの判別領域を設定し、Graph Cut領域分割方法を用いて、判別領域を対象領域と背景領域とに分割する。なお、Graph Cut領域分割方法は、「Yuri Y. Boykov, Marie-Pierre Jolly, “Interactive Graph Cuts for Optimal Boundary and Region Segmentation of Objects in N-D images”, Proceedings of “International Conference on Computer Vision”, Vancouver, Canada, July 2001 vol.I, p.105-112.」に記載されている。

Graph Cut領域分割方法においては、まず、図６に示すように、判別領域内の各画素を表すノードＮ_ij、対象領域（血管領域）を表すノードＳ、背景領域を表すノードＴ、隣接する画素のノード同士が同一の領域であることの確からしさを表すＮ−ｌｉｎｋ、各画素を表すノードＮ_ijと対象領域を表すノードＳとをつなぐＳ−ｌｉｎｋ、およびノードＮ_ijと背景領域を表すノードＴとをつなぐリンクであるＴ−ｌｉｎｋから構成されるグラフを作成する。なお、説明を簡単にするために、図６においては判別領域を３×３の２次元の領域としている。

ここで、各画素を表すノードＮ_ijと対象領域を表すノードＳとをつなぐＳ−ｌｉｎｋは、各画素が対象領域に含まれる画素である確からしさをリンクの太さ（値の大きさ）で表すものであり、各画素を表すノードＮ_ijと背景領域を表すノードＴとをつなぐＴ−ｌｉｎｋは、各画素が背景領域に含まれる画素である確からしさをリンクの太さ（値の大きさ）で表すものである。Ｓ−ｌｉｎｋは、判別部３０が算出した線状構造らしさの評価値Ｌ０が大きいほど大きい値を設定する（リンクを太くする）。Ｔ−ｌｉｎｋは、判別部３０が算出した面状構造らしさの評価値Ｐ０が大きいほど大きい値を設定する（リンクを太くする）。本実施形態においては、式（５）、（６）に示す最終項exp（−Ｓ_1st／２ｄ²）、exp（−Ｓ_1st／２ｈ²）の値が小さいほど、線状構造らしさおよび面状構造らしさの評価値Ｌ０，Ｐ０は小さくなるため、対象領域のＴ−ｌｉｎｋが太くなったり、背景領域のＳ−ｌｉｎｋが太くなったりすることがなくなる。

なお、Ｓ−ｌｉｎｋおよびＴ−ｌｉｎｋの値は、評価値Ｌ０，Ｐ０のみならず画素の輝度値を用いて下記の式（１２）、（１３）に示すように設定してもよい。なお、式（１２）において、ｇ１（）は線状構造らしさの評価値Ｌ０が大きいほど値が大きくなる関数、ｇ２（）は３次元画像Ｍ０の各画素の輝度値（ＣＴ値）Ｉが血管の輝度値に統計的に近いほどＳ−ｌｉｎｋの値を大きくする関数である。また、式（１３）において、ｇ３（）は面状構造らしさの評価値Ｐ０が大きいほど値が大きくなる関数、ｇ４（）は輝度値Ｉが骨の輝度値に統計的に近いほどＴ−ｌｉｎｋの値を大きくする関数である。

Ｓ−ｌｉｎｋ＝ｇ１（Ｌ０）×ｇ２（Ｉ） （１２）
Ｔ−ｌｉｎｋ＝ｇ３（Ｐ０）×ｇ４（Ｉ） （１３）
また、Ｎ−ｌｉｎｋは、隣接する画素が同一領域の画素である確からしさをリンクの太さ（値の大きさ）で表すものであり、式（５）の最終項exp（−Ｓ_1st／２ｄ²）を除いた線状構造らしさの評価値Ｌ０′が大きいほど大きい値を設定する。なお、式（５）の最終項exp（−Ｓ_1st／２ｄ²）の重み付けを小さくして評価値Ｌ０′を算出してもよい。ここで、本実施形態が対象とする血管領域は、Ｎ−ｌｉｎｋにより血管領域全体が同一領域であることを反映できるため、血管領域内に設定されたＳ−ｌｉｎｋに従い、３次元画像Ｍ０において血管領域全体を他の領域と分割できる。とくに、ノイズを有する画像においては、画素間の輝度差が血管領域内においても大きくなり、血管領域内の類似度が小さくなってＮ−ｌｉｎｋも小さくなるが、式（５）の最終項exp（−Ｓ_1st／２ｄ²）を除いた線状構造らしさの評価値（Lineness）を用いることにより、Ｎ−ｌｉｎｋの値を大きくすることができる。

なお、Ｎ−ｌｉｎｋの値についても評価値Ｌ０′のみならず、Ｎ−ｌｉｎｋを設定する２つの画素の輝度値を用いて下記の式（１４）に示すように設定してもよい。なお、式（１４）において、ｇ５（）は評価値Ｌ０′が大きいほど値が大きくなる関数、ｇ６（）は２つの画素の輝度値Ｉ１，Ｉ２の平均値が血管の輝度値および骨の輝度値の平均値に統計的に近いほどＮ−ｌｉｎｋの値を大きくする関数、ｇ７（）は２つの画素の輝度値Ｉ１，Ｉ２の差が小さいほどＮ−ｌｉｎｋの値を大きくする関数である。

Ｎ−ｌｉｎｋ＝ｇ５（Ｌ０′）×ｇ６（（Ｉ１＋Ｉ２）／２）×ｇ７（Ｉ１−Ｉ２） （１４）
ここで、図６においては、ノードＮ₁₁，Ｎ₁₂，Ｎ₂₁，Ｎ₂₂，Ｎ₃₁により表される画素においては、線状構造らしさの評価値Ｌ０が大きいため、各ノードＮ₁₁，Ｎ₁₂，Ｎ₂₁，Ｎ₂₂，Ｎ₃₁とノードＳとをつなぐＳ−ｌｉｎｋは太くなる。また、ノードＮ₁₁，Ｎ₁₂，Ｎ₂₁，Ｎ₂₂，Ｎ₃₁により表される画素においては、式（５）の最終項exp（−Ｓ_1st／２ｄ²）を除いた線状構造らしさの評価値Ｌ０′が大きいため、各ノードＮ₁₁，Ｎ₁₂，Ｎ₂₁，Ｎ₂₂，Ｎ₃₁をつなぐＮ−ｌｉｎｋは太くなる。一方、また、ノードＮ₁₃，Ｎ₂₃，Ｎ₃₂，Ｎ₃₃により表される画素においては、面状構造らしさの評価値Ｐ０が大きいため、各ノードＮ₁₃，Ｎ₂₃，Ｎ₃₂，Ｎ₃₃とノードＴとをつなぐＴ−ｌｉｎｋは太くなる。

そして、血管領域と血管以外の骨を含む領域とは互いに排他的な領域であるので、例えば図６に破線で示すように、Ｓ−ｌｉｎｋ、Ｔ−ｌｉｎｋおよびＮ−ｌｉｎｋのうち適当なリンクを切断してノードＳをノードＴから切り離すことにより、判別領域を対象領域と背景領域とに分割できる。ここで、切断するすべてのＳ−ｌｉｎｋ、Ｔ−ｌｉｎｋおよびＮ−ｌｉｎｋの値の合計が最も小さくなるような切断を行うことにより、最適な領域分割をすることができる。

表示部５０は、２次元画像または３次元画像等を表示するモニタ、ＣＲＴ画面等である。本実施形態においては、表示部５０に、対象領域として分割された線状構造および面状構造をボリュームレンダリング表示することにより線状構造や面状構造の全体を概観し、その連続性を可視化することができる。

入力部６０は、キーボード、マウス等である。

次いで、本実施形態において行われる処理について説明する。図７は本実施形態において行われる処理を示すフローチャートである。まず、画像取得部１０が、Ｘ線ＣＴ装置２により撮像された２次元画像から３次元画像Ｍ０を生成する（ステップＳＴ１）。次に、検出領域設定部２０が、３次元画像Ｍ０を等方化するとともに多重解像度変換して解像度が異なる複数の３次元多重解像度画像Ｍｓｉを生成する（ステップＳＴ２）。

次いで、判別部３０のフィルタリング部３２が、３次元多重解像度画像Ｍｓｉ毎に、ガウシアンカーネルを用いたフィルタリングを行い、一階偏微分値および二階偏微分行列を算出する（ステップＳＴ３）。次いで、評価部３４が、線状構造らしさの評価値Ｌ０および面状構造らしさの評価値Ｐ０を算出する（ステップＳＴ４）。

そして、分割部４０が上述したGraph Cut領域分割方法を用いて、３次元画像Ｍ０を対象領域（血管領域）と背景領域とに分割する（ステップＳＴ５）。そして、表示部５０が分割された対象領域および背景領域をボリュームレンダリング表示し（ステップＳＴ６）、処理を終了する。

このように、本実施形態においては、３次元画像Ｍ０の各画素位置における画素値の二階偏微分行列および一階偏微分値を算出し、二階偏微分行列の固有値に基づいて、画素位置における線状構造らしさおよび面状構造らしさの評価値Ｌ０，Ｐ０を算出するに際して、一階偏微分値が大きいほど評価値Ｌ０，Ｐ０を小さくするようにしたものである。このため、一階偏微分値が大きく、理想的な線状構造および面状構造から離れた形状となるほど、評価値Ｌ０，Ｐ０が小さくなることから、本実施形態によれば、３次元画像Ｍ０に含まれる線状構造物と面状構造物との誤判別を防止して、３次元画像Ｍ０に含まれる血管等の線状構造物と骨等の面状構造物とを精度良く判別できる。

また、一階偏微分値は、二階偏微分行列を算出する際に必ず算出される一次微分カーネルを用いて算出されるため、従来のヘッセ解析を行う場合と比較して演算量が増大することもない。

また、評価値Ｌ０，Ｐ０に基づいて、３次元画像Ｍ０の各画素に対して、対象領域らしさを表すＳ−ｌｉｎｋ、背景領域らしさを表すＴ−ｌｉｎｋ、および隣接する画素が同一の領域であることの確からしさを表すＮ−ｌｉｎｋを設定して、３次元画像Ｍ０を対象領域および背景領域に分割することにより、対象領域および背景領域を精度良く分割できる。

なお、上記実施形態においては、３次元画像Ｍ０に含まれる線状構造物および面状構造物を判別しているが、２次元画像に含まれる線状構造物および面状構造物を判別するようにしてもよい。

また、上記実施形態においては、線状構造物として血管の判別を例に挙げたが、気管支等の他の線状構造物の判別にも適用できる。また、また骨のみならず、皮膚、葉間胸膜等の面状構造物の判別にも利用することができる。

また、上記実施形態においては、Graph Cut領域分割方法を用いて３次元画像Ｍ０に含まれる線状構造物と面状構造物とを分割しているが、Watershedアルゴリズム等の他の領域分割手法を用いてもよいことはもちろんである。ここで、Watershedアルゴリズムは、画像の画素値情報を高度とみなした地形において水を満たしていく際に、異なるくぼみに溜まる水の間で境界ができるように画像を分割する手法である。このため、３次元画像Ｍ０上の評価値Ｌ０，Ｐ０に対し、適当な平滑化を行った後でWatershedアルゴリズムを実行することによって、線状構造物および面状構造物を分割することが可能である。

また、上記実施形態においては、３次元画像Ｍ０を多重解像度変換して、複数の３次元多重解像度画像Ｍｓｉを作成し、一つのサイズのフィルタを用いて一階偏微分値および二階偏微分行列を算出することにより、各種サイズの線状構造物および面状構造物を判別しているが、３次元画像Ｍ０に対して、サイズが異なる複数のフィルタを用いて一階偏微分値および二階偏微分行列を算出するようにしてもよい。

また、上記実施形態においては、線状構造らしさの評価値Ｌ０および面状構造らしさの評価値Ｐ０の双方を算出しているが、線状構造らしさの評価値Ｌ０および面状構造らしさの評価値Ｐ０のいずれか一方のみを算出するようにしてもよい。

また、上記実施形態においては、ヘッセ行列を固有値分解することにより取得した固有値を用いて線状構造らしさの評価値Ｌ０および面状構造らしさの評価値Ｐ０を算出しているが、固有値を得ることなく、ヘッセ行列の要素をそのまま用いて線状構造らしさの評価値および面状構造らしさの評価値を算出するようにしてもよい。

ところで、理想的な線状構造においては、上述したヘッセ行列を固有値分解することにより得られる３つの固有値λ１，λ２，λ３のうち、線状構造が延在する方向に直交する方向の２つの固有値λ２，λ３は、式（３）に示すように略等しくなる。しかしながら、例えば図１０に示す心臓１１０の辺縁を走行する冠動脈１１１のように、面状構造物（心臓１１０の壁面）の近傍に線状構造物（冠動脈１１１）が存在する場合、線状構造物から面状構造物の方向、すなわち面状構造物の法線方向における固有値が大きくなる。例えば、面状構造物（心臓１１０の壁面）の法線方向が固有ベクトルｅ３の方向と一致する場合、その方向の固有値λ３は、これに直交する方向の固有値λ２よりも大きくなる。これは、心臓の内側に対して外側は肺が存在することから、心臓の内側と外側とでＣＴ値の差が大きくなるためである。このため、線状構造らしさの評価値Ｌ０が低くなり、その結果、線状構造物の判別に失敗するおそれがある。

以下、この問題を解決するための手法について本発明の第２の実施形態として説明する。上述したように算出したＸ方向、Ｙ方向およびＺ方向の一階偏微分値ρｘ，ρｙ，ρｚは、固有値λ１，λ２，λ３の固有ベクトルｅ１，ｅ２，ｅ３の方向とずれているため、第２の実施形態においては、評価部３４は、固有ベクトルｅ１，ｅ２，ｅ３の方向に対応した一階偏微分値ρ１，ρ２，ρ３を下記の式（１５）により算出する。

ここで、固有値λ２，λ３のうち、値が大きい方の固有値をλ３とすると、評価部３４は、固有値λ３を下記の式（１５）により補正する。すなわち、固有ベクトルｅ２，ｅ３の方向における一階偏微分値ρ２，ρ３の大きさに応じて、固有値λ３を補正して、補正された固有値λ３′を得る。

λ３′＝λ２＋（λ３−λ２）×ｆ（｜ρ２−ρ３｜） （１６）
式（１６）において、ｆ（）は０〜１の値を出力し、｜ρ２−ρ３｜の値が大きいほど小さい値を出力する関数である。これにより、｜ρ２−ρ３｜の値が大きいほど、固有値λ３′は固有値λ２の値に近い値となるため、式（３）に示した固有値λ２と固有値λ３とが略等しくなるという関係を満たすこととなる。

そして、評価部３４は、Ｒ_A，Ｒ_B，Ｒ_Cの各値を算出する式（７）〜（９）において、固有値λ３に代えて固有値λ３′を使用し、これにより算出した値Ｒ_A，Ｒ_B，Ｒ_Cを用いて式（５）により線状構造らしさの評価値Ｌ０を算出する。

このように、第２の実施形態においては、線状構造らしさの評価値Ｌ０を算出するに際し、固有値λ２，λ３の固有ベクトルｅ２，ｅ３の方向に一致する一階偏微分値ρ２，ρ３の大きさに応じて固有値λ３を補正し、補正された固有値λ３′に基づいて評価値Ｌ０を算出するようにしたものである。このため、面状構造物の近傍に線状構造物が存在する場合であっても、面状構造物が存在する方向およびこれに直交する方向における固有値λ２，λ３を略等しくすることができ、これにより、線状構造らしさの評価値Ｌ０が低くなることを防止でき、その結果、線状構造物を精度よく判別することができることとなる。

１ 画像処理装置
２ Ｘ線ＣＴ装置
１０ 画像取得部
２０ 検出領域設定部
３０ 判別部
３２ フィルタリング部
３４ 評価部
４０ 分割部
５０ 表示部
６０ 入力部

Claims (19)

1. 画像中の各画素位置における画素値について、前記画像の次元に応じた複数方向のそれぞれにおける一階偏微分値および二階偏微分値を算出し、かつ該二階偏微分値を要素とする二階偏微分行列を算出する微分値算出手段と、
   前記二階偏微分行列を固有値分解することにより前記画像の次元に応じた固有値を算出し、該固有値に基づいて、前記画素位置における線状構造らしさの評価値および／または面状構造らしさの評価値を算出する評価手段とを備え、
   前記評価手段は、前記複数の一階偏微分値のパワーを算出し、該パワーが大きいほど前記線状構造らしさの評価値および／または面状構造らしさの評価値を小さくするよう該各評価値を算出する手段であることを特徴とする画像処理装置。
2. 前記微分値算出手段は、複数サイズのフィルタを用いて前記一階偏微分値および前記二階偏微分値を算出するに際し、該二階偏微分値を算出するフィルタのサイズよりも大きいサイズのフィルタにより、前記一階偏微分値を算出する手段であることを特徴とする請求項１記載の画像処理装置。
3. 前記微分値算出手段は、前記画像を多重解像度変換して解像度が異なる複数の解像度画像を取得し、該複数の解像度画像の対応する画素位置において所定サイズのフィルタを用いて前記一階偏微分値および前記二階偏微分値を算出するに際し、該二階偏微分値を算出する解像度画像よりも低い解像度の解像度画像により、前記一階偏微分値を算出する手段であることを特徴とする請求項１記載の画像処理装置。
4. 前記微分値算出手段は、１次元の基本ガウシアンカーネル、該基本ガウシアンカーネルを一階微分した一次微分カーネル、および該基本ガウシアンカーネルを二階微分した二次微分カーネルを用いて、前記一階偏微分値および前記二階偏微分値を算出する手段であることを特徴とする請求項１から３のいずれか１項記載の画像処理装置。
5. 前記評価値に基づいて、前記画像の前記各画素に対して、対象領域らしさ、背景領域らしさ、および隣接する画素が同一の領域であることの確からしさを設定して、該対象領域および該背景領域を分割する分割手段をさらに備えたことを特徴とする請求項１から４のいずれか１項記載の画像処理装置。
6. 前記分割手段は、前記対象領域らしさを前記線状構造らしさの評価値に基づいて設定し、前記背景領域らしさを前記面状構造らしさの評価値に基づいて設定する手段であることを特徴とする請求項５記載の画像処理装置。
7. 前記分割手段は、前記同一の領域であることの確からしさを、前記固有値にのみ基づいて算出された、前記画素位置における線状構造らしさの評価値および面状構造らしさの評価値、または前記線状構造らしさの評価値および前記面状構造らしさの評価値における、前記一階偏微分値の影響を小さくした評価値に基づいて設定する手段であることを特徴とする請求項５または６記載の画像処理装置。
8. 画像中の各画素位置における画素値について、前記画像の次元に応じた複数方向のそれぞれにおける一階偏微分値および二階偏微分値を算出し、かつ該二階偏微分値を要素とする二階偏微分行列を算出し、
   前記二階偏微分行列を固有値分解することにより前記画像の次元に応じた固有値を算出し、該固有値に基づいて、前記画素位置における線状構造らしさの評価値および／または面状構造らしさの評価値を算出するに際し、前記複数の一階偏微分値のパワーを算出し、該パワーが大きいほど前記線状構造らしさの評価値および／または面状構造らしさの評価値を小さくするよう該各評価値を算出することを特徴とする画像処理方法。
9. 画像中の各画素位置における画素値について、前記画像の次元に応じた複数方向のそれぞれにおける一階偏微分値および二階偏微分値を算出し、かつ該二階偏微分値を要素とする二階偏微分行列を算出する手順と、
   前記二階偏微分行列を固有値分解することにより前記画像の次元に応じた固有値を算出し、該固有値に基づいて、前記画素位置における線状構造らしさの評価値および／または面状構造らしさの評価値を算出するに際し、前記複数の一階偏微分値のパワーを算出し、該パワーが大きいほど前記線状構造らしさの評価値および／または面状構造らしさの評価値を小さくするよう該各評価値を算出する手順とを有することを特徴とする画像処理方法をコンピュータに実行させるためのプログラム。
10. 画像中の各画素位置における画素値について、前記画像の次元に応じた複数方向のそれぞれにおける一階偏微分値および二階偏微分値を算出し、かつ該二階偏微分値を要素とする二階偏微分行列を算出する微分値算出手段と、
    前記二階偏微分行列を固有値分解することにより前記画像の次元に応じた固有値を算出し、該固有値に基づいて、前記画素位置における線状構造らしさの評価値を算出する評価手段とを備え、
    前記評価手段は、前記評価手段は、前記複数の一階偏微分値のパワーを算出し、該パワーが大きいほど前記評価値を小さくするよう該評価値を算出する手段であることを特徴とする画像処理装置。
11. 前記評価手段は、前記二階偏微分方向に一致する前記一階偏微分値の大きさに応じて前記固有値を補正し、該補正された固有値に基づいて、前記評価値を算出する手段であることを特徴とする請求項１０記載の画像処理装置。
12. 前記微分値算出手段は、複数サイズのフィルタを用いて前記一階偏微分値および前記二階偏微分値を算出するに際し、該二階偏微分値を算出するフィルタのサイズよりも大きいサイズのフィルタにより、前記一階偏微分値を算出する手段であることを特徴とする請求項１０または１１記載の画像処理装置。
13. 前記微分値算出手段は、前記画像を多重解像度変換して解像度が異なる複数の解像度画像を取得し、該複数の解像度画像の対応する画素位置において所定サイズのフィルタを用いて前記一階偏微分値および前記二階偏微分値を算出するに際し、該二階偏微分値を算出する解像度画像よりも低い解像度の解像度画像により、前記一階偏微分値を算出する手段であることを特徴とする請求項１０または１１記載の画像処理装置。
14. 前記微分値算出手段は、１次元の基本ガウシアンカーネル、該基本ガウシアンカーネルを一階微分した一次微分カーネル、および該基本ガウシアンカーネルを二階微分した二次微分カーネルを用いて、前記一階偏微分値および前記二階偏微分値を算出する手段であることを特徴とする請求項１０から１３のいずれか１項記載の画像処理装置。
15. 前記評価値に基づいて、前記画像の前記各画素に対して、対象領域らしさ、背景領域らしさ、および隣接する画素が同一の領域であることの確からしさを設定して、該対象領域および該背景領域を分割する分割手段をさらに備えたことを特徴とする請求項１０から１４のいずれか１項記載の画像処理装置。
16. 前記分割手段は、前記対象領域らしさを前記線状構造らしさの評価値に基づいて設定する手段であることを特徴とする請求項１５記載の画像処理装置。
17. 前記分割手段は、前記同一の領域であることの確からしさを、前記固有値にのみ基づいて算出された、前記画素位置における線状構造らしさの評価値、または前記線状構造らしさの評価値における、前記一階偏微分値の影響を小さくした評価値に基づいて設定する手段であることを特徴とする請求項１５または１６記載の画像処理装置。
18. 画像中の各画素位置における画素値について、前記画像の次元に応じた複数方向のそれぞれにおける一階偏微分値および二階偏微分値を算出し、かつ該二階偏微分値を要素とする二階偏微分行列を算出し、
    前記二階偏微分行列を固有値分解することにより前記画像の次元に応じた固有値を算出し、該固有値に基づいて、前記画素位置における線状構造らしさの評価値を算出するに際し、前記複数の一階偏微分値のパワーを算出し、該パワーが大きいほど前記評価値を小さくするよう該評価値を算出することを特徴とする画像処理方法。
19. 画像中の各画素位置における画素値について、前記画像の次元に応じた複数方向のそれぞれにおける一階偏微分値および二階偏微分値を算出し、かつ該二階偏微分値を要素とする二階偏微分行列を算出する手順と、
    前記二階偏微分行列を固有値分解することにより前記画像の次元に応じた固有値を算出し、該固有値に基づいて、前記画素位置における線状構造らしさの評価値を算出するに際し、前記複数の一階偏微分値のパワーを算出し、該パワーが大きいほど前記評価値を小さくするよう該評価値を算出する手順とを有することを特徴とする画像処理方法をコンピュータに実行させるためのプログラム。

Images