(57)【要約】【課題】過剰分割や過剰統合の問題を回避できる、カーネル形状の決定を可能とする。【解決手段】任意のＮ次元ベクトルの注目サンプルと、γ次元の位置及びη次元の特徴量を持つＮ次元空間内のサンプル集合中の複数の候補サンプルとを結ぶパスを設定し、このパス上で前記特徴量を積分して積分値を求め、この積分値を予め設定されたカーネル関数の変数値として注目サンプルに対応するカーネル関数値を算出する。【選択図】  図２
2008021175.tif 000002
【特許請求の範囲】
【請求項１】
任意のＮ次元ベクトル（Ｎは自然数）の第１画素と、γ次元の位置及びη次元の特徴量を持つＮ次元空間（γは２以上の自然数，ηは１以上の自然数）内の複数の第２画素とを結ぶパスを設定するステップと；
前記パス上で前記特徴量を積分して積分値を得るステップと；
前記積分値を予め設定されたカーネル関数の変数値として前記第１画素に対応するカーネル関数値を算出するステップと；
前記カーネル関数値と前記第２画素の特徴量との積和を前記カーネル関数値の総和で割ることにより平均値シフトベクトルを算出するステップと；
前記第１画素を前記平均値シフトベクトルに従って移動させるステップと；
前記第１画素の特徴量を移動後の第１画素の特徴量で置き換えるステップと；
を具備する画像のフィルタリング方法。
【請求項２】
前記積分値は、前記第１画素の特徴量と前記パス上に設定された複数の位置に存在する各画素の特徴量との差分の総和である請求項１記載の方法。
【請求項３】
前記積分値は、前記パス上に設定された複数の位置に存在する隣接画素の特徴量の差分の総和である請求項１記載の方法。
【請求項４】
前記カーネル関数は、正規分布関数である請求項１記載の方法。
【請求項５】
前記カーネル関数は、前記変数値が予め定められた値を超える場合は０を出力し、そうでない場合は予め定められた定数を出力する関数である請求項１記載の方法。
【請求項６】
任意のＮ次元ベクトル（Ｎは自然数）の第１画素と、γ次元の位置及びη次元の特徴量を持つＮ次元空間（γは２以上の自然数，ηは１以上の自然数）内の複数の第２画素とを結ぶパスを設定するステップと；
前記パス上で前記特徴量を積分して積分値を得るステップと；
前記積分値を予め設定されたカーネル関数の変数値として前記第１画素に対応するカーネル関数値を算出するステップと；
前記カーネル関数値と前記第２画素の特徴量との積和を前記カーネル関数値の総和で割ることにより平均値シフトベクトルを算出するステップと；
前記第１画素を前記平均値シフトベクトルに従って移動させるステップと；
前記第１画素の特徴量を移動後の第１画素の特徴量で置き換えるステップと；
前記移動後の第１画素のうち特徴量の類似度が予め定められた閾値以上の複数の画素に対して同一のクラスタを割り当てるステップと；を具備する画像のクラスタリング方法。
【請求項７】
前記積分値は、前記第１画素の特徴量と前記パス上に設定された複数の位置に存在する各画素の特徴量との差分の総和である請求項６記載の方法。
【請求項８】
前記積分値は、前記パス上に設定された複数の位置に存在する隣接画素の特徴量の差分の総和である請求項６記載の方法。
【請求項９】
前記カーネル関数は、正規分布関数である請求項６記載の方法。
【請求項１０】
前記カーネル関数は、前記変数値が予め定められた値を超える場合は０を出力し、そうでない場合は予め定められた定数を出力する関数である請求項６記載の方法。
【請求項１１】
２次元の位置と１次元の特徴量を持つ複数の画素からなる画像を取得するステップと；
前記画像に対して被写体領域及び背景領域をそれぞれ概略的に表す初期第１領域及び初期第２領域を設定するステップと；
任意のＮ次元ベクトル（Ｎは自然数）の第１画素と、前記画像中の複数の第２画素とを結ぶパスを設定するステップと；
前記パス上で前記特徴量を積分して積分値を得るステップと；
前記積分値を予め設定されたカーネル関数の変数値として前記第１画素に対応するカーネル関数値を算出するステップと；
前記カーネル関数値のうち前記初期第１領域に属する各画素における第１カーネル関数値を加算して第１確率密度関数値を得るステップと；
前記カーネル関数値のうち前記初期第２領域に属する各画素における第２カーネル関数値を加算して第２確率密度関数値を得るステップと；
前記第１確率密度関数値及び第２確率密度関数値を比較して前記画像中の少なくとも一部の画素に設定された注目画素が前記被写体領域及び背景領域のいずれに属するかを表す領域情報を得るステップと；
前記領域情報を出力するステップと；を具備する画像の領域フィッティング方法。
【請求項１２】
前記積分値は、前記第１画素の特徴量と前記パス上に設定された複数の位置に存在する各画素の特徴量との差分の総和である請求項１１記載の方法。
【請求項１３】
前記積分値は、前記パス上に設定された複数の位置に存在する隣接画素の特徴量の差分の総和である請求項１１記載の方法。
【請求項１４】
前記カーネル関数は、正規分布関数である請求項１１記載の方法。
【請求項１５】
前記カーネル関数は、前記変数値が予め定められた値を超える場合は０を出力し、そうでない場合は予め定められた定数を出力する関数である請求項１１記載の方法。
【請求項１６】
任意のＮ次元ベクトル（Ｎは自然数）の第１画素と、γ次元の位置及びη次元の特徴量を持つＮ次元空間（γは２以上の自然数，ηは１以上の自然数）内の複数の第２画素とを結ぶパスを設定する設定部と；
前記パス上で前記特徴量を積分して積分値を得る積分部と；
前記積分値を予め設定されたカーネル関数の変数値として前記第１画素に対応するカーネル関数値を算出する第１算出部と；
前記第２サンプルの特徴量を抽出する抽出部と；
前記特徴量を記憶する記憶部と；
前記カーネル関数値と前記第２サンプルの特徴量との積和を前記カーネル関数値の総和で割ることにより平均値シフトベクトルを算出する第２算出部と；
前記第１サンプルを前記平均値シフトベクトルに従って移動させ、前記第１サンプルの特徴量を移動後の第１サンプルの特徴量で置き換える置き換え部と；
を具備する画像のための平均値シフトフィルタ。
【請求項１７】
２次元の位置と１次元の特徴量を持つ複数の画素からなる画像を取得する取得部と；
前記画像に対して被写体領域及び背景領域をそれぞれ概略的に表す初期第１領域及び初期第２領域を設定する第１設定部と；
任意のＮ次元ベクトル（Ｎは自然数）の第１画素と、前記画像中の複数の第２画素とを結ぶパスを設定する第２設定部と；
前記パス上で前記特徴量を積分して積分値を得る積分部と；
前記積分値を予め設定されたカーネル関数の変数値として前記第１画素に対応するカーネル関数値を算出する第１算出部と；
前記カーネル関数値のうち前記画像の被写体領域を概略的に表す初期第１領域に属する各画素における第１カーネル関数値を加算して第１確率密度関数値を得る第１加算部と；
前記カーネル関数値のうち前記画像の背景領域を概略的に表す初期第２領域に属する各画素における第２カーネル関数値を加算して第２確率密度関数値を得る第２加算部と；
前記第１確率密度関数値及び第２確率密度関数値を比較して前記画像中の少なくとも一部の画素に設定された注目画素が前記被写体領域及び背景領域のいずれに属するかを表す領域情報を得る比較部と；
前記領域情報を出力する出力部と；を具備する領域フィッティング装置。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、カーネル関数値を用いた、画像の平均値シフトによるフィルタリングとクラスタリング及び領域フィッティングの方法及び装置に関する。
【背景技術】
【０００２】
画像空間のような対象空間内に散らばるサンプルを均一な特徴量を持つクラスタに分けるクラスタリングの手法として、非特許文献１に開示された平均値シフト（mean shift)が知られている。平均値シフトとは、カーネル内の各サンプルの特徴量の平均で注目サンプルの特徴量を置き換える手法である。非特許文献１では、注目サンプルを中心とするある半径の円、球または楕円の固定カーネルを用いる。平均値シフトは、対象空間内に散らばる複数のサンプルに対して確率密度関数の推定（一般に、これをカーネル密度推定という）を行い、推定された確率密度関数の勾配の極大位置を求める操作によって実現される。
【０００３】
カーネル密度推定については、例えば非特許文献２のpp.13-15 及び pp.75-76で解説されている。カーネル関数の例は、非特許文献２のp.43に示されている。カーネル関数やプロファイルは、カーネル密度推定や平均値シフトにおいて個々のサンプルに対する重み値として働く。
【０００４】
カーネル密度推定によって得られる確率密度関数の極大位置ごとにクラスタを割り当てれば、クラスタリングが実現できる。具体的には、対象空間内の各サンプルを注目サンプルとして平均値シフトを繰り返し、クラスタ内の各サンプルの特徴量を当該クラスタ内の代表的な特徴量に収束させる。収束後、単一の特徴量を持つ部分ごとにクラスタを割り当てる。適当な距離またはそれに準ずるパラメータが定義されている対象空間については、このようなクラスタリングが可能である。対象空間が特に画像空間であり、位置や特徴量（例えば、輝度のベクトル、またはＲＧＢのような色信号値のベクトル）に対してクラスタリングを適用する操作は、領域分割とも呼ばれる。
【非特許文献１】D. Comaniciu et al., “Mean Shift: A Robust Approach Toward Feature Space Analysis,” IEEE Trans. Pattern Analysis and Machine Intelligence, vol.24, no.5, May 2002.
【非特許文献２】B. W. Silverman, “Density Estimation for Statistics and Data Analysis,” Chapman & HALL/CRC, ISBN 0-412-24620-1, 1986.
【発明の開示】
【発明が解決しようとする課題】
【０００５】
平均値シフトでは、注目サンプルを特徴量の平均位置に移動させる操作を行う。対象空間内のサンプルの分布が偏っていなければ、このような移動は起こらない。対象空間が２次元の位置（ｘ座標及びｙ座標）と１次元の特徴量（輝度）で定義される３次元の画像空間である場合、各画素を画像空間内のサンプルと考えれば、平均値シフトによるクラスタリングが可能である。このとき、例えば画像中に複数の固定カーネルがそれぞれ収まるような複数の部分領域により構成される一つの領域が存在すると、固定カーネルを用いた平均値シフトでは部分領域毎に収束点が生成される。その結果、単一領域に対して複数のクラスタが割り当てられてしまう（過剰分割と呼ぶ）。これは一つと見なせる領域に一つのクラスタを割り当てる、というクラスタリングの目的に反する。カーネルの半径を増やせば、単一領域に複数のクラスタが割り当てられることは起こりにくくなる。
【０００６】
対象空間内に複数のサンプル集合が存在し、クラスタリングによって各サンプル集合に対してクラスタを割り当てることを期待することを考える。このとき、カーネルの半径が小さければ、平均値シフトにより個々のサンプル集合内のサンプルはそれぞれの集合内の一箇所に収束し、各集合に対してそれぞれクラスタが割り当てられる。しかし、カーネルの半径が大きすぎると、平均値シフトにより全てのサンプルが一箇所に収束してしまう。この結果、本来分かれるべき複数のクラスタが一つに統合されてしまう（過剰統合と呼ぶ）。
【０００７】
さらに具体的に、例えば画像検索のために正しい輪郭のパターンを準備し、一般化ハフ変換（例えば、D. H. Ballard, Generalizing the Hough Transform to Detect Arbitrary Shapes, Pattern Recognition 1981, vol.13, no.2, pp.111-122 参照）により領域境界がエッジパターンと一致するものを探す場合を考える。このとき、領域境界が失われれば（過剰統合）、対象物の検出に失敗し、不要な領域境界（過剰分割）があれば、対象物でないものが誤って検出されてしまうという問題がある。
【０００８】
本発明の主な目的は、過剰分割や過剰統合の問題を回避するような、カーネル形状の決定を可能とすることにある。また、このようにして決定されるカーネル形状に従って平均値シフトフィルタリング、平均値シフトクラスタリング及び領域フィッティング等の画像処理を可能とすることにある。
【課題を解決するための手段】
【０００９】
本発明の第１の態様によると、任意のＮ次元ベクトル（Ｎは自然数）の第１画素と、γ次元の位置及びη次元の特徴量を持つＮ次元空間（γは２以上の自然数，ηは１以上の自然数）内の複数の第２画素とを結ぶパスを設定するステップと；前記パス上で前記特徴量を積分して積分値を得るステップと；前記積分値を予め設定されたカーネル関数の変数値として前記第１画素に対応するカーネル関数値を算出するステップと；前記カーネル関数値と前記第２画素の特徴量との積和を前記カーネル関数値の総和で割ることにより平均値シフトベクトルを算出するステップと；前記第１画素を前記平均値シフトベクトルに従って移動させるステップと；前記第１画素の特徴量を移動後の第１画素の特徴量で置き換えるステップと；を具備する画像のフィルタリング方法を提供する。
【００１０】
本発明の第２の態様によると、任意のＮ次元ベクトル（Ｎは自然数）の第１画素と、γ次元の位置及びη次元の特徴量を持つＮ次元空間（γは２以上の自然数，ηは１以上の自然数）内の複数の第２画素とを結ぶパスを設定するステップと；前記パス上で前記特徴量を積分して積分値を得るステップと；前記積分値を予め設定されたカーネル関数の変数値として前記第１画素に対応するカーネル関数値を算出するステップと；前記カーネル関数値と前記第２画素の特徴量との積和を前記カーネル関数値の総和で割ることにより平均値シフトベクトルを算出するステップと；前記第１画素を前記平均値シフトベクトルに従って移動させるステップと；前記第１画素の特徴量を移動後の第１画素の特徴量で置き換えるステップと；前記移動後の第１画素のうち特徴量の類似度が予め定められた閾値以上の複数の画素に対して同一のクラスタを割り当てるステップと；を具備する画像のクラスタリング方法を提供する。
【００１１】
本発明の第３の態様によると、２次元の位置と１次元の特徴量を持つ複数の画素からなる画像を取得するステップと；前記画像に対して被写体領域及び背景領域をそれぞれ概略的に表す初期第１領域及び初期第２領域を設定するステップと；任意のＮ次元ベクトル（Ｎは自然数）の第１画素と、前記画像中の複数の第２画素とを結ぶパスを設定するステップと；前記パス上で前記特徴量を積分して積分値を得るステップと；前記積分値を予め設定されたカーネル関数の変数値として前記第１画素に対応するカーネル関数値を算出するステップと；前記カーネル関数値のうち前記初期第１領域に属する各画素における第１カーネル関数値を加算して第１確率密度関数値を得るステップと；前記カーネル関数値のうち前記初期第２領域に属する各画素における第２カーネル関数値を加算して第２確率密度関数値を得るステップと；前記第１確率密度関数値及び第２確率密度関数値を比較して前記画像中の少なくとも一部の画素に設定された注目画素が前記被写体領域及び背景領域のいずれに属するかを表す領域情報を得るステップと；前記領域情報を出力するステップと；を具備する画像の領域フィッティング方法を提供する。
【００１２】
本発明の第４の態様によると、任意のＮ次元ベクトル（Ｎは自然数）の第１画素と、γ次元の位置及びη次元の特徴量を持つＮ次元空間（γは２以上の自然数，ηは１以上の自然数）内の複数の第２画素とを結ぶパスを設定する設定部と；前記パス上で前記特徴量を積分して積分値を得る積分部と；前記積分値を予め設定されたカーネル関数の変数値として前記第１画素に対応するカーネル関数値を算出する第１算出部と；前記第２サンプルの特徴量を抽出する抽出部と；前記特徴量を記憶する記憶部と；前記カーネル関数値と前記第２サンプルの特徴量との積和を前記カーネル関数値の総和で割ることにより平均値シフトベクトルを算出する第２算出部と；前記第１サンプルを前記平均値シフトベクトルに従って移動させ、前記第１サンプルの特徴量を移動後の第１サンプルの特徴量で置き換える置き換え部と；を具備する画像のための平均値シフトフィルタを提供する。
【００１３】
本発明の第５の態様によると、２次元の位置と１次元の特徴量を持つ複数の画素からなる画像を取得する取得部と；前記画像に対して被写体領域及び背景領域をそれぞれ概略的に表す初期第１領域及び初期第２領域を設定する第１設定部と；任意のＮ次元ベクトル（Ｎは自然数）の第１画素と、前記画像中の複数の第２画素とを結ぶパスを設定する第２設定部と；前記パス上で前記特徴量を積分して積分値を得る積分部と；前記積分値を予め設定されたカーネル関数の変数値として前記第１画素に対応するカーネル関数値を算出する第１算出部と；前記カーネル関数値のうち前記画像の被写体領域を概略的に表す初期第１領域に属する各画素における第１カーネル関数値を加算して第１確率密度関数値を得る第１加算部と；前記カーネル関数値のうち前記画像の背景領域を概略的に表す初期第２領域に属する各画素における第２カーネル関数値を加算して第２確率密度関数値を得る第２加算部と；前記第１確率密度関数値及び第２確率密度関数値を比較して前記画像中の少なくとも一部の画素に設定された注目画素が前記被写体領域及び背景領域のいずれに属するかを表す領域情報を得る比較部と；前記領域情報を出力する出力部と；を具備する領域フィッティング装置を提供する。
【発明の効果】
【００１４】
本発明に従ってカーネル関数値を算出すると、空間内に散らばるサンプルの分布やカーネルの位置に適したカーネル形状を設定することができるため、平均値シフトを行った場合、一つと見なせる領域に単一の収束点が形成されることによって、過剰分割や過剰統合の問題を避けることができる。
【発明を実施するための最良の形態】
【００１５】
以下、図面を参照して本発明の実施形態について説明する。以下の説明では、まずカーネル関数値の算出（カーネル形状の決定）について述べ、次いでカーネル関数値算出の応用である平均値シフト（平均値シフトフィルタ）、カーネル密度推定（確率密度関数の推定）及び領域フィッティングについて述べることとする。また、以下に説明する処理はコンピュータを用いてソフトウェアにより行うことができる。すなわち、以下の処理をコンピュータに行わせるためのプログラムを提供することができる。
【００１６】
＜カーネル関数値算出＞  
図１を参照すると、本発明の一実施形態に従うカーネル関数値算出装置は、データメモリ１０、パス設定部１１、積分部１２及びカーネル関数算出部１３を有する。このカーネル関数値算出装置は、γ次元の位置とη次元の特徴量を持つＮ次元空間（γ、η、Ｎは自然数）内のサンプル集合中の複数の候補サンプルに対応するカーネル関数値を算出する。従来ではカーネルに円、球または楕円のような固定された形状を用いていたのに対し、本実施形態ではカーネルの形状を予め定められた基準に従い、空間内に散らばるサンプルに合わせて適応的に決める。ここで、カーネルの形状を決めることは、本実施形態においてＮ次元空間内の個々の候補サンプルについてカーネル関数値（候補サンプルに与える重み値）を求めることに相当する。
【００１７】
次に、図２を用いてカーネル関数値算出の処理手順を説明する。図２のフローチャートは、Ｎ次元空間内の任意の注目サンプル（空間内に実在するサンプルでなくても構わない）と複数の候補サンプル（これは空間内に存在するサンプルである）を入力として与えたときに、カーネル関数値を算出するアルゴリズムの流れを示している。Ｎ次元空間は、これに限られないが、例えば３次元の画像空間である。画像空間の場合、サンプルとは画素である。画像空間は、２次元の位置ｘ, ｙと１次元の特徴量（輝度Ｉ(x, y)）を持つ。
【００１８】
まず、予めカーネル関数ｋ(・)を設定する（ステップＳ１０１）。カーネル関数の微分であるプロファイルも、ここでいうカーネル関数の一種である。カーネル関数は、予め固定されていてもよいし、選択可能としてもよい。例えば、コンピュータプログラムにおいて、ベクトルを入力するとカーネル関数値を返すサブルーチンへのポインタを入力とすれば、カーネル関数値の算出時にカーネル関数を選択することが容易である。カーネル関数は、具体的には正規分布関数、あるいは、変数値が予め定められた値を超える場合は０を出力し、そうでない場合は予め定められた定数を出力するような関数が用いられる。
【００１９】
ステップＳ１０１によりカーネル関数ｋ(・)を設定した後、Ｎ次元空間内のサンプル集合を入力し、このサンプル集合の各サンプル（サンプルデータ）をデータメモリ１０に保存する（ステップＳ１０２）。
【００２０】
次に、それぞれがＮ次元ベクトルからなる１個の注目サンプルＸと１個の候補サンプルＸ’を入力する（Ｓ１０３）。候補サンプルＸ’は、データメモリ１０に保存されたサンプル集合中の一つのサンプルである。このためステップＳ１０３では、データメモリ１０内の候補サンプルＸ’の保存位置を示す候補サンプルアドレスを入力してもよいし、候補サンプルのベクトルを入力してもよい。注目サンプルＸがデータメモリ１０に保存されたサンプル集合中の一つのサンプルである場合、ステップＳ１０３ではデータメモリ１０内の注目サンプルＸの保存位置を示す注目サンプルアドレスを入力してもよいし、注目サンプルのベクトルを入力してもよい。
【００２１】
ステップＳ１０３により入力された注目サンプルＸ及び候補サンプルＸ’は、パス設定部１１に与えられる。パス設定部１１は、注目サンプルＸと候補サンプルＸ’とを結ぶパスを設定する（Ｓ１０４）。パスの設定は、必要に応じデータメモリ１０内のサンプルデータを参照しながら行う。
【００２２】
ステップＳ１０４により設定されたパスは、積分部１２に渡される。積分部１２は、設定されたパスを積分区間として、必要に応じてデータメモリ１０内のサンプルデータを参照しつつサンプルの特徴量（例えば、サンプルデータが画像の場合は輝度）について定積分を行い、積分値ｒ(X, X’)を算出する（ステップＳ１０５）。積分値ｒ(X, X’)は、例えば次式で与えられる。
【数１】
2008021175.tif 000003
【００２３】
ここで、Ｒ(I(v))はスキャン中のサンプルの特徴量を注目サンプルＸに対応するカーネルに含めるかどうかを判定するための関数である。Ｒ(I(v))の具体的な形については、後述する平均値シフトに関する実施形態で詳しく述べる。なお、XおよびX’はベクトルであり、数式（１）の積分記号は線積分を表す。
【００２４】
最後に、カーネル関数値算出部１３により積分値ｒ(X, X’)をカーネル関数ｋ(・)の変数値として用いて、注目サンプルＸに対応するカーネルのカーネル関数値ｋ(r(X, X’)/h2)を算出し（ステップＳ１０６）、このカーネル関数値ｋ(r(X, X’)/h2)を出力する（ステップＳ１０７）。ここで、ｈは平均値シフトでいう半径に相当し、予め定められた値とされる。以上の処理をＮ次元空間内のサンプル集合中の各注目サンプルについて行ってそれぞれカーネル関数値を算出することにより、カーネル形状を決定することができる。カーネルの形状がどのように設定されるかについては、後に詳しく説明する。
【００２５】
なお、画像空間内の位置の表現方法はｘ，ｙの２次元に限定されず、例えばｘ，ｙと時刻ｔを用いた３次元（時空間画像）や、ｘ，ｙと奥行きｚを用いた３次元、あるいはｘ，ｙ，ｚ，ｔを用いた４次元でもよい。特徴量は１次元の輝度に限定されず、例えばカラー画像の色やＸ線画像におけるＸ線検出量でもよい。カラー画像の色を表現する色空間は任意だが、例えばＲＧＢ空間、Ｌ＊ａ＊ｂ＊空間、Ｌ＊ｕ＊ｖ＊空間、あるいはＣＭＹＫ空間を用いることができる。また、特徴量としてテクスチャを表す特徴量（例えば点毎に求めた自己相関係数、ウェーブレット係数、あるいはGaborウェーブレット係数）を用いてもよい。なお、特に断らない限り、ηが２以上の場合、本実施例中の特徴量I(x,y)をベクトルとして扱い同じ式を用いる。この場合、Ｒ(I(v))はη次元の値に対応する１次元の値を出力する関数として扱う。
【００２６】
図１のカーネル関数値算出装置及び図２に示す手順は、後述する平均値シフトやカーネル密度推定に適用できる。図１のカーネル関数値算出装置を実現するハードウェアや図２のフローチャートに従ったプログラムを設計しておけば、これらを平均値シフトやカーネル密度推定に流用することができるため、開発コストを節約できる。
【００２７】
＜平均値シフトを利用した領域分割＞  
次に、上述したカーネル関数値算出の一つの応用例として、平均値シフトを利用して画像の領域分割を行う実施形態について説明する。まず、カーネルの形状を固定した場合の問題点について述べ、次に本発明を領域分割に適用する際のカーネルの形状の決定方法について述べる。領域分割の手順についてはその後で述べる。
【００２８】
前述したように、画像は２次元の位置x, yと１次元の特徴量（輝度Ｉ(x, y)）を持つ３次元空間内のサンプル集合とみなせる。図３に示すような、ある色を持つサンプル集合の領域Ｒ１が他の色を持つ背景Ｒ２中に存在するような画像を考える。図３の画像に対して、非特許文献１に従って図４に示すような円形のカーネルＫを用いて平均値シフトを行うと、カーネルＫの内部の点は図５の矢印で示す平均値シフトベクトルに従って移動する。すなわち、図５の左側の部分領域Ｒ１１内の点はＲ１１の中心方向に移動し、右側の部分領域Ｒ１２内の点はＲ１２の中心方向に移動する。この結果、平均値シフトを繰り返すと図６に示すように部分領域Ｒ１１，Ｒ１２毎に収束点Ｐ１１，Ｐ１２が形成される。すなわち、単一領域Ｒ１に対して複数の収束点Ｐ１１，Ｐ１２が形成されてしまう。従って、平均値シフトを用いてクラスタリングを行うと、単一の領域Ｒ１に対して複数のクラスタが割り当てられてしまう。これが前述の過剰分割であり、画像に対するクラスタリングである領域分割にとって望ましくない現象である。クラスタリングの目的は、単一領域に一つのクラスタを割り当てることであるからである。図４に示したカーネルＫの半径ｒを大きくすれば、このような過剰分割は起こりにくくなる。
【００２９】
ここで、例えば図７のように同じ色を持つサンプル集合の複数の領域Ｒ１Ａ，Ｒ１Ｂが存在する画像を考える。この場合、クラスタリングによって各領域Ｒ１Ａ，Ｒ１Ｂが別々のクラスタに分けられることが期待される。カーネルＫの半径ｒが小さければ、平均値シフトにより領域Ｒ１Ａの各サンプルと領域Ｒ１Ｂ内の各サンプルはそれぞれ別の点に収束する結果、領域Ｒ１Ａと領域Ｒ１Ｂに別々のクラスタが割り当てられる。しかし、図８のようにカーネルＫの半径ｒが大きすぎると、図９の矢印で示すように平均値シフトベクトルが形成されるため、全てのサンプルが一点に収束する。その結果、本来分かれるべき２つのクラスタが一つに統合されてしまう。これが前述した過剰統合である。
【００３０】
このように平均値シフトの結果は、カーネルの半径ｒにより大きく左右される。半径ｒをチューリングすればよいが、それを人手で行うには相応の経験が必要とされ、チューニングに手間がかかる。また、最適な半径ｒは画像内の位置によって異なる。さらに、サンプルの分布が例えば空間内でのサンプルの分布が円や球から大きく外れて場合などにおいては、半径ｒをどのように設定しても過剰分割と過剰統合を同時に抑制できないことがある。複数の半径ｒで試してみて最適な半径を選択したり、カーネルの形状を円や球以外の形状（例えば、楕円あるいは超楕円）としたりすることも考えられる。しかし、サンプルが予め想定した形状のカーネルから大きく外れて分布している場合に過剰分割と過剰統合の少なくとも一方が起こるという問題は解決できない。
【００３１】
本発明の一実施形態では、カーネルの形状を円形のような固定された特定の形状に限定せず、カーネルの位置毎に位置に適した形状を設定することによって、例えば図１０に示すようにカーネルＫの形状を領域Ｒ１の全体にわたって広げる。このとき平均値シフトベクトルは、図１１に示すように領域Ｒ１内の全ての点がカーネルＫの中心を向く。このため図１２のように単一の領域Ｒ１内に一つのみの収束点Ｐが形成され、収束点に割り当てられるラベルはただ一つとなるので、上述した過剰分割や過剰統合の問題は避けられる。
【００３２】
以下に、本発明の一実施形態に従って図１０のようにカーネルＫの形状を領域Ｒ１の全体にわたって広げる方法の具体例を示す。図１３に示すように、各注目サンプルＸについて、注目サンプルＸの特徴量（例えば輝度）と同一の特徴量を持つ候補サンプルＸ’の領域（同一特徴量領域という）を探索する。図１３の矢印は、注目サンプルＸの周囲複数の方向（図１３では、上、下、左、右、左上、右上、左下、右下の８方向）における同一特徴量領域を示し、矢印の長さｒは各方向における同一特徴量領域の半径を示している。次に、図１４のように図１３の同一特徴量領域の先端を結ぶ線分（点線矢印で示す）の内部をカーネルと定める。カーネルがガウシアン関数のように無限の長さを持つ場合には、図１４の点線矢印の内部（カーネル内部）でカーネル関数値が大きくなるようにする。
【００３３】
関数の値がｘが１未満であれば１、それ以外であれば０となる一様なカーネルには名前がつけられている。ｋ(x)が一様なカーネルのとき、k(x)はUniform kernelと呼ばれる。また、カーネル関数の微分であるｇ(x)が一様なカーネルのとき、g(x)はEpanechnikovプロファイルと呼ばれる（非特許文献１参照）。ここでは、平均値シフトの算出に用いるカーネル関数の微分であるｇ(x)がEpanechnikovカーネルの場合を考える。このとき、候補サンプルＸ’が注目サンプルＸの特徴量と同一の特徴量を持つかどうかを判断する方法として、例えば以下の方法を用いることができる。候補サンプルの位置をＸ＝(x, y)、注目サンプルＸの特徴量をＩ(X)、候補サンプルＸ’の特徴量をＩ(X’)としたときに、注目サンプルＸから候補サンプルＸ’までのパスｖを考える。式（１）に示した積分値ｒ(X, X’)が次式（２）を満たす範囲内の候補サンプルは、カーネルＫの内部のサンプルであるとする。
【数２】
2008021175.tif 000004
【００３４】
ここで、ｈは平均値シフトでいう半径に相当し、予め定められた値とされる。
【００３５】
すなわち、式（２）は積分値ｒ(X, X’)をｈ2で割った値が閾値（ここでは１）に満たない候補サンプルはカーネルＫの内部のサンプルであることを表している。従って、式（１）の積分値をカーネル関数の変数値として用いることにより、カーネルの形状を決定することができる。
【００３６】
式（１）の積分値ｒ(X, X’)は、例えば次式（３）に示すように、注目サンプルＸの特徴量Ｉ(X)とパスｖ上の各サンプルの特徴量Ｉ(v)との差分（ここでは絶対値差）の積分を行うことで求められる。
【数３】
2008021175.tif 000005
【００３７】
図１５のように注目サンプルＸを積分開始点、候補サンプルＸ’を積分終了点とし、ＸからＸ’までの間のパスｖ上に、積分に利用されるサンプルが等間隔で離散的に存在していれば、式（３）は次式（４）のように変形される。
【数４】
2008021175.tif 000006
【００３８】
ただし、Ｘi（ｉ＝１，２，…，ｎ）はパスｖ上のサンプルを表す。
【００３９】
積分に利用されるサンプルが等間隔でない場合は、式（３）においてＩ(X)とＩ(v)の差をサンプル間の距離で割ってもよいし、またＩ(X)とＩ(v)の差をパスｖの長さで正規化してもよい。
【００４０】
積分値ｒ(X, X’)としては、例えば次式（５）に示すような、パスｖ上のサンプルの特徴量Ｉ(v)の変化量を用いることもできる。
【数５】
2008021175.tif 000007
【００４１】
図１５のようにＸとＸ’間のパスｖ上の積分に利用されるサンプルが等間隔で離散的に存在していれば、式（５）は次式のように変形される。
【数６】
2008021175.tif 000008
【００４２】
積分に利用されるサンプルが等間隔でない場合は、パスｖ上のサンプルの特徴量Ｉ(v)の変化量を間隔で割った式を用いてもよく、また特徴量Ｉ(v)の変化量をパスｖの長さで正規化してもよい。
【００４３】
上述した方法でカーネルの形状を決めると、カーネルは、注目サンプルの特徴量と近い特徴量を持ち、かつ注目サンプルと連結されている領域（パス中に別の領域とみなすべき部分が含まれていない領域）に一致しやすくなる。ここで、単にスキャンを行って式（１）（２）を満たす領域をカーネル内部とする方法ではなく、上述のような方法を用いることの利点は、ノイズに対する耐性を高くできるためである。
【００４４】
カーネル関数がガウシアンのように無限の長さを持つ場合は、ｒ(X, X’)／ｈ2をカーネル関数値とすることで、Epanechnikovプロファイル（‖x‖が一定値未満であればｋ（‖ｘ‖）＝１、そうでなければｋ（‖ｘ‖）＝０とするプロファイル）を用いた先の処理と同じように、注目サンプルＸの特徴量と同一の特徴量を持つと判断できる範囲内でのカーネル関数値を増やすことができる。
【００４５】
本発明の一実施形態に従うカーネル関数を用いた領域分割は、後述する平均値シフトフィルタの最後にクラスタの割り当てを行うことで実現できる。詳細は平均値シフトフィルタの説明で述べるが、その概要は次のとおりである。
【００４６】
ステップ１：画像上のすべての画素について、後述する平均値シフトフィルタを適用して（γ＋η）次元空間上で画素の位置を移動させる。このとき、各画素の移動元と移動先の対応関係がわかるように、移動元を記録しておく。
【００４７】
ステップ２：ラベル値を保持するメモリを確保し、すべての移動後の画素に対して不明なラベル値を割り当てる。
【００４８】
ステップ３：不明なラベル値を持つ移動後の画素を１つ選び、その画素に新しいラベル値を割り当てる。
【００４９】
ステップ４：ステップ３で選んだ画素に対して（γ＋η）次元空間上で近い位置にあり、かつ不明なラベル値を持つ画素をすべて探し、それらの画素にstep 3で選んだ画素と同じラベル値を割り当てる。ここで、ステップ３で選んだ画素と近い画素を探すには、例えば、ステップ３で選んだ画素との（γ＋η）次元空間上でのユークリッド距離がしきい値よりも小さい画素を探せばよい。
【００５０】
ステップ５：不明なラベル値を持つ画素がなくなるまで、ステップ３〜ステップ４を繰り返す。
【００５１】
ステップ６：対応する移動元の各画素に対しラベル値を割り当てたラベル画像を出力する。ステップ４でユークリッド距離によるしきい値判定を行うのは、実画像ではノイズの影響により、同じ点への収束が期待される画素であっても、収束点が完全に一致しないためである。この判定に用いる距離はユークリッド距離（Ｌ２距離）でなくてもよい。例えば、Ｌ１距離やＬｎ距離が予め定めた閾値以下の点を同一の点と見なすようにしても良い。
【００５２】
以下、上述した領域分割の効果について説明する。図１６Ａはオリジナル画像であり、自動車が道路上を走っている様子を後方から撮影して得られた画像である。図１６Ｂは、図１６Ａの画像に対して従来の平均値シフトによる領域分割を行った結果を示し、図１６Ｃは本実施形態に基づく領域分割を行った結果を示している。図１６Ｂでは、１６１で示されるように道路上の白線の一部が消滅し、また１６２で示されるように道路表面が不必要に分割されている。１６１で示す現象は過剰統合であり、１６２で示す現象は過剰分割である。これに対して本実施形態の領域分割によると、図１６Ｃの１６３で示すように白線が消えることはなく過剰統合が効果的に抑制されており、また１６４で示すように過剰分割も抑制されている。
【００５３】
＜平均値シフトフィルタ＞  
次に、以上述べたカーネル関数値算出手順を利用する平均値シフトフィルタについて説明する。図１７に示されるように、本発明の一実施形態に従う平均値シフトフィルタは、図１に示したカーネル関数値算出装置と同様の要素である、データメモリ（ここでは画像メモリ）１０、パス設定部１１、積分部１２及びカーネル関数値算出部１３に加えて、注目アドレス記憶部２１、候補アドレス記憶部２２、注目アドレス更新部２３、候補アドレス更新部２４、特徴量抽出部２５、積和演算部２６、積和中間値記憶部２７、特徴量画像更新部２８、特徴量画像メモリ２９及びシーケンス制御部３０を有する。平均値シフトフィルタはＮ次元空間が画像空間以外の場合にも適用が可能であるが、ここでは画像空間の場合を例にとり説明する。
【００５４】
以下、図１８を用いて３次元（２次元の位置及び１次元の輝度）の画像空間内で平均値シフトフィルタの処理を行う手順について説明する。図１７の処理は、シーケンス制御部３０により制御される。
【００５５】
まず、処理対象の入力画像データの各画素（サンプル）の輝度値を入力して画像メモリ１０に保存する（ステップＳ２０１）。入力画像データの画素のうち、次に処理すべき注目画素（注目サンプル）を設定する（ステップＳ２０２）。具体的には、予め定めた画素の処理順序（例えば、ラスタスキャン順）の先頭の画素を注目画素とし、注目画素を示すアドレスを設定し、注目アドレス記憶部２１に記憶させる。
【００５６】
次に、注目画素のためのカーネルを考えた場合の各候補画素に与える重み、すなわちカーネル関数値ｇ(r(X, X’)/h2)をステップＳ２０３〜Ｓ０６により求める。まず、注目画素に対して、重みを求めるべき候補画素のアドレスを設定し（ステップＳ２０３）、候補アドレス記憶部２２に記憶する。候補アドレス記憶部２２に記憶された候補画素のアドレスは、シーケンス制御部３０による制御のもとで候補アドレス更新部２４によって更新される。
【００５７】
ステップＳ２０４〜Ｓ２０６の処理は、図２中のステップＳ１０４〜Ｓ１０６と同様である。すなわち、パス設定部１１によって注目画素と候補画素とを結ぶパスを設定する（Ｓ２０４）。パスの設定は、必要に応じ画像メモリ１０内の画像データを参照しながら行う。ステップＳ２０４により設定されたパスは、積分部１２に渡される。積分部１２は、設定されたパスを積分区間として、必要に応じて画像メモリ１０内の画像データを参照しつつ輝度について定積分を行い、積分値ｒ(X, X’)を算出する（ステップＳ２０５）。次に、カーネル関数値算出部１３により積分値ｒ(X, X’)をカーネル関数の変数値として用いて、注目画素に対応するカーネルのカーネル関数値ｇ(r(X, X’)/h2)を算出する（ステップＳ２０６）。
【００５８】
次に、カーネル関数値ｇ(r(X, X’)/h2)を用いてカーネル内で重み付き平均を行って平均値シフトベクトルを求め（ステップＳ２０７）、さらに平均値シフトベクトルに従って注目画素の位置を移動させる（ステップＳ２０８）。言い換えれば、ステップＳ２０８では注目画素の位置及び特徴量を重み付き平均に変更する。これは特徴量画像更新部２８が特徴量画像メモリ２９の内容を更新することで行う。このとき特徴量抽出部２５により候補画素の特徴量（例えば、２次元の位置と１次元の輝度からなる３次元のベクトル）を抽出し、積和演算部２６が積和中間値記憶部２７の内容を参照しつつ各注目画素に対応するカーネル関数値ｇ(r(X, X’)/h2)と特徴量との積和を求め、この積和を全注目画素に対応するカーネル関数値の総和で割ることにより、平均値シフトベクトルが求まる。
【００５９】
ステップＳ２０３〜Ｓ２０８の処理をシーケンス制御部３０によりステップＳ２０８で所定回数繰り返されたと判断されるまで繰り返し行う。所定回数繰り返された時点でステップＳ２０８により最終的に得られる移動後の注目画素を保存する（Ｓ２１０）。すなわち、移動後の注目画素のアドレスを注目アドレス記憶部２１において記憶する。
【００６０】
ステップＳ２１１において全ての画素が注目画素として処理されるまで、注目画素を次の画素にずらしてステップＳ２０２〜Ｓ２１０の処理を繰り返す。注目画素のアドレスは注目アドレス記憶部２１に記憶されており、注目アドレス更新部２３により更新される。
【００６１】
注目画素毎に特徴量をステップＳ２０１〜Ｓ２１１で注目画素を移動した後の特徴量、すなわち重み付き平均値で置き換え（Ｓ２１２）、得られた特徴量画像を特徴量画像メモリ２９から出力する（Ｓ２１３）。ステップＳ２１２の処理は、積和演算部２６により求められた平均値シフトベクトルに従って特徴量画像更新部２８が特徴量画像メモリ２９の内容を更新することで行う。
【００６２】
なお、ステップＳ２０９で繰り返しを行わず、ステップＳ２０１〜Ｓ２１１の処理全体を繰り返し実行してもよい。ステップＳ２０９で繰り返しを行い、その後ステップＳ２０２〜Ｓ２１１の処理全体を繰り返し実行してもよい。この例では各注目画素が移動する前の位置や画素値を利用して注目画素の移動先を求めているが、注目画素の移動後の位置や画素値を利用してもよい。
【００６３】
図１８のような手順により、各注目画素は３次元空間内の幾つかの点に収束する。移動元の各点の位置における輝度値を収束後の輝度値で置き換えることで、ノイズ除去とエッジ鮮鋭を同時に行う平均値シフトフィルタを実現することができる。また、先に述べたように、各収束点にラベルを割り当てれば、画像の領域分割（対象が画像でない場合はクラスタリング）ができる。
【００６４】
注目画素に対応するカーネルを求める方法は種々考えられる。例えば、短い処理時間で良好なカーネルを求めたい場合は、次の方法を用いることができる。
【００６５】
まず、注目画素を開始点として、周囲４方向（９０°間隔の方向）または周囲８方向（４５°間隔の方向）に直線を引く。これらの各直線上で画素を終了点としたときに式（２）が満たされなくなる最初の位置を求めると、カーネルの境界となるべき点が４個あるいは８個得られる。これらカーネルの境界となるべき点を線分で結び、カーネルを形成する。より処理時間を短くしたければ、各開始点と終了点を結んだ線分上の画素のみをカーネルとしてもよい。処理時間よりも精度を重視したければ、注目画素を開始点として引く直線の本数を増やしてもよい。逆に、直線による近似を行いたくない場合には、例えば次の手順を用いる。
【００６６】
ステップＡ１：注目画素を未処理画素セット（画素位置の配列を記録できるバッファ）に所属させる。処理済画素セット（画素位置の配列を記録できるバッファ）を空に初期化する。画素毎に積分値を記録できる積分値バッファを準備する。
【００６７】
ステップＡ２：未処理画素セットから一つの画素を取り出し（当該画素を読み込んで、未処理画素セットから削除）、当該画素（処理中画素という）を処理済画素セットに追加する。
【００６８】
ステップＡ３：処理中画素の４近傍あるいは８近傍の画素（近傍画素という）のそれぞれに対し、処理済画素セットに記録された画素における積分値と近傍画素における被積分関数値との複数の加算値を求める。当該加算値の最小値を各近傍画素における積分値として積分値バッファに記録する。もし、積分値バッファに既に記録された値があれば、上書きをする。被積分関数値は、例えば次式（７）または（８）により求められる。
【数７】
2008021175.tif 000009
【００６９】
【数８】
2008021175.tif 000010
【００７０】
ステップＡ４：４つまたは８つの近傍画素を用いて次の操作をする。算出された積分値（これは近傍画素におけるｒ(X, X’)/h2で、単調増加である）がｇ(r(X, X’)/h2)を０とする値（すなわち、カーネルの半径を超えた値）でなく、かつ処理済画素セットに含まれていなければ、近傍画素を未処理画素セットに追加する。そうでない場合は、近傍画素に対しては何もしない。
【００７１】
ステップＡ５：未処理画素セットが空であれば終了する。そうでなければステップＡ２に戻る。
【００７２】
画像空間内で位置の次元γ＝２の場合の例として、注目画素を開始点として周囲４方向あるいは周囲８方向に対して操作を行う例を示した。この操作は２次元に限らず、３次元以上（γが３以上）の場合に拡張することもできる。例えば、３次元空間では、注目画素を中心として９０°間隔の方向に対し操作する場合は周囲６方向への操作に、４５°間隔の方向に対し操作する場合は周囲２６方向への操作になる。
【００７３】
画像空間は、一般的なレンズとカメラを用いた撮像装置により得られる画像空間のような、格子状にサンプリングされた均一な空間である必要はない。例えば、画素が６角形であってもかまわないし、歪みのあるレンズ（例えば、全方位センサ）を用いた撮像装置により得られる画像であってもよい。その場合には、４近傍や８近傍といった「近傍」の定義を画素の形状に合わせて変更する必要がある。例えば、注目画素に隣接する全ての画素を「近傍」とする。
【００７４】
ところで、上記実施形態では任意形状のカーネルを平均値シフトに用いた例について述べたが、任意形状のカーネルの応用例はこれに限られない。背景技術などでも説明したように、平均値シフトとはカーネル密度推定で得られる密度関数の勾配の極大値を求めるための演算の一手法にすぎない。
【００７５】
説明の便宜上、平均値シフトを先に説明したが、上記のようなカーネル形状の決定方法はカーネル密度推定にも適用できる。例えば、先のクラスタリングの例では、平均値シフトは勾配しか評価しないため、平均値シフトでは生成されたクラスタが真のクラスタであるか、ノイズにより生成された偽のクラスタであるかはわからない。しかし、カーネル密度推定によって収束位置における確率密度関数値を求めれば、真のクラスタであれば確率密度関数値は大きく、偽のクラスタであれば確率密度関数値は小さくなる。従って、例えば予め定められた閾値と確率密度関数値とを比較することで、クラスタリングを行うことができる。具体的にカーネル密度推定を行う方法については、後述する。
【００７６】
以上の説明では、γ次元の位置及びη次元の特徴量を持つＮ次元空間の例として、主に２次元の位置と１次元の特徴量（輝度）を有する画像空間を例に述べた。これらの次元数γ及びηが別の値でよいことは既に述べた通りであり、またＮ次元空間は特に画像空間である必要もない。距離あるいはそれに準ずる量（例えば、半距離）が定義され、空間上の点がＮ次元（Ｎは２以上）のベクトルで表される空間（Ｎ次元空間）があるとき、そのＮ次元空間の次元を位置と特徴量を表す２つに分ける。例えば、先の説明はＮ＝３、γ＝２（位置が２次元）、η＝１（特徴量が輝度である１次元）の場合に相当する。このＮ次元空間の特徴量をＩ(v)とし、式（１）に示したようにパスｖを積分区間としてＸからＸ’までの積分を行い、得られた積分値をカーネル関数の変数値として代入することにより、カーネル関数値が求められる。このカーネル関数値を用いて、先の説明と同様に平均値シフトやカーネル密度推定を行うことができる。
【００７７】
＜カーネル密度推定＞  
次に、カーネル密度推定について説明する。本発明の一実施形態に従うカーネル密度推定装置は、図１９に示されるように図１に示したカーネル関数値算出装置の要素である、データメモリ１０、パス設定部１１、積分部１２及びカーネル関数値算出部１３に加えて、候補アドレス記憶部３１、候補アドレス更新部３２、積和演算部３３、積和中間値記憶部３４及びシーケンス制御部３５を有する。
【００７８】
以下、図２０を参照してγ次元の位置とη次元の特徴量を有するＮ次元空間に対して、注目サンプルＸに対応するカーネル密度推定を行う手順について述べる。図２０の手順は、シーケンス制御部３５により制御される。
【００７９】
まず、入力データとしてＮ次元空間内の全サンプルのデータを入力し、データメモリ１０に保存する（ステップＳ３０１）。一方、確率密度関数を求めたい位置及び特徴量を持つ複数の注目サンプルＸを設定する（ステップＳ３０２）。注目サンプルＸは、Ｎ次元空間内の全サンプルのうちの一つであってもよいし、Ｎ次元空間内のサンプルに属さない任意のベクトルであってもよい。図１９では、データメモリ１０に保存されたサンプルに属するサンプルのみが注目サンプルＸとして設定されるときの例として、１つの注目サンプルＸのアドレス（注目アドレス）をパス設定部１１に与えている。注目アドレスに代えて、注目サンプルのベクトルそのものをパス設定部１１に与えてもよい。
【００８０】
注目サンプルＸのためのカーネルを考え、各サンプルに与える重みであるカーネル関数値ｋ(r(X, X’)/h2を「カーネル関数値算出」の項で述べた手法に従って算出する。すなわち、まず注目サンプルＸと各候補サンプルＸ’を結ぶパスを設定する（ステップＳ３０４）。パスの設定は、パス設定部１１において必要に応じてデータメモリ１０内のサンプルのデータを参照しながら行われる。候補サンプルＸ’は、例えばＮ次元空間内の全サンプルとされる。処理の高速化のために、候補サンプルＸ’をＮ次元空間内の全サンプルに対して均等の確率でランダムサンプリングを行って得られる複数のサンプル、あるいはγ次元空間で注目サンプルＸとの距離が予め設定した適当な半径ｒ以内に存在するサンプルとしてもよい。
【００８１】
次に、ステップＳ３０４により設定されたパスを積分区間として、積分部１２により必要に応じてデータメモリ１０内のデータを参照しつつ定積分を行い、積分値ｒ(X, X’)を算出する（ステップＳ３０５）。引き続きカーネル関数値算出部１３において、予め設定されたカーネル関数(・)の変数値として積分値ｒ(X, X’)を用いることにより、カーネル関数値ｋ(r(X, X’)/h2を算出する（ステップＳ３０６）。
【００８２】
次に、カーネル関数値ｋ(r(X, X’)/h2)を用いてカーネル密度推定を行い、確率密度関数ｆ(X)を算出する（ステップＳ３０７）。対象空間全体に対するサンプルの生起確率を１とするための正規化定数をＺとおくと、確率密度関数ｆ(X)は次式で求められる。
【数９】
2008021175.tif 000011
【００８３】
カーネルの形状が注目サンプルＸに応じて変化するため、正規化定数Ｚの算出は平均値シフトと違い容易ではない。しかし、勾配や確率密度の比較などではＺは特に意味を持たず、求める必要がないことも多い。また、先に述べたクラスタの閾値判定においては、正規化定数Ｚを求めずに、与えられた各サンプルの収束点に対して次式の値Ｚf(X)を求め、Ｚf(X)を平均値で割った値を閾値との比較に利用してもよい。
【数１０】
2008021175.tif 000012
【００８４】
ステップＳ３０５〜Ｓ３０７の処理は、例えば正規化定数Ｚが必要でない場合（すなわち、式（１０）のＺf(X)を算出すればよい場合）、積和演算部３３によって積和中間値記憶部３４の内容を参照しつつ積和演算を行うことで実現される。
【００８５】
なお、推定密度関数の勾配が必要な場合であっても、例えば平均値シフトを用いずにステップＳ３０１〜Ｓ３０５の処理をＮ次元空間内の注目サンプルＸの周辺の数点に対して適用し、得られた点と最小２乗法によって推定密度関数の勾配を求めてもよい。
【００８６】
連続空間で記述された積分値を離散的な観測データから求める操作は、次のように実現される。例えば、積分値ｒ(X, X’)が式（３）に従って算出される場合、γ次元で定義されたパス上に存在するη次元のサンプルを求める。それらη次元のサンプルの隣接サンプルの特徴量の差分の大きさを求め、それらの差分の大きさを合計する。
【００８７】
一方、積分値ｒ(X, X’)が式（５）に従って算出される場合は、γ次元で定義されたパス上に存在するη次元のサンプルを求める。それらη次元のサンプルの隣接サンプルの特徴量の差分の大きさをサンプル間のパスの長さで割った値を求め、それらの値を合計する。サンプル間のパスの長さは、Ｌｎ距離（ｎは距離の次数）で求められる。言い換えれば、サンプル間のパスの長さは例えばＬ１距離（市街地距離）、Ｌ２距離（ユークリッド距離）、あるいはＬ∞距離（チェス板距離）で求められる。
【００８８】
与えられた各サンプルをＮ次元空間内の点とすると、積分値の算出に利用されるのが望ましいサンプルがパスと交点を持たないことがある。各サンプルが予め定められた形状、例えばγ次元における正方形（超正方形を含む）、長方形（超長方形を含む）、あるいは円（超円を含む）と大きさを持っているとした上でパス上の点を求めると、この問題は起こりにくくなる。
【００８９】
図２のステップＳ１０４、図１８のステップＳ２０４及び図２０のステップＳ３０４におけるパスの設定は、他の結果に依存しないため、各パスの算出を並列に行うことができる。図２のステップＳ１０５、図１８のステップＳ２０５及び図２０のステップＳ３０５における積分値の算出において、パス中の次の点における積分値はパス中の直前の積分値との差分を加算した値である。このため逐次的な処理によって積分を行えば、パス中の直前の積分値を再利用できるので、計算量を減らすことができる。特に、次の条件下では計算量を大幅に減らすことができる。
【００９０】
（条件１）対象空間を画像空間、または画像と同等に扱うことのできる空間（位置の次元においてサンプルが格子点に配置されている空間）であるとする。
【００９１】
（条件２）平均値シフトでは、そのプロファイルg(x)をEpanechnikovプロファイルに限定する。カーネル密度推定では、そのカーネル関数ｋ(x)をUniformカーネル（）に限定する。Uniformカーネルの形を次式に示すが、Epanechnikovプロファイルも全く同じ形である。
【数１１】
2008021175.tif 000013
【００９２】
なお、式（１１）中、２ヶ所の定数１は別の定数であってもよい。
【００９３】
（条件３）カーネルの形状を決めるためのパスは、与えられた点を開始点とする直線であるとする。
【００９４】
（条件４）積分値は、直線上の隣接サンプル（画素）間での特徴量の差分の大きさの総和、すなわち式（５）で与えられるとする。
【００９５】
これらの４つの条件を満たすとき、本発明の一実施形態に基づく平均値シフトやカーネル密度推定、に対応するカーネル関数値の算出において、積分値の算出を高速化できる。特に、あるスキャン方向において隣接する直前の画素で用いた積分値が存在する場合、注目画素における積分値を高速に算出できる。具体的には、例えば次の方法により注目画素における積分値を順次算出する。
【００９６】
ステップＢ１：前述した「近傍」の定義に合わせて、直線方向のセットを初期化する。例えば幅Ｗ、高さＨ、位置の座標が（ｘ，ｙ）で表される画像の４近傍であれば、直線方向セットを｛（１，０）,（０，１）,（−１，０）,（０，−１）｝とする。
【００９７】
ステップＢ２：ステップＢで初期化された直線方向セットから一つの直線方向を取り出して、次に処理する直線方向とする。
【００９８】
ステップＢ３：処理中の直線方向に対してスキャンを行ったときに、画像上で終端となる点の組（端点セットという）を全て求める。例えば直線方向が（１，０）方向のとき端点セット、つまり開始点と終了点の組は｛（０、ｊ），（Ｗ−１、ｊ）｝のように求められる（ｊは整数）。
【００９９】
ステップＢ４：端点セットから１つの組を取り出す。以下、その開始点をＡ、終了点をＢとおく。
【０１００】
ステップＢ５：注目画素をＡとする。また、カーネルの端点（画素位置）をＣ、Ｄとし、ＣとＤにＡを設定する。ここでＣはカーネルの開始点側の端点、Ｄはカーネルの終了点側の端点を表す。Ｃ及びＤの積分値を０とする。
【０１０１】
ステップＢ６：注目画素Ａの次の画素（処理中の直線方向における隣接画素、以下Ａnextと記載する）に対して、ＡからＡnextまでをパスとしたときの積分値、つまりＡnextにおける、式（７）または（８）に示した被積分関数値を求める。Ｃ，Ｄの積分値からＡnextの被積分関数値を引く。
【０１０２】
ステップＢ７：Ｄの次の画素（Ｄnextと記載）に対して、Ｄnextにおける被積分関数値を求める。もしＤの積分値が正であり、Ｄの積分値にＤnextにおける被積分関数値を加算した値を変数値とするカーネルの値を０とする閾値を超えるなら、次のステップＢ８に進む。そうでなければ、ＤをＤnextで置き換え、Ｄの積分値にＤnextにおける被積分関数値を加算してステップＢ７に戻る。なお、ＤnextがＢを超える場合はステップＢ７を飛ばす。
【０１０３】
ステップＢ８：もしＣの積分値が負でないか、Ｃの積分値を変数値とするカーネルの値が０でなければステップＢ９に進む。そうでなければ、Ｃの次の画素（Ｃnextと記載する）の被積分関数値を求め、Ｃの積分値から求めた値を引く。その後ＣをＣnextで置き換え、ステップＢ８を繰り返す。
【０１０４】
ステップＢ９：注目画素Ａ、処理中の直線方向に対するUniformカーネルの範囲としてＣとＤを記録する。
【０１０５】
ステップＢ１０：注目画素がＢに到達するまで、ステップＢ６〜ステップＢ９の処理を繰り返す。
【０１０６】
ステップＢ１１：端点セットが空になるまで、ステップＢ４〜ステップＢ１０を繰り返す。
【０１０７】
ステップＢ１２：直線方向セットが空になるまで、ステップＢ２〜ステップＢ１１を繰り返す。
【０１０８】
以上の手順により、画像内の全画素に対するUniformカーネルの範囲を高速に求めることができる（画像内の全画素がＡとなるように、ステップＢ３で端点を設定した場合）。もし画像内の全画素ではなく、関心領域内の画素におけるUniformカーネルの範囲を求めたければ、次のようにする。すなわち、ステップＢ５において、処理中の直線方向に対しＡとＤの値として関心領域内で最初の画素を設定し、ステップＢ１０において、Ｂではなく処理中の直線方向に対し関心領域内で最後の画素に到達するまでステップＢ６〜ステップＢ９の処理を繰り返すように変更する。
【０１０９】
以上で説明した方法により、クラスタリングや領域分割において、従来の平均値シフトで起こっていた過剰統合や過剰分割を抑制できる。また、カーネル密度推定により確率密度関数の値を推定すれば、クラスタとノイズとの閾値判定が可能となり、ノイズ除去もできる。
【０１１０】
＜領域フィッティング＞  
本実施形態による任意形状のカーネルを用いたカーネル密度推定は、クラスタリングや領域分割以外に、画像に対する領域フィッティングという応用にも適用が可能である。例えば、文字のような被写体を撮影して得られた画像があり、画像中の被写体の存在する領域（被写体領域）とそれ以外の領域（背景領域）が概略的にわかっているときに、被写体領域と背景領域を分離する操作を領域フィッティングという。ここで、概略的に表される被写体領域及び背景領域をそれぞれ初期被写体領域（初期第１領域）及び初期背景領域（初期第２領域）と呼び、これらの領域を区別するための情報を初期領域情報と呼ぶことにする。
【０１１１】
領域フィッティングでは、画像とこれに対応する初期領域情報を入力とし、補正された（すなわちフィッティングされた）領域情報を出力する。ここで、初期領域情報及び補正領域情報は、いずれも画像の各画素に当該画素が被写体領域であれば１、背景領域あれば０というバイナリ値をそれぞれ割り当てるバイナリマスク（初期バイナリマスク及び補正後バイナリマスクという）で表現されるものとする。
【０１１２】
先の領域分割と同様に、位置と画素値（輝度値）を並べたベクトルを考える。従来のカーネル密度推定、及び本発明の一実施形態に従う任意形状のカーネルを用いたカーネル密度推定のいずれによっても、画像中に注目画素を設定したとき、次の３ステップの手順で注目画素が被写体に属するか背景に属するかを推定できる。
【０１１３】
ステップＣ1：初期領域情報において被写体領域とされた画素のみが存在する空間で、注目画素に対しカーネル密度推定を適用する。このとき注目画素の位置と画素値を入力とし、確率密度を出力とする確率密度関数が得られる。この確率密度関数に注目画素の位置と画素値を入力すると、被写体領域に属する確率密度関数値の推定値が得られる。
【０１１４】
ステップＣ２：初期領域情報において背景領域とされた画素のみが存在する空間を考えて、ステップＣ１と同様の処理を行うことにより、背景領域に属する確率密度関数値の推定値が得られる。
【０１１５】
ステップＣ3：ステップＣ１及びＣ２で得られる２つの確率密度関数値の推定値を比較し、推定値が大きい方の領域に注目画素は属すると推定する。すなわち、被写体領域に属する確率密度関数値の推定値の方が高ければ注目画素は被写体領域に属し、拝啓領域に属する確率密度関数値の推定値の方が高ければ注目画素は背景領域に属すると推定する。
【０１１６】
上記の３つのステップＣ１，Ｃ２及びＣ３の処理を全画面あるいは画像の一部の全ての画素について、それらの画素を注目画素として適用すれば、領域フィッティングが可能である。この場合、確率密度関数値の推定値は相対値のみに意味があるため、正規化定数Ｚは考えなくてよく、例えばＺ＝１とすればよい。
【０１１７】
図２１は、上述した原理に基づく一実施形態に従う領域フィッティング装置を示し、図１に示したカーネル関数値算出装置と同様の要素である、データメモリ（ここでは画像メモリ）１０、パス設定部１１、積分部１２及びカーネル関数値算出部１３に加えて、初期バイナリマスク４０、注目アドレス記憶部４１、候補アドレス記憶部４２、注目アドレス更新部４３、候補アドレス更新部４４、特徴量抽出部４５、マスク付き加算部４６、中間値記憶部４７、確率密度関数値記憶部４８、比較部４９、及び出力バイナリマスクメモリ５０を有する。
【０１１８】
以下、図２２を参照して領域フィッティングを行う手順について述べる。まず、画像データ及び初期バイナリマスクを入力し（ステップＳ４０１）、それぞれ画像メモリ１０及びバイナリマスクメモリ４０に記憶する。画像中に領域フィッティングの処理対象領域を決める（ステップＳ４０２）。処理対象領域は画像の全画面でもよいし、例えば初期バイナリマスクにおいて０と１の値が変化する画素の近傍のみでもよい。処理対象領域は、アドレス更新部４３が注目アドレス記憶部４１に記憶されている注目画素のアドレス（注目アドレス）を更新する際に利用される。処理対象領域が画像の全画面の場合は、ラスタ順で注目アドレスを更新する。処理対象領域が画面の一部分である場合は、アドレス更新部４３と結合された図示しない処理対象領域記録部の内容を入力の画像データに従って設定し、処理対象領域記録部の内容を注目アドレスの更新の際に参照する。
【０１１９】
注目画素を処理対象領域の最初の画素（例えば、処理対象領域においてラスタ順で初めの画素）に設定する（ステップＳ４０３）。この設定は、アドレス更新部４３によって注目アドレス記憶部４１に記憶されている注目アドレス（注目画素のアドレス）を更新することで行う。次いでアドレス更新部４４により候補アドレス記憶部４２に記憶されている候補アドレス（候補画素のアドレス）を更新しつつ、パス設定部１１、積分部１２及びカーネル関数値算出部１３を用いて、注目画素のためのカーネル形状の算出、すなわちカーネル関数値の算出（ステップＳ４０４）を行う。カーネル関数値の算出は、全画素に対して、つまり被写体領域のバイナリマスクを使うことなく、先に「カーネル関数値算出」の項で説明した手順により候補アドレスを更新しつつ行う。
【０１２０】
次に、確率密度関数値の推定を行う（ステップＳ４０５）。ステップＳ４０５では、マスク付き加算部４６及び中間値記憶部４７により、初期バイナリマスク４０に記憶された初期バイナリマスクを参照して、マスク付きの加算、すなわち被写体領域毎及び背景領域毎のカーネル関数値により形状が決められたカーネル内の特徴量の加算を行うことにより、被写体領域に属する確率密度関数値及び背景領域に属する確率密度関数値（推定値）を算出する。算出された被写体領域に属する確率密度関数値及び背景領域に属する確率密度関数値は、記憶部４８に記憶される。
【０１２１】
次に、確率密度関数値記憶部４８に記憶された被写体領域に属する確率密度関数値及び背景領域に属する確率密度関数値を比較部４９で比較することによって、注目画素が被写体領域に属するか背景領域に属するかの所属推定を行い、注目画素が何れの領域に属するかを示す領域情報を出力する（ステップＳ４０６）。ステップＳ４０４〜Ｓ４０６の処理を処理対象領域内の全ての画素について完了するまで、注目画素を例えばラスタ順で更新しながら繰り返す（ステップＳ４０７〜Ｓ４０８）。領域情報は、出力バイナリマスクメモリ５０に記憶される。処理対象領域内の全ての画素について処理が終了し、全画素について領域情報が得られたら、補正後バイナリマスクメモリ５０を介して領域情報をバイナリマスクとして出力する（ステップＳ４０９）。
【０１２２】
本実施形態に従う領域フィッティングにおいて従来の円形カーネルを用いた場合に比較して優位な例を説明する。図２３に示すような２次元の位置及び１次元の特徴量（画素値）を持つ被写体２３１と背景２３２を有する画像と、図２３に対応して図２４に示すような初期被写体領域２４１及び初期背景領域２４２とを境界線２４３により区別するための初期領域情報が入力された場合を考える。また、図２３の画像中の注目画素は図２５の２５１あるいは図２６の２６１であるとする。図２５及び図２６においては、図２４の境界線２４３に対応して境界線２５３及び２６３が描かれている。図２５及び図２６に示されるように、注目画素２５１及び２６１は背景領域と判定されることを期待しているものとする。
【０１２３】
従来の技術によると、図２５に示すように注目画素２５１に対して円形カーネル２５２（図２５では２次元であるが、実際は３次元の球である）が設定される。ここで、円形カーネル２５２内の白色領域は、図２４中の初期被写体領域２４１とされた部分の方が注目画素に近く、かつその面積も広い。従って、確率密度関数値の推定値においても被写体領域の確率密度関数値の方が大きな値になるため、注目画素２５１は背景領域であるにもかかわらず被写体領域と判定されてしまう。
【０１２４】
一方、本発明の一実施形態に従うと、図２６に示すように注目画素２６１に対して図２４において初期背景領域２４２とされた部分の方が広い面積を持つようなカーネル２６２が設定される。従って、注目画素２６１は正しく背景領域にあると判定される。このように、本実施形態によると領域フィッティングを行う場合に領域の判定誤りを減らすことができ、信頼性の高い結果が得られる。これはカーネルの形状が固定でなく、画像の被写体領域及び背景領域の形状を反映するように適応的に設定されているためである。
【図面の簡単な説明】
【０１２５】
【図１】第１の実施形態に従うカーネル関数値算出装置を示すブロック図
【図２】図１の装置の処理手順を示すフローチャート
【図３】ある色を持つ領域が背景中に存在する画像の例を示す図
【図４】円形カーネルについて説明する図
【図５】図３の画像に対して図４の円形カーネルを用いた場合の平均値シフトベクトルを示す図
【図６】図５の平均値シフトベクトルに従って複数の収束点が形成される様子を示す図
【図７】同じ色を持つサンプル集合の複数の領域が存在する画像の例を示す図
【図８】図７の画像に対して半径の大きい円形カーネルを用いる様子を示す図
【図９】図７の画像に対して図８のカーネルを用いた場合の平均値シフトベクトルを示す図
【図１０】図３の画像に対して一実施形態に従って形成されるカーネルの形状の例を示す図
【図１１】図３の画像に対して図１０のカーネルを用いた場合の平均値シフトベクトルを示す図
【図１２】図１１の平均値シフトベクトルに従って単一の領域に一つの収束点が形成される様子を示す図
【図１３】注目サンプルについて注目サンプルと同一特徴量を持つ候補サンプルの存在する領域を探索する様子を示す図
【図１４】図１３の探索により求まるカーネルの例を示す図
【図１５】第１の実施形態に従ってカーネルの形状を算出するための積分動作について説明する図
【図１６Ａ】領域分割の対象のオリジナル画像の例を示す図
【図１６Ｂ】図１６Ａの画像に対して従来の平均値シフトによる領域分割を行った結果を示す図
【図１６Ｃ】図１６Ａの画像に対して第１の実施形態に従う平均値シフトによる領域分割を行った結果を示す図
【図１７】第２の実施形態に従う平均値シフトフィルタを示すブロック図
【図１８】図１７のフィルタの処理手順を示すフローチャート
【図１９】第３の実施形態に従う確率密度関数値算出装置を示すブロック図
【図２０】図１９の装置の処理手順を示すフローチャート
【図２１】第４の実施形態に従う領域フィッティング装置を示すブロック図
【図２２】図２２の装置の処理手順を示すフローチャート
【図２３】図２２の装置に入力される画像の例を示す図
【図２４】図２２の装置に図２３の画像に対応して入力される初期領域情報について説明する図
【図２５】図２３の画像に対して注目画素に対応する円形カーネルを適用した例を示す図
【図２６】図２３の画像に対して注目画素に対応する一実施形態に従って形成されるカーネルを適用した例を示す図
【符号の説明】
【０１２６】
１０・・・データメモリ（画像メモリ）
１１・・・パス設定部
１２・・・積分部
１３・・・カーネル関数値算出部
２１・・・注目アドレス記憶部
２２・・・候補アドレス記憶部
２３・・・注目アドレス更新部
２４・・・候補アドレス更新部
２５・・・特徴量抽出部
２６・・・積和演算部
２７・・・積和中間値記憶部
２８・・・特徴量画像更新部
２９・・・特徴量画像メモリ
３０・・・シーケンス制御部
３１・・・候補アドレス記憶部
３２・・・候補アドレス更新部
３３・・・積和演算部
３４・・・積和中間値記憶部
４１・・・注目アドレス記憶部
４２・・・候補アドレス記憶部
４３・・・注目アドレス更新部
４４・・・候補アドレス更新部
４５・・・特徴量抽出部
４６・・・マスク付き加算部
４７・・・積和中間値記憶部
４８・・・確率密度関数値記憶部
４９・・・比較部
５０・・・補正後バイナリマスクメモリ
【図１】
2008021175.tif 000014
【図２】
2008021175.tif 000015
【図３】
2008021175.tif 000016
【図４】
2008021175.tif 000017
【図５】
2008021175.tif 000018
【図６】
2008021175.tif 000019
【図７】
2008021175.tif 000020
【図８】
2008021175.tif 000021
【図９】
2008021175.tif 000022
【図１０】
2008021175.tif 000023
【図１１】
2008021175.tif 000024
【図１２】
2008021175.tif 000025
【図１３】
2008021175.tif 000026
【図１４】
2008021175.tif 000027
【図１５】
2008021175.tif 000028
【図１６Ａ】
2008021175.tif 000029
【図１６Ｂ】
2008021175.tif 000030
【図１６Ｃ】
2008021175.tif 000031
【図１７】
2008021175.tif 000032
【図１８】
2008021175.tif 000033
【図１９】
2008021175.tif 000034
【図２０】
2008021175.tif 000035
【図２１】
2008021175.tif 000036
【図２２】
2008021175.tif 000037
【図２３】
2008021175.tif 000038
【図２４】
2008021175.tif 000039
【図２５】
2008021175.tif 000040
【図２６】
2008021175.tif 000041