以下、本発明の実施の形態に係る画像処理装置、画像処理方法、及び画像処理プログラムについて、図面を参照しながら説明する。なお、これら実施の形態によって本発明が限定されるものではない。また、各図面の記載において、同一部分には同一の符号を付して示している。
(実施の形態1)
図1は、本発明の実施の形態1に係る画像処理装置を示すブロック図である。本実施の形態1に係る画像処理装置1は、一例として、カプセル型内視鏡によって生体の管腔内を撮像することにより取得された管腔内画像(以下、単に画像ともいう)に対し、粘膜表面の微細構造である柔毛が膨張した病変部(異常部)の有無を判別する画像処理を施す装置である。管腔内画像は、通常、各画素位置においてR(赤)、G(緑)、B(青)の波長成分(色成分)に対して、所定の(例えば256階調の)画素レベル(画素値)を有するカラー画像である。
図1に示すように、画像処理装置1は、該画像処理装置1全体の動作を制御する制御部10と、内視鏡によって撮像された画像に対応する画像データを取得する画像取得部20と、外部から入力される入力信号を受け付ける入力部30と、各種表示を行う表示部40と、画像取得部20によって取得された画像データや種々のプログラムを格納する記録部50と、画像データに対して所定の画像処理を実行する演算部100とを備える。
制御部10は、CPU等のハードウェアによって実現され、記録部50に記録された各種プログラムを読み込むことにより、画像取得部20から入力される画像データや入力部30から入力される操作信号等に従って、画像処理装置1を構成する各部への指示やデータの転送等を行い、画像処理装置1全体の動作を統括的に制御する。
画像取得部20は、内視鏡を含むシステムの態様に応じて適宜構成される。例えば、カプセル型内視鏡との間の画像データの受け渡しに可搬型の記録媒体が使用される場合、画像取得部20は、この記録媒体を着脱自在に装着し、記録された画像の画像データを読み出すリーダ装置で構成される。また、内視鏡によって撮像された画像の画像データを保存しておくサーバを設置する場合、画像取得部20は、サーバと接続される通信装置等で構成され、サーバとデータ通信を行って画像データを取得する。或いは、画像取得部20を、内視鏡からケーブルを介して画像信号を入力するインターフェース装置等で構成しても良い。
入力部30は、例えばキーボードやマウス、タッチパネル、各種スイッチ等の入力デバイスによって実現され、受け付けた入力信号を制御部10に出力する。
表示部40は、LCDやELディスプレイ等の表示装置によって実現され、制御部10の制御の下で、管腔内画像を含む各種画面を表示する。
記録部50は、更新記録可能なフラッシュメモリ等のROMやRAMといった各種ICメモリ、内蔵若しくはデータ通信端子で接続されたハードディスク、又は、CD−ROM等の情報記録媒体及びその読取装置等によって実現される。記録部50は、画像取得部20によって取得された画像データの他、画像処理装置1を動作させると共に、種々の機能を画像処理装置1に実行させるためのプログラムや、このプログラムの実行中に使用されるデータ等を格納する。具体的には、記録部50は、粘膜表面における柔毛の異常部の有無を判別する画像処理プログラム51や、このプログラムの実行中に使用される種々の情報等を格納する。
演算部100は、CPU等のハードウェアによって実現され、画像処理プログラム51を読み込むことにより、管腔内画像に対して画像処理を施し、粘膜表面における柔毛の異常部の有無を判別するための種々の演算処理を実行する。
次に、演算部100の構成について説明する。
図1に示すように、演算部100は、処理対象の管腔内画像に注目領域を設定する注目領域設定部110と、周囲よりも画素値が高い画素が所定の画素数以上連続する線形状の領域(以下、線形状領域ともいう)を抽出する線形状凸領域抽出部120と、線形状領域に沿った1つ以上の円弧の曲率に基づく特徴量(以下、曲率特徴量ともいう)の分布を算出する領域内曲率特徴量算出部130と、曲率特徴量に基づき注目領域内に異常部が存在するか否かを判別する異常判別部140を備える。
このうち、注目領域設定部110は、画像から残渣や暗部といった粘膜以外の領域を除外して粘膜領域を抽出する粘膜領域抽出部111を備え、抽出された粘膜領域を注目領域に設定する。
また、線形状凸領域抽出部120は、凸形状高周波成分算出部121と、孤立点除外部122と、細線化部123とを備える。
凸形状高周波成分算出部121は、周囲よりも画素値が高い画素領域に対し、空間周波数が所定値以上である成分(以下、高周波成分という)の強度を算出する。なお、以下において、周囲よりも画素値が高い画素領域のことを、凸形状領域又は単に凸形状ともいい、凸形状であって、且つ高周波成分の強度が所定値以上である画素領域のことを、凸形状高周波領域ともいう。また、凸形状高周波領域における高周波成分の強度のことを、凸形状高周波成分ともいう。
孤立点除外部122は、凸形状高周波成分算出部121により算出された凸形状高周波成分に基づき、孤立点を正常な柔毛として除外する。ここで、孤立点とは、凸形状高周波領域のうち、該凸形状高周波領域の周囲の全方向において、連続する画素数が所定の閾値よりも少ない領域のことである。
細線化部123は、凸形状高周波領域を細線化する。
また、領域内曲率特徴量算出部130は、線形状領域に沿った1つ以上の円弧の曲率及び該円弧までの撮像距離に対応する距離情報を、円弧のサイズ特徴量として算出するサイズ特徴量算出部131と、円弧のサイズ特徴量の頻度分布を作成する頻度分布作成部132とを備える。このうち、サイズ特徴量算出部131は、線形状領域から1つ以上の円弧の曲率を算出する曲率算出部131aと、算出された1つ以上の円弧の曲率から代表値を算出する曲率代表値算出部131bと、管腔内画像の撮像位置(即ち、カプセル型内視鏡の位置)から線形状をなす凸領域までの距離情報を算出する距離情報算出部131cとを備える。
次に、画像処理装置1の動作について説明する。図2は、画像処理装置1の動作を示すフローチャートである。
まず、ステップS10において、演算部100は、記録部50に記録された画像データを読み込むことにより、処理対象の管腔内画像を取得する。
続くステップS20において、注目領域設定部110は、管腔内画像に対して注目領域を設定する。図3は、注目領域設定部110が実行する処理を詳細に示すフローチャートである。
ステップS201において、粘膜領域抽出部111は、管腔内画像内の各画素のエッジ強度に基づいて、管腔内画像を複数の小領域に分割する。この管腔内画像を複数の小領域に分割する処理を、図4を参照しながら詳しく説明する。まず、ステップS2011において、粘膜領域抽出部111は、管腔内画像内の各画素の画素値のうちのG成分を用いて、G成分画像を作成する。ここでG成分を用いる理由は、G成分は血液の吸光帯域に近いことから、粘膜の隆起や残渣の境界といった被写体の構造を最もよく表すからである。なお、本実施の形態1においてはG成分を用いるが、その代わりに、他の色成分や、画素値(R、G、B各成分)を変換して得られる輝度値、色差(YCbCr変換)、色相、彩度、明度(HSI変換)、色比等を用いても良い。
ステップS2012において、粘膜領域抽出部111は、このG成分画像にフィルタ処理を施してエッジ強度を算出し、エッジ強度を各画素の画素値とするエッジ画像を作成する。エッジ強度を算出するためのフィルタ処理としては、ブリューウイットフィルタ、ソーベルフィルタといった1次微分フィルタや、ラプラシアンフィルタ、LOG(Laplacian of Gaussian)フィルタといった2次微分フィルタ等が用いられる(参考:CG-ARTS協会、「ディジタル画像処理」、第114〜121頁)。
さらに、粘膜領域抽出部111は、エッジ画像におけるエッジ強度をもとに、領域分割を行う。領域分割の手法の一例として、本実施の形態1においては、国際公開第2006/080239号に開示された方法を適用する。即ち、ステップS2013において、エッジ画像に対し、ノイズを除去するために平滑化処理を施す。続くステップS2014において、エッジ画像内の各画素位置において、エッジ強度の勾配が最大となる方向(最大勾配方向)を取得する。そして、ステップS2015において、各画素(出発画素)から最大勾配方向に沿ってエッジ画像を探索し、画素値(エッジ強度)が極値となる画素位置を取得する。さらに、ステップS2016において、同一の画素位置において極値となる出発画素同士を1つの領域としてラベリングしたラベル画像を作成する。ここで、ラベル画像とは、ラベリングにより得られたラベル番号(1〜n,nはラベル数)を各画素の画素値とする画像のことである。このラベル画像において同一の画素値を有する画素領域が、管腔内画像を領域分割した各小領域に相当する。その後、処理はメインルーチンに戻る。
なお、管腔内画像を領域分割する方法としては、上述した方法以外の公知の方法を適用しても良い。具体的には、分水嶺(watershed)アルゴリズム(参考:Luc Vincent and Pierre Soille,“Watersheds in digital spaces:An efficient algorithm based on immersion simulations”,IEEE Transactions on Pattern Analysis and Machine Intelligence,Vol.13, No.6,pp.583−598,June 1991)を用いても良い。分水嶺アルゴリズムは、画像の画素値情報を高度とみなした地形において水を満たしていく際に、異なるくぼみに溜まる水の間で境界ができるように画像を分割する方法である。
また、本実施の形態1においては、孤立的な画素ノイズによる影響を軽減したり、粘膜や残渣等の境界に沿った領域分割により後段における異常判別を実行し易くする目的で、エッジに基づく領域分割を行ったが、他の領域分割の手法として、管腔内画像を所定サイズの矩形領域に分割しても良い。この場合、エッジ画像の作成及びエッジ強度に基づく領域分割といった処理が不要となり、処理時間を短縮することができる。或いは、領域分割を行わずに、各画素を以下に説明する各ステップにおける処理単位としても良い。この場合、ステップS201の処理を省略することができる。
ステップS201に続くステップS202において、粘膜領域抽出部111は、領域分割により得られた各小領域に対し、特徴量の平均値を算出する。より詳細には、粘膜領域抽出部111は、小領域内の各画素の画素値から、R、G、B各成分の平均値を算出する。そして、これらの平均値に対してHSI変換(参考:CG-ARTS協会、「ディジタル画像処理」、第64頁(HSI変換と逆変換))を施すことにより、小領域の色相平均値を算出する。ここで、色相平均値を算出するのは、粘膜や管腔内の内容物等の被写体は、それらの成分である血液や胆汁等の吸光特性の違いから、色相等の色特徴量によってある程度判別できるからである。基本的に、胃粘膜、腸粘膜、内容物の色相は、この順に、赤系の色相から黄系に向かって変化する。
なお、各小領域の色相平均値は、小領域内の各画素の画素値(R、G、B各成分)に対してHSI変換を施して色相を算出した後、各画素の色相の平均値を算出することにより求めても良い。また、色相の代わりに、例えば、色差、色比(G/R、B/G)等の特徴量の平均値を算出しても良い。
ステップS203において、粘膜領域抽出部111は、各小領域に対し、特徴量の平均値を基に暗部領域であるか否かを判別する。より詳細には、粘膜領域抽出部111は、まず、各小領域の判別結果リストを作成する。ここで、判別結果リストとは、管腔内画像を複数の小領域に分割した際のラベル画像におけるラベル番号と、各ラベル番号に対応するフラグ情報とを関連付けたリストのことである。なお、判別結果リストのサイズは、ラベル数nに対応する。続いて、粘膜領域抽出部111は、判別結果リストを初期化して、全ての小領域に対して粘膜領域のフラグ(0:粘膜領域)を付与する。そして、粘膜の隆起や残渣の境界といった被写体の構造に最もよく表すG成分を特徴量として用い、該特徴量の平均値が所定の閾値以下となる小領域に対して暗部領域のフラグ(1:暗部領域)を付与する。なお、この際の閾値は、明度の変化に対して小領域の色相平均値等の変化が線形性を保つ範囲内で決定される。これは、暗部領域においては、ノイズ等の影響により線形性が崩れるためである。
ステップS204において、粘膜領域抽出部111は、ステップS202において算出した特徴量の平均値(色相平均値)の分布をクラスタリングする。ここで、クラスタリングとは、特徴空間内におけるデータ分布を、データ間の類似性を基にクラスタと呼ばれる塊に分ける手法であり、階層法やk-means法等の公知の様々な手法で実行することができる(参考:CG-ARTS協会:「ディジタル画像処理」、第231〜232頁(クラスタリング))。本実施の形態1においては、後段で色相空間におけるデータのヒストグラムを用いる処理を行うため、特徴空間を色相の1次元空間とし、データのヒストグラムに基づくクラスタリングの手順を示す。なお、クラスタリングは、これに限定する必要はなく、他の手法を用いても良い。
以下、色相平均値の分布をクラスタリングする処理を、図5〜図7を参照しながら詳しく説明する。まず、横軸を色相平均値の階級区間とし、縦軸を各階級区間に相当する色相平均値の度数とするヒストグラムを作成する。この際、色相平均値は、判別結果リストにおいて粘膜領域のフラグ(0:粘膜領域)が付与された小領域から取得する。また、階級区間に関しては、管腔内画像における色相分布範囲を十分に含む黄系から赤系の色相範囲を、事前に定められた所定の間隔で分割することにより設定する。図5は、ヒストグラムの作成例である。なお、図5においては、度数の総和が1になるように正規化されている。
続いて、色相平均値の分布を、ヒストグラムの谷を境界として複数のクラスタに分ける。図6は、クラスタリング方法を説明する模式図である。図6に示す度数データのグラフは、ヒストグラムの度数変化を折れ線で簡易的に示したものである。また、ここでは、ヒストグラムの横軸に示す階級区間を便宜的に座標で示している。なお、クラスタリングされたクラスタ数を制限するために、ヒストグラムの階級区間の幅を変更したり、クラスタリングの前にヒストグラムに平滑化処理を施したりしても良い。
このような度数データに対し、各座標における勾配方向を求める。ここで、勾配方向とは、注目座標における度数データと、注目座標に隣接する座標における度数データとの差に基づいて決まる方向であり、度数の値が増加する方向を矢印で示している。なお「極」と表示された座標は、隣接する何れの座標よりも度数が高い極値(極大値)座標である。
さらに、各座標から勾配方向に沿って極値座標を探索する。図6においては、極値座標を探索した際の座標の変化を、出発座標(n=0における各座標)の位置から下方向に向かって座標の数値で示している。例えば、座標1から探索を開始する場合、座標1においては勾配方向が右方向となっているため、1回目の探索ステップ(n=1)では右隣の座標2が得られる。続いて、座標2においても勾配方向は右方向となっているため、2回目の探索ステップ(n=2)では、座標2の右隣の座標3が得られる。その後も順次、勾配方向に沿って探索を続けると、最終的に極値座標である座標5に辿り着く(n=4)。次に、座標2から探索を開始すると、やはり、最終的には極値座標である座標5の極値座標に辿り着く(n=3)。同様にして、全ての座標から同様の探索を行うと、座標1〜9を出発座標とする場合には、座標5に辿り着き、座標10〜15を出発座標とする場合には、座標12に辿り着く。
このような探索の結果、極値座標として座標5が得られた座標1〜9を第1のクラスタとし、極値座標として座標12が得られた座標10〜15を第2のクラスタとすることにより、ヒストグラムの谷を境界としてデータ分布をクラスタに分けることができる。
図7は、上述したクラスタリング方法を図5に示すヒストグラムに適用した結果を示すグラフである。図5の場合には、ヒストグラムの谷を境界として、色相平均値をクラスタ1〜3の3つに分類することができる。
ステップS204に続くステップS205において、粘膜領域抽出部111は、特徴量平均値の分布を基に残渣等の不要領域を判別する。ここで、管腔内においては一般に、胆汁と血液との吸光特性の違いから、残渣等の不要領域は黄系の色になる。そこで、ステップS204において得られたクラスタのうち、色相が黄系寄りのクラスタを、不要領域である小領域の色相平均値データからなるクラスタ(以下、不要領域クラスタという)と推定し、不要領域クラスタを特定する。
より詳細には、まず、各クラスタa(aはクラスタ番号であり、図7の場合、a=1、2、3)の重心ω
aバー(バーはωの上に平均を示す記号「−」が記されていることを示す)を、クラスタa内の色相平均値の総和ω
a、及びクラスタaに属するデータの数(即ち、小領域の数)m
aを用いて、次式(1)によって算出する。
式(1)において、
また、記号iは小領域のラベル番号を示す。記号nはラベル数、即ち、画像内の小領域の数を示す。記号h
iは、小領域iの色相平均値を示す。記号a_minは、クラスタaにおける色相平均値の最小値を示し、記号a_maxは、クラスタaにおける色相平均値の最大値を示す。さらに、
図7に示す色相平均値のヒストグラムにおいては、色相軸(横軸)の左に行くほど、黄系が強くなる。そこで、予め取得した教師データから、残渣の色相の範囲を示す残渣判定閾値を算出しておき、この残渣判定閾値よりも重心がグラフの左側に位置するクラスタを、不要領域クラスタと判別する。例えば図7においては、クラスタ1が不要領域クラスタとして判別される。
さらに、粘膜領域抽出部111は、判別結果リストにおいて粘膜領域のフラグ(0:粘膜領域)が付与された小領域のうち、不要領域クラスタに属する小領域に対し、不要領域のフラグ(2:不要領域)を付与する。
なお、不要領域の判別にあたっては、色相以外の特徴量の平均値を用いて同様の処理を行っても良い。
続くステップS206において、注目領域設定部110は、ラベル画像及び判別結果リストのフラグ情報をもとに、画像内の粘膜領域(不要領域以外の領域)の座標情報を取得し、該粘膜領域の座標情報を注目領域として出力する。その後、処理はメインルーチンに戻る。
ステップS20に続くステップS30において、線形状凸領域抽出部120は、画像内の注目領域から、周囲よりも画素値が高い画素が所定の画素数以上連続する線形状領域を抽出する。図8は、線形状凸領域抽出部120が実行する処理を詳細に示すフローチャートである。
ステップS301において、凸形状高周波成分算出部121は、凸形状高周波成分の強度を算出する。より詳細には、まず、凸形状高周波成分算出部121は、管腔内画像内の各画素の画素値のうちのG成分を用いて、G成分画像を作成する。G成分を用いる理由は、上述したように、G成分は血液の吸光帯域に近いため、被写体の構造を最もよく表すからである。なお、本実施の形態1においてはG成分を用いるが、その代わりに、他の色成分や、画素値(R、G、B各成分)を変換して得られる輝度値、色差(YCbCr変換)、色相、彩度、明度(HSI変換)、色比等を用いても良い。
続いて、凸形状高周波成分算出部121は、G成分画像に対して濃淡モフォロジ(モルフォロジー)によるトップハット変換(参考:コロナ社、「モルフォロジー」、第103〜104頁)を実行する。ここで、トップハット変換とは、図9に示すように、元の信号出力y=f(x)に対してローパスフィルタ処理を施した出力y=fg(x)を、元の信号出力f(x)から引き去った処理のことであり、ハイパスフィルタ処理の出力y=f(x)−fg(x)を与える。また、この際のローパスフィルタ処理としては、モフォロジのオープニング処理が実行される。オープニング処理とは、図9に示すような球形の構造要素によって、与えられた信号f(x)を下からトレースする処理に相当する。この際、信号が緩やかに変化しているところでは信号出力y=f(x)は正確になぞられるが、正のパルスがあるところでは球(構造要素)が入りきれないためにトレース可能な範囲から外れる。そのため、オープニング処理により、正のパルスが除かれた出力y=fg(x)が得られる。この出力y=fg(x)を元の信号出力y=f(x)から差し引くことにより、結果として、正のパルスの出力である凸形状領域における高周波成分を取得することができる。
なお、この他にも、フーリエ空間によるハイパスフィルタ処理や、DoG(Difference of Gaussian)等を用いて凸形状高周波成分の強度を算出しても良い。
続くステップS302において、孤立点除外部122は、ステップS301において算出された凸形状高周波成分の強度に基づいて孤立点を抽出し、この孤立点を正常な柔毛として除外する。
詳細には、孤立点除外部122はまず、凸形状高周波成分の強度に対して閾値処理を施すことにより、二値化画像を作成する。この際、孤立点除外部122は、凸形状高周波成分の強度が所定の閾値以下の画素には0を与え、該閾値以上の画素には1を与える。なお、閾値としては、例えば、注目領域内の凸形状高周波成分の強度から算出された1つの代表値(最大値又は平均値)に対する所定の割合の値を設定する。或いは、閾値として、0より大きい固定値を設定しても良い。後者の場合、凸形状高周波成分に対する代表値の算出及び該代表値からの閾値の算出処理を省略することができる。なお、この処理により作成された二値化画像は、凸形状高周波領域が抽出された画像となる。
続いて、孤立点除外部122は、粒子解析を実行して、二値化画像から孤立点を除外する。本実施の形態1においては、粒子解析として、二値化画像内の対象領域(粒子)の面積に基づいて解析を行う。具体的には、孤立点除外部122は、二値化画像に対して画素値1を有する領域をラベリングすることにより、ラベル番号1〜s(sはラベル数)を割り与えたラベル画像を作成する。そして、ラベル画像を走査してラベル番号ごとに画素数をカウントし、面積リストを作成する。ここで、面積リストは、ラベル番号と各ラベル番号に対応する面積とによって構成され、面積リストのサイズはラベル数sに対応する。さらに、孤立点除外部122は、各ラベル番号に対応するラベル領域の面積を所定の閾値で閾値処理し、面積が閾値以下であるラベル番号を有するラベル領域を孤立点と判定する。孤立点と判定されたラベル領域のラベル番号は0に設定される。
なお、孤立点の除外処理を行う際には、対象領域の面積以外にも、対象領域の周囲長、フェレ径(最大フェレ径)、対象領域の輪郭線上の任意の2点間の距離の最大値を示す絶対最大長、対象領域と等しい面積をもつ円の直径を示す円相当径等を用いた粒子解析を行っても良い。
ステップS303において、細線化部123は、凸形状高周波領域を細線化する。より詳細には、細線化部123は、孤立点除外部122により孤立点を除外されたラベル画像において、1以上のラベル番号を有する領域に対して細線化処理を実行する。そのために、細線化部123は、まず、ラベル画像に対し、ラベル番号が0の領域に画素値0を与え、ラベル番号が1以上の領域に画素値1を与えることにより、二値化画像を作成する。続いて、細線化部123は、二値化画像に対し、細線化フィルタ処理を施す。細線化フィルタ処理とは、図10に示すような3×3画素からなる局所領域のパターンM1〜M8を二値化画像内から順次探索して、中心画素を削除する処理のことである(参考:東京大学出版会、「画像解析ハンドブック」、第577〜578頁(細線化と縮退化))。具体的には、まず、パターンM1を二値化画像内から探索して、中心画素の画素値を0に変更する。なお、図10に示す*は考慮しなくて良い画素を示す。次に、パターンM2を二値化画像内から探索して、中心画素の画素値を0に変更する。このような処理を、パターンM8まで繰り返して1サイクルとする。このサイクルを、削除する点がなくなるまで反復する。
さらに、ステップS304において、線形状凸領域抽出部120は、細線化処理が施された二値化画像(以下、線形状画像ともいう)を出力する。その後、処理はメインルーチンに戻る。
ステップS30に続くステップS40において、領域内曲率特徴量算出部130は、線形状領域に沿った円弧の曲率に基づく曲率特徴量を算出する。図11は、領域内曲率特徴量算出部130が実行する処理を詳細に示すフローチャートである。
ステップS401において、サイズ特徴量算出部131は、円弧のサイズの特徴量として、線形状画像において観察される線形状領域に沿った円弧の曲率、及び撮像位置(カプセル型内視鏡)から円弧までの撮像距離に対応する距離情報を算出する。このサイズ特徴量算出部131が実行する処理について、図12を参照しながら詳細に説明する。
ステップS4011において、サイズ特徴量算出部131は、線形状画像内の連続する線形状領域をラベリングする。本実施の形態1においては、8連結の連結成分に関するラベル付けを行う。具体的には、サイズ特徴量算出部131は、まず、線形状画像に対してラスタスキャンを行い、画素値1を有する画素のうち、ラベル番号が付与されていない画素を探索する。そして、ラベル番号が付与されていない画素を発見した際に、当該画素を注目画素とする。
サイズ特徴量算出部131は、注目画素の上隣又は左上の画素がラベル番号を有するとき、当該上隣又は左上の画素のいずれかのラベル番号を注目画素に付与する。その後、注目画素の左隣の画素のラベル番号が注目画素のラベル番号と異なる場合には、予め用意されたルックアップテーブルに、注目画素及びその左隣の画素のラベル番号が同一連結成分に属することを記録する。
また、注目画素の上隣及び左上の画素のラベル番号が0(ラベル番号なし)で、左隣の画素がラベル番号を有するとき、左隣の画素のラベル番号を注目画素に付与する。
さらに、サイズ特徴量算出部131は、注目画素の上隣、左隣、及び左上の画素のいずれにもラベル番号が付与されていない場合、当該注目画素に新たなラベル番号を付与する。
サイズ特徴量算出部131は、これらの処理を、線形状画像内の全画素に対して実行する。そして、最後に再度ラスタスキャンを行い、ルックアップテーブルを参照しながら、同一の連結成分に属する画素群に対し、それらの画素群に付与されたラベル番号のうち、最も小さいラベル番号を選んで付け直す(参考:CG-ARTS協会、「ディジタル画像処理」、第181〜182頁)。
ステップS4012において、曲率算出部131aは、図13に示すように、線形状画像G1内の線形状領域Lの端点Pt、及び/又は線形状領域L同士の交点Pcによって区切られる線形状領域Lの各区間(即ち、端点間、交点間、又は端点と交点との間)から、1つ以上の円弧の曲率を算出する。
より詳細には、曲率算出部131aは、線形状画像を走査して、線形状領域の端点Pt又は線形状領域同士の交点Pcの画素を検出する。そして、図14に示すように、検出された端点Pt又は交点Pcの画素を開始点画素M(xs,ys)に設定する。さらに、開始点画素M(xs,ys)から線形状領域L上(即ち、同じラベル番号が付与された画素上)を所定の画素数Δsだけ辿って到着した画素を、終点画素N(xe,ye)に設定する。
続いて、曲率算出部131aは、開始点画素Mの接線の傾き、及び終点画素Nの接線の傾きを算出する。本実施の形態においては、各接線の傾きを、x軸とのなす角度DEG
s、DEG
eで表す。各接線とx軸とのなす角度DEG
i(i=s、e)は、線形状領域に沿って開始点画素M(又は終点画素N)から前後に所定の画素数だけ離れた2つの画素(座標(x
a,y
a)、(x
b,y
b))を設定し、次式(2)を用いて算出される。
なお、線形状画像の左上を座標(0,0)に設定し、左上から各行を右方向に走査する場合、座標(xa,ya)は開始点画素M(又は終点画素N)の前方座標とし、座標(xb,yb)は開始点画素M(又は終点画素N)の後方座標とする。ただし、開始点画素M(又は終点画素N)から所定の画素数だけ離れると線形状領域から外れてしまう場合(端点又は交点を越えてしまう場合)には、線形状領域の端点又は交点の座標を上記座標(xa,ya)、(xb,yb)のいずれかに設定する。
続いて、曲率算出部131aは、開始点画素Mと終点画素Nとの間の円弧の中心角αを算出する。この中心角αは、角度DEGsと角度DEGe角度との差DEGe−DEGsに相当する。
さらに、曲率算出部131aは、開始点画素Mと終点画素Nとの間の円弧の曲率κを算出する。ここで、曲率κは円弧の曲率半径Rの逆数であり、円弧の長さに対応する画素数ΔsはR×α×(π/180)と近似することができるから、曲率κは、次式(3)によって与えられる。
曲率算出部131aは、図15に示すように、線形状領域上の開始点画素M(M1、M2、…)及び終点画素N(N1、N2、…)をずらしつつ、曲率κ(κ1、κ2、…)の算出処理を繰り返す。そして、曲率算出部131aは、終点画素Nが線形状領域Lの端点Pt、又は線形状領域L同士の交点Pcを越えると曲率αの算出処理を終了し、曲率αを算出した線形状領域に対して、1から順にインクリメントされたラベル番号を付与する。
曲率算出部131aは、このような曲率κの算出及びラベル番号の付与を、図13に示す線形状画像G1から検出された全て線形状領域Lに対して実行する。図13に示す括弧付きの数字は、各線形状領域Lに付与されたラベル番号を示している。
続くステップS4013において、曲率代表値算出部131bは、曲率算出部131aによって同一のラベル番号を付与された線形状領域ごとに、円弧の曲率の代表値を算出する。なお、本実施の形態1においては、代表値として中央値を算出するが、曲率の平均値、最大値、最小値等を代表値として算出しても良い。
続くステップS4014において、距離情報算出部131cは、撮像位置から線形状領域までの奥行き方向における距離情報を算出する。具体的には、距離情報算出部131cは、元の管腔内画像内の注目領域に対応する部分から、各画素の画素値のうちでヘモグロビンによる吸光が少なく、粘膜の表層の形状を最もよく表すR成分値を取得してR成分画像を作成する。
なお、本ステップS4014においては、粘膜上の柔毛の微細構造に及ぼす影響が少ない画像であれば、他の手法により作成された画像から距離情報を取得しても良い。例えば、凸形状高周波成分算出部121の処理において説明した濃淡モフォロジによるオープニング処理を施した画像を用いても良い。
続いて、距離情報算出部131cは、曲率代表値算出部131bによって曲率が算出された各円弧の位置における距離情報を算出する。具体的には、円弧上の画素位置におけるR成分値を取得し、円弧ごとに、R成分値の平均値、中央値、最小値、最大値といった代表値を取得する。本実施の形態1においては、このR成分値の代表値を、撮像距離に対応する距離情報として用いる。なお、この場合、距離情報の値が大きいほど撮像距離が短く、距離情報の値が小さいほど撮像距離が長くなる。また、円弧の曲率中心の画素位置や、円弧の曲率中心と曲率半径とからなる円形状領域や扇状領域の内部領域におけるR成分値を、当該円弧の距離情報として用いても良い。
さらに、ステップS4015において、サイズ特徴量算出部131は、円弧の曲率及び距離情報を、円弧のサイズの特徴量として出力する。その後、処理はメインルーチンに戻る。
ステップS401に続くステップS402において、頻度分布作成部132は、サイズ特徴量算出部131から出力された特徴量の頻度分布を作成する。より詳細には、円弧の曲率及び距離情報の2軸からなる頻度分布を作成する。この際、頻度分布作成部132は、作成した頻度分布から分散等の統計量も併せて算出する。
さらに、ステップS403において、領域内曲率特徴量算出部130は、円弧の曲率及び距離情報の2軸からなる頻度分布を曲率特徴量として出力する。その後、処理はメインルーチンに戻る。
ステップS40に続くS50において、異常判別部140は、曲率特徴量に基づき注目領域に異常部が存在するか否かを判別する。
ここで、管腔内画像に見られる円弧形状の領域としては、粘膜表面の柔毛が膨張した異常部、管腔内の液体が泡状になった泡領域、粘膜の溝といった構造物の輪郭(以下、粘膜輪郭という)等が挙げられる。以下において、膨張した柔毛のことを異常柔毛という。
図16Aは、正常な柔毛モデルを示す模式図である。正常な柔毛は、通常、粘膜表面の法線方向から見た平面形状が円形で、全体として中心部が突出した形状をなすため、管腔内画像においては、柔毛の中心部及びその近傍における輝度が最も高くなる。このため、正常な柔毛が映った管腔内画像から、周囲の画素よりも輝度値が高い凸領域を抽出すると、図16Bに示すように、柔毛の中心部に対応する点状の凸領域PR1が得られる。このような点状の凸領域PR1は、孤立点として予め線形状画像から除外されている。
図17Aは、異常柔毛モデルを示す模式図である。異常柔毛を粘膜表面の法線方向から見た平面形状は、楕円のような扁平状をなすため、管腔内画像においては、撮像方向に対する柔毛の向きに応じた輪郭部分(曲率が高い輪郭部分、又は曲率が低い輪郭部分)の輝度が高くなる。従って、異常柔毛が映った管腔内画像から凸領域を抽出すると、図17Bに示すように、撮像方向に応じた曲率の円弧状をなす線形状領域PR2が得られる。このため、線形状領域PR2に沿った円弧の曲率の分布は、図17Cに示すように、ある程度の狭い範囲に収まる。また、各線形状領域PR2に沿った円弧の曲率代表値として曲率の最小値を採用する場合、曲率は値が小さい側に偏って分布する(実線参照)。一方、曲率代表値として曲率の最大値を採用する場合には、曲率は値が大きい側に偏って分布する(破線参照)。
図18Aは、泡モデルを示す模式図である。泡領域において、1つひとつの泡はほぼ球形状をなす。また、泡の表面全体は高反射領域となっているため、管腔内画像においては、泡の輪郭部分の輝度が高くなる。従って、泡領域が映った管腔内画像から凸領域を抽出すると、図18Bに示すように、ほぼ円形状をなす様々な大きさの線形状領域PR3が得られる。このため、線形状領域PR3に沿った円弧の曲率の分布は、図18Cに示すように、大きくばらつく。
また、粘膜輪郭に沿った円弧の曲率は、一般に、異常柔毛や泡領域と比較して小さい値を取る。
図19は、円弧の曲率及び距離情報の2軸からなる頻度分布を示す模式図である。粘膜輪郭を示す線形状領域に沿った円弧の曲率は、概ね、撮像距離にはよらず、値が小さく、且つ狭い範囲に分布する。また、泡領域を示す線形状領域に沿った円弧の曲率は、撮像距離によらず、広い範囲に分布する。それに対して、異常柔毛を示す線形状領域に沿った円弧の曲率の分布範囲は、撮像距離に応じて変化する。具体的には、撮像距離が長くなると(撮像情報の値が小さくなると)、曲率は値が大きい方(円弧のカーブが強い方)に分布し、撮像距離が短くなると(撮像情報の値が大きくなると)、曲率は値が小さい方(円弧のカーブが緩い方)に分布する。
そこで、異常判別部140は、以下のように判別を行う。
まず、異常判別部140は、ある距離情報Larbにおける曲率の分散が、距離情報Larbに応じて予め定められた所定の閾値(分散閾値)TH1よりも大きい場合、当該注目領域には泡領域が映っており、異常部は存在しないと判別する。
また、異常判別部140は、距離情報Larbにおける曲率の分散が上記分散閾値TH1以下である場合において、曲率が所定の閾値(曲率閾値)TH2よりも小さい範囲に分布する場合、当該注目領域には粘膜輪郭が映っており、異常部は存在しないと判別する。
一方、異常判別部140は、距離情報Larbにおける曲率の分散が上記分散閾値TH1以下である場合において、曲率が上記曲率閾値TH2以上の範囲に分布する場合、当該注目領域には異常柔毛が映っており、異常部が存在すると判別する。
続くステップS60において、演算部100は、ステップS50においてなされた異常部の有無の判別結果を出力する。これに応じて、制御部10は、判別結果を表示部40に表示させると共に、処理対象の管腔内画像の画像データと関連付けて記録部50に記録する。
以上説明したように、実施の形態1によれば、管腔内画像から抽出された粘膜領域を注目領域として設定し、この注目領域から、周囲の画素よりも画素値が高い画素群からなる凸領域のうちの高周波成分を抽出する。そして、該高周波成分を細線化した線形状領域に対し、曲率の特徴量として曲率及び距離情報を算出し、これらの特徴量の頻度分布に基づいて注目領域に映った被写体を識別する。従って、粘膜輪郭や泡等の誤検出(過検出)を抑制し、柔毛の微細構造の撮像距離に応じたサイズの変化や、撮像方向に対する柔毛の向きに応じた形状の変化にも対応して、異常柔毛を精度良く識別して、異常部の有無を判別することが可能となる。
(変形例1−1)
上記実施の形態1においては、管腔内画像から不要領域を除去した残りの領域全体を1つの注目領域として設定したが、残りの領域を分割した複数の領域を注目領域としてそれぞれ設定しても良い。この場合、管腔内画像内のどの位置に異常部が存在するかを特定することが可能となる。分割方法は特に限定されず、例えば、不要領域を除去した残りの領域を単にマトリクス状に分割しても良い。或いは、不用領域を除去した残りの領域内の各画素の画素値から撮像距離に対応する距離情報を取得し、距離情報を階層分けした階層ごとに領域を分割しても良い。この場合、図12に示すステップS4014の処理を省略することができる。
(変形例1−2)
次に、実施の形態1の変形例1−2について説明する。
図20は、変形例1−1に係る画像処理装置を示すブロック図である。図20に示すように、変形例1−1に係る画像処理装置1−2は、図1に示す演算部100の代わりに、演算部100−2を備える。
演算部100−2は、図1に示す線形状凸領域抽出部120の代わりに、線形状凸領域抽出部150を備える。この線形状凸領域抽出部150は、画素値が尾根形状に変化する画素群を尾根形状の領域として抽出する尾根形状抽出部151と、抽出された尾根形状の領域から、周囲に対して孤立している画素群を孤立点として除外する孤立点除外部152と、該孤立点を除外した後の尾根形状の領域を細線化して線形状凸領域を抽出する細線化部153とを有する。
演算部100−2の動作は全体として図2に示すものと同様であり、ステップS30において線形状凸領域抽出部150が実行する処理の内容のみが異なる。図21は、線形状凸領域抽出部150が実行する処理の詳細を示すフローチャートである。
まず、ステップS311において、尾根形状抽出部151は、注目領域から尾根形状の領域を抽出する。より詳細には、尾根形状抽出部151は、注目領域内において注目画素を順次設定し、各注目画素における画素値の勾配の最大方向(以下、最大勾配方向という)及び最小方向(以下、最小勾配方向という)を求める。そして、最大勾配方向及び最小勾配方向の各々に対し、画素値が変化する形状の曲率を算出し、これらの曲率の比をもとに尾根領域を検出する。
そのために、尾根形状抽出部151は、まず、注目画素に対し、次式(4)に示すヘッセ行列Hによる固有方程式を解くことにより、固有値である最大勾配方向における曲率k1、及び、最小勾配方向における曲率k2を算出する。
det(H−λE)=0 …(4)
式(4)において、ヘッセ行列Hは、次式(5)によって与えられる。
また、記号Eは単位行列であり、λEは次式(6)のとおりである。
続いて、曲率k1、k2から、次式(7−1)及び(7−2)によって与えられるガウス曲率KG及び平均曲率Akを算出する。
KG=k1×k2 …(7−1)
Ak=(k1+k2)/2 …(7−2)
このとき、Ak<0、且つ、KG≒0である画素群、即ち、絶対値|KG|が所定の閾値以下である画素群が、尾根形状の領域であると判別される。尾根形状抽出部151は、尾根形状の領域と判別された画素群に対して画素値1を付与し、その他の画素に画素値0を付与することにより、二値化画像を作成する。
なお、曲率k1、k2を算出した後、例えばShape index and Curvednessのように、任意の3次元表面形状を認識する公知技術を用いて尾根形状の領域の判別を行っても良い(参考:Chitra Dorai,Anil K Jain,“COSMOS-A Representation Scheme for Free-Form Surfaces”)。
続くステップS312において、孤立点除去部152は、ステップS311において抽出された尾根形状の領域から孤立点を抽出し、該孤立点を正常な柔毛として除外する。より詳細には、ステップS311において作成された二値化画像に対し、連続する画素の数が所定値以下である画素群(即ち、所定値以下のサイズの領域)を孤立点として抽出する。そして、孤立点として抽出した画素群に対し、画素値0を付与する。
続くステップS313において、細線化部153は、孤立点を除外した後の尾根形状の領域を細線化する。より詳細には、細線化部153は、ステップS312において作成された二値化画像に対し、図10に示すパターンM1〜M8を用いた細線化フィルタ処理を施す。なお、細線化フィルタ処理の詳細な内容については、実施の形態1と同様である。
さらに、ステップS314において、線形状凸領域抽出部150は、細線化処理が施された二値化画像を出力する。その後、処理はメインルーチンに戻る。
(変形例1−3)
実施の形態1において、注目領域から抽出された凸形状高周波領域を細線化する処理(図8のステップS303参照)としては、上述した方法以外にも、種々の公知技術を適用することができる。以下、公知技術の1つであるHilditch細線化アルゴリズムについて説明する。
まず、第1ステップとして、凸形状高周波領域が抽出された二値化画像に対し、細線化対象である画素Pkのうち、以下に示す6つの条件を満たす画素を境界画素として逐次削除する。ここで、kは二値化画像内の画素の画素番号である(kは自然数)。また、画素Pkの画素値をB(Pk)(B(Pk)=1)と示す。この場合、第1ステップは、画素値B(Pk)を1から−1に置き換える処理に相当する。なお、細線化の非対象画素の画素値はB(Pk)=0である。
条件1:注目画素が細線化の対象画素である。即ち、次式(a1)を満たす。
B(Pk)=1 …(a1)
条件2:注目画素に対し、縦、横で隣接する画素のうちいずれか1つが0である。即ち、次式(a2)を満たす。
条件4:孤立点でない。即ち、次式(a4)を満たす。ここでは、逐次処理により得られる結果をB(P
k)、1つ前のラスタ操作完了時の結果をB
*(P
k)として区別して示している。なお、画素値B
*(P
k)は−1を取らない。
条件5:連結性を保持する。即ち、図22に示すパターンM11〜M14のいずれかに当てはまる。なお、図22において、※で示す画素は、画素値1又は−1を取る。この条件を条件式で表すと、次式(a5)となる。
条件6:線幅が2である線分に対して、その片方のみを削除する。即ち、1つ前のラスタ操作完了時の結果B
*(P
k)が、図23に示すパターンM21〜M24のいずれかに当てはまる。この条件を条件式で表すと、次式(a6)となる。
次に第2ステップとして、境界画素として逐次削除された画素(即ち、画素値B(Pk)が1から−1に置き換えられた画素)の画素値を、非対象領域の画素値B(Pk)=0に置き換える。
これらの第1ステップ及び第2ステップを、非対象領域の画素値への置換えがなくなるまで繰り返す。それにより、対象領域の細線化がなされる。
(実施の形態2)
次に、本発明の実施の形態2について説明する。
図24は、実施の形態2に係る画像処理装置の構成を示すブロック図である。図24に示すように、実施の形態2に係る画像処理装置2は、図1に示す演算部100の代わりに演算部200を備える。この演算部200は、注目領域設定部110と、線形状凸領域抽出部120と、領域内曲率特徴量算出部210と、異常判別部220とを備える。なお、注目領域設定部110及び線形状凸領域抽出部120の構成及び動作については、実施の形態1と同様である。また、演算部200以外の画像処理装置の各部の構成及び動作についても、実施の形態1と同様である。
領域内曲率特徴量算出部210は、線形状凸領域抽出部120によって抽出された線形状領域に沿った1つ以上の円弧の曲率を算出すると共に、それらの曲率のばらつきを特徴量として算出する形状特徴量算出部211と、特徴量の頻度分布を作成する頻度分布作成部212とを備える。より詳細には、形状特徴量算出部211は、線形状領域の端点及び/又は線形状領域同士の交点によって区切られる線形状領域の各区間から1つ以上の円弧の曲率を算出する曲率算出部211aと、各区間から算出された1つ以上の円弧の曲率の標準偏差を算出する曲率標準偏差算出部211bとを有する。
異常判別部220は、上記曲率のばらつきに基づいて異常部の有無を判別する。
次に、画像処理装置2の動作について説明する。
画像処理装置2の動作は全体として図2に示すものと同様であり、ステップS40において領域内曲率特徴量算出部210が実行する処理、及びステップS50において異常判別部220が実行する処理の内容が実施の形態1とは異なる。
図25は、領域内曲率特徴量算出部210が実行する処理を示すフローチャートである。まず、ステップS421において、形状特徴量算出部211は、ステップS30において抽出された線形状領域の形状を示す特徴量として、線形状領域に沿った円弧の曲率のばらつきを算出する。この形状特徴量算出部211が実行する処理について、図26を参照しながら詳細に説明する。なお、図26に示すステップS4211及びS4212は、図12に示すステップS4011及びS4012と対応している。このうち、ステップS4012の処理は、曲率算出部131a(図1参照)の代わりに、曲率算出部211aが実行する。
ステップS4212に続くステップS4213において、曲率標準偏差算出部211bは、線形状領域の区間ごとに、1つ以上の円弧の曲率の標準偏差を算出する。具体的には、曲率標準偏差算出部211bは、曲率算出部211aによって付与された同一のラベル番号を有する円弧に対して算出されたm個の曲率κ
i(i=1〜m)から、次式(8)によって与えられる標準偏差σを算出する。
式(8)において、κの上にバーが附されたκバーは、曲率κ
iの平均値であり、次式(9)によって与えられる。
続くステップS4214において、形状特徴量算出部211は、線形状領域の各区間について算出された標準偏差を、線形状領域の形状を表す特徴量として出力する。その後、処理はメインルーチンに戻る。
ステップS421に続くステップS422において、頻度分布作成部212は、特徴量の頻度分布として、形状特徴量算出部211から出力された標準偏差の頻度分布を作成する。
さらに、ステップS423において、領域内曲率特徴量算出部210は、標準偏差の頻度分布を曲率特徴量として出力する。その後、処理はメインルーチンに戻る。
次に、異常判別部220が実行する処理について説明する。
異常判別部220は、曲率特徴量である標準偏差の頻度分布に基づき注目領域に異常部が存在するか否かを判別する。
ここで、注目領域に異常柔毛が存在する場合、図27Aに示すように、当該注目領域から作成された線形状画像には、扁平な線形状領域PR4が現れる。図27Bに示すように、このような線形状領域PR4に沿った円弧の曲率は、線形状領域PR4内の位置に応じて様々な値を取る。このため、曲率のばらつきは比較的大きくなる。なお、図27Bに示すベクトルvi(i=1、2、…)は、円弧の曲率の大きさを表すベクトルである。従って、図27Cに示すように、線形状領域PR4について算出された曲率の標準偏差は、比較的、値が大きい範囲に分布する。
一方、注目領域に泡領域が存在する場合、図28Aに示すように、当該注目領域から作成された線形状画像には、円形の線形状領域PR5が現れる。図28Bに示すように、このような線形状領域PR5において、円弧の曲率は、ほぼ一定値となり、曲率のばらつきは小さくなる。従って、図28Cに示すように、線形状領域PR5について算出された曲率の標準偏差は、値が小さい範囲(0又はその近傍)に分布する。
そこで、異常判別部220は、標準偏差の頻度分布に基づいて、異常部の有無を判別する。具体的には、標準偏差の頻度分布が、所定の閾値よりも大きい範囲に偏っている場合、当該注目領域に異常部が存在すると判別する。
以上説明したように、実施の形態2によれば、管腔内画像から抽出された粘膜領域を注目領域として設定し、この注目領域から、周囲の画素よりも画素値が高い画素群からなる凸領域のうちの高周波成分を抽出する。そして、該高周波成分を細線化した線形状領域に対し、線形状領域に沿った円弧の曲率の標準偏差の頻度分布を算出し、この標準偏差の頻度分布に基づいて、注目領域における異常部の有無を判別する。従って、異常柔毛と泡領域とを精度良く識別して、異常部の判別精度を向上させることが可能となる。
また、実施の形態2によれば、注目領域内の被写体までの撮像距離(距離情報)を用いることなく判別を行うので、演算処理を簡略化することが可能となる。
(実施の形態3)
次に、本発明の実施の形態3について説明する。
図29は、実施の形態3に係る画像処理装置の構成を示すブロック図である。図29に示すように、実施の形態3に係る画像処理装置3は、図1に示す演算部100の代わりに演算部300を備える。この演算部300は、注目領域設定部110と、線形状凸領域抽出部120と、領域内曲率特徴量算出部310と、異常判別部320とを備える。なお、注目領域設定部110及び線形状凸領域抽出部120の構成及び動作については、実施の形態1と同様である。また、演算部300以外の画像処理装置3の各部の構成及び動作についても、実施の形態1と同様である。
領域内曲率特徴量算出部310は、線形状凸領域抽出部120によって抽出された線形状領域に沿った1つ以上の円弧の中心方向を示す曲率特徴量を算出する方向特徴量算出部311と、曲率特徴量の頻度分布を作成する頻度分布作成部312とを備える。このうち、方向特徴量算出部311は、各円弧から曲率中心に向かう方向(以下、曲率中心方向という)を算出する曲率中心方向算出部311aと、各円弧の位置における被写体の勾配方向を算出する勾配方向算出部311bとを有する。ここで、被写体の勾配方向とは、画像の奥行き方向に被写体(具体的には粘膜構造)が勾配する方向のことである。
異常判別部320は、粘膜構造の勾配の方向と、曲率中心方向の分布とに基づいて、異常部の有無を判別する。
次に、画像処理装置3の動作について説明する。
画像処理装置3の動作は全体として図2に示すものと同様であり、ステップS40において領域内曲率特徴量算出部310が実行する処理、及びS50において異常判別部220が実行する処理の内容が実施の形態1とは異なる。
図30は、領域内曲率特徴量算出部310が実行する処理を示すフローチャートである。まず、ステップS431において、方向特徴量算出部311は、円弧の中心方向を示す特徴量として、円弧の曲率中心方向及び勾配方向を算出する。この方向特徴量算出部311が実行する処理について、図31を参照しながら詳細に説明する。
ステップS4311において、曲率中心方向算出部311aは、線形状領域に沿った円弧の曲率中心方向を算出する。具体的には、まず、図12に示すステップS4012と同様にして、各円弧の曲率を算出し、さらに、曲率の逆数である曲率半径Rを算出する。
続いて、曲率中心方向算出部311aは、当該円弧上の位置座標(x,y)及び曲率半径Rを用いてハフ(hough)変換を実行することにより近似円を作成する(参考:CG−ARTS協会、「ディジタル画像処理」、第213〜214頁)。ここで、ハフ変換による近似円の検出は、円弧上の位置座標(x,y)を通る半径Rの円を、円の中心座標(a,b)及び半径Rからなるパラメータ空間に投票し、その投票結果を調べることにより作成することができる。
さらに、曲率中心方向算出部311aは、円弧上の任意の点(例えば、中心点)から近似円の中心に向かう方向ベクトルを算出し、この方向ベクトルを曲率中心方向として出力する。
続くステップS4312において、勾配方向算出部311bは、線形状領域に沿った円弧の位置における被写体の勾配方向を算出する。具体的には、勾配方向算出部311bは、元の管腔内画像のうち処理対象である注目領域に対応する部分から、各画素の画素値のうちでヘモグロビンによる吸光が少なく、粘膜の表層の形状に最も近いR成分値を取得してR成分画像を作成する。
なお、本ステップS4312においては、粘膜上の柔毛の微細構造に及ぼす影響が少ない画像であれば、他の手法により作成された画像から勾配方向を取得しても良い。例えば、凸形状高周波成分算出部121の処理において説明した濃淡モフォロジによるオープニング処理により画像を用いても良い。
続いて、勾配方向算出部311bは、R成分画像から、曲率中心方向算出部311aにより曲率中心方向を算出した円弧の位置における勾配方向を算出する。具体的には、エッジ強度を算出するためのフィルタ処理(参考:CG-ARTS協会、「ディジタル画像処理、第114〜117頁)である1次微分フィルタ処理(プリューウイットフィルタ、ソーベルフィルタ等)を実行する。
ステップS4313において、方向特徴量算出部311は、円弧の曲率中心方向及び勾配方向を特徴量として出力する。その後、処理はメインルーチンに戻る。
ステップS431に続くステップS432において、頻度分布作成部312は、方向特徴量算出部311から出力された特徴量の頻度分布を作成する。具体的には、円弧の勾配方向別に、円弧の曲率中心方向の頻度分布を作成する。
ステップS433において、領域内曲率特徴量算出部310は、円弧の曲率中心方向の頻度分布を曲率特徴量として出力する。
次に、異常判別部320が実行する処理を説明する。
異常判別部320は、曲率特徴量である曲率中心方向の頻度分布に基づき、注目領域に異常部が存在するか否かを判別する。
ここで、注目領域に異常柔毛が存在する場合、図32Aに示すように、被写体(粘膜構造)の勾配に対して1つの撮像方向から撮像を行うと、図32Bに示すように、線形状領域PR6がほぼ1方向を向いた線形状画像が得られる。この場合、図32Cに示すように、線形状領域PR6における曲率中心方向dは概ね一定の方向を向く。従って、図32Dに示すように、勾配方向別の曲率中心方向の分布には偏りが生じる。
一方、注目領域に泡領域が存在する場合、図33Aに示すように、被写体(粘膜構造)の勾配に対して1つの撮像方向から撮像を行うと、図33Bに示すように、線形状領域PR7が泡領域の輪郭に対応する円形状をなす線形状画像が得られる。この場合、図33Cに示すように、線形状領域PR7における曲率中心方向dはあらゆる方向を向く。従って、図33Dに示すように、勾配方向別の曲率中心方向は、幅広く、且つ概ね均一に分布する。
そこで、異常判別部320は、勾配方向別の曲率中心方向の頻度分布において、曲率中心方向が偏っている場合、即ち、曲率中心方向の分散が所定の閾値以下である場合、当該注目領域に異常部が存在すると判別する。
以上説明したように、実施の形態3によれば、管腔内画像から抽出された粘膜領域を注目領域として設定し、この注目領域から、周囲の画素よりも画素値が高い画素群からなる凸領域のうちの高周波成分を抽出する。そして、該高周波成分を細線化した線形状領域に対し、線形状領域に沿った円弧の曲率中心方向及び勾配方向を算出し、勾配方向別の曲率中心方向の頻度分布に基づいて、注目領域における異常部の有無を判別する。従って、撮像方向によって変化する柔毛の向きに応じて、異常柔毛と泡領域とを精度良く識別し、異常部の判別精度を向上させることが可能となる。
(実施の形態4)
次に、本発明の実施の形態4について説明する。
図34は、実施の形態4に係る画像処理装置の構成を示すブロック図である。図34に示すように、実施の形態4に係る画像処理装置4は、図1に示す演算部100の代わりに、演算部400を備える。この演算部400は、注目領域設定部410と、線形状凸領域抽出部120と、領域内曲率特徴量算出部420と、異常判別部430とを備える。なお、線形状凸領域抽出部120の構成及び動作については、実施の形態1と同様である。また、演算部400以外の画像処理装置4の各部の構成及び動作についても、実施の形態1と同様である。
注目領域設定部410は、処理対象の管腔内画像から残渣や暗部といった粘膜以外の領域を除外して、粘膜領域を抽出する粘膜領域抽出部111と、抽出された粘膜領域をさらに複数の領域に分割する領域分割部411とを備え、分割された各領域を注目領域に設定する。なお、粘膜領域抽出部111の動作は実施の形態1と同様である。
領域内曲率特徴量算出部420は、曲率算出部131a及び曲率代表値算出部131bを備えるサイズ特徴量算出部421と、頻度分布作成部132とを備える。なお、曲率算出部131a、曲率代表値算出部131b、及び頻度分布作成部132の動作は実施の形態1と同様である。
異常判別部430は、領域内曲率特徴量算出部420によって算出された曲率特徴量に基づいて、各注目領域に異常部が存在するか否かを判別する。
次に、画像処理装置4の動作について説明する。
画像処理装置4の動作は全体として図2に示すものと同様であり、ステップS20〜S60の各々において実行される処理の内容が実施の形態1とは異なる。
図35は、ステップS20において注目領域設定部410が実行する処理を示すフローチャートである。なお、図35におけるステップS201〜S205は、図3と対応している。
ステップS205に続くステップS241において、注目領域設定部410は、不要領域を除去した後に残った粘膜領域を、各々が所定サイズ以下の複数の領域に分割する。分割方法は特に限定されず、本実施の形態4においては、粘膜領域を矩形状の領域に分割する。また、分割される1つの領域のサイズは、分割された領域内で被写体までの撮像距離の差が所定の範囲に収まるサイズとなるように、予め設定されている。
続くステップS242において、注目領域設定部410は、分割された各領域の座標情報を注目領域として出力する。
ステップS30において、線形状凸領域抽出部120は、ステップS20において出力された注目領域の各々から線形状領域を抽出する。
ステップS40において、領域内曲率特徴量算出部420は、ステップS20において出力された注目領域ごとに、線形状領域に沿った円弧の曲率に基づく曲率特徴量を算出する。なお、曲率特徴量としては、実施の形態1と同様に、線形状領域の各区間から算出された1つ以上の円弧の曲率の代表値の頻度分布が算出されるが、実施の形態1とは異なり、撮像情報は算出されない。
ステップS50において、異常判別部430は、注目領域ごとに、曲率特徴量である円弧の曲率の代表値の頻度分布に基づいて、当該注目領域に異常部が存在するか否かを判別する。具体的には、まず、曲率が所定の閾値(曲率閾値)よりも小さい場合、当該注目領域には粘膜輪郭が映っており、異常部は存在しないと判別する。
また、異常判別部430は、曲率が上記曲率閾値より大きく、且つ曲率の分散が所定の閾値(分散閾値)よりも大きい場合、当該注目領域には泡領域が映っており、異常部は存在しないと判別する。これは、撮像距離が一定であるとしても、本来、泡領域には様々な曲率の泡が含まれるからである。
一方、異常判別部430は、曲率が上記曲率閾値より大きく、且つ曲率の分散が所定の閾値(分散閾値)以下である場合、当該注目領域に異常柔毛が映っており、異常部が存在すると判別する。これは、本来、近傍領域内における異常柔毛は互いに類似した形状を有しているため、注目領域のサイズが小さく、1つの注目領域内における撮像距離の差を無視できる場合には、異常柔毛の輪郭に対応する円弧の曲率が揃うからである。
ステップS60において、演算部400は、ステップS50において各注目領域に対してなされた判別の結果を、注目領域の座標情報と共に出力する。これに応じて、制御部10は、判別結果を表示部に表示させると共に、処理対象の管腔内画像の画像データと関連付けて記録部50に記録する。この際、制御部10は、表示部40に表示された管腔内画像に対し、異常部が存在すると判別された注目領域の位置を表すマーク等を重ねて表示させても良い。
以上説明したように、実施の形態4によれば、管腔内画像から抽出された粘膜領域を分割することにより複数の注目領域を設定し、各注目領域から、周囲の画素よりも画素値が高い画素群からなる凸領域のうちの高周波成分を抽出する。そして、該高周波成分を細線化した線形状領域に対し、曲率の特徴量として曲率の頻度分布を算出し、この特徴量に基づいて、注目領域に映った被写体を判別する。従って、異常柔毛と泡領域及び粘膜輪郭とを精度良く識別して、異常部の有無を判別することが可能となる。
また、実施の形態4によれば、管腔内画像から抽出された粘膜領域を、被写体までの撮像距離の差を無視できる程度に分割し、分割された各領域を注目領域に設定してから異常部の有無を判別するための各処理を行うので、距離情報の算出処理及び距離情報に応じた異常部の判別処理を省略することができ、演算処理を簡単にすることができる。
また、実施の形態4によれば、注目領域ごとに異常部の有無の判別を行うので、管腔内画像における異常部の存在範囲を特定することが可能となる。
(変形例4−1)
次に、実施の形態4の変形例4−1について説明する。
注目領域設定部410が設定する注目領域のサイズは可変としても良い。以下、注目領域のサイズを可変とする場合の画像処理装置4の動作について説明する。図36は、変形例4−1における画像処理装置4の動作を示すフローチャートである。
ステップS20において、注目領域設定部410は、不要領域を除去した残りの粘膜領域を1つの注目領域として設定する。
続くステップS30〜S50における処理内容は、実施の形態4と同様である。
ステップS50に続くステップS71において、演算部400は、設定された注目領域のいずれかにおいて異常部が存在するとの判別結果が得られたか否かを判定する。なお、1巡目においては、注目領域は1つである。
異常部が存在するとの判別結果が得られた場合(ステップS71:Yes)、処理はステップS60に移行する。なお、ステップS60における処理内容は、実施の形態4と同様である。
一方、異常部が存在するとの判別結果が得られていない場合(ステップS71:No)、演算部400は、現在設定されている注目領域のサイズが所定サイズ以下であるか否かを判定する(ステップS72)。なお、この際の所定サイズとしては、注目領域内の被写体に対する撮像距離の差を無視できる程度のサイズが設定される。
注目領域のサイズが上記所定サイズ以下である場合(ステップS72:Yes)、処理はステップS60に移行する。
一方、注目領域のサイズが上記所定サイズよりも大きい場合(ステップS72:No)、注目領域設定部410は、現在設定されている注目領域よりもサイズを小さくし、粘膜領域に対して注目領域を再設定する(ステップS73)。その後、処理はステップS30に移行する。
以上説明したように、本変形例4−1によれば、注目領域のサイズを徐々に小さくしながら、各注目領域における異常部の有無の判別処理を行う。ここで、注目領域のサイズを小さくすることにより、各注目領域内における被写体までの撮像距離の差が小さくなるので、ステップS30〜ステップS73の処理を一巡するごとに、各注目領域に対する異常部の判別精度が上がることになる。このように管腔内画像に対する異常部の判別精度を低い方から高い方に変化させることにより、異常部の判別処理を効率化しつつ、異常部の判別漏れを抑制することが可能となる。
以上説明した実施の形態1〜4並びにこれらの変形例に係る画像処理装置は、記録媒体に記録された画像処理プログラムをパーソナルコンピュータやワークステーション等のコンピュータシステムで実行することによって実現することができる。また、このようなコンピュータシステムを、ローカルエリアネットワーク、広域エリアネットワーク(LAN/WAN)、又は、インターネット等の公衆回線を介して、他のコンピュータシステムやサーバ等の機器に接続して使用しても良い。この場合、実施の形態1〜4に係る画像処理装置は、これらのネットワークを介して管腔内画像の画像データを取得したり、これらのネットワークを介して接続された種々の出力機器(ビュアーやプリンタ等)に画像処理結果を出力したり、これらのネットワークを介して接続された記憶装置(記録媒体及びその読取装置等)に画像処理結果を格納するようにしても良い。
なお、本発明は、実施の形態1〜4及びこれらの変形例に限定されるものではなく、各実施の形態や変形例に開示されている複数の構成要素を適宜組み合わせることによって、種々の発明を形成できる。例えば、各実施の形態や変形例に示される全構成要素からいくつかの構成要素を除外して形成しても良いし、異なる実施の形態や変形例に示した構成要素を適宜組み合わせて形成しても良い。