JP4142966B2 - 画像解析装置、画像解析プログラム、及び画像解析方法 - Google Patents
画像解析装置、画像解析プログラム、及び画像解析方法 Download PDFInfo
- Publication number
- JP4142966B2 JP4142966B2 JP2003076294A JP2003076294A JP4142966B2 JP 4142966 B2 JP4142966 B2 JP 4142966B2 JP 2003076294 A JP2003076294 A JP 2003076294A JP 2003076294 A JP2003076294 A JP 2003076294A JP 4142966 B2 JP4142966 B2 JP 4142966B2
- Authority
- JP
- Japan
- Prior art keywords
- reflectance
- background
- analysis
- independent component
- wavelength
- 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.)
- Expired - Fee Related
Links
- 238000010191 image analysis Methods 0.000 title claims description 18
- 238000003703 image analysis method Methods 0.000 title claims description 6
- 238000004458 analytical method Methods 0.000 claims description 45
- 238000012880 independent component analysis Methods 0.000 claims description 41
- 239000011159 matrix material Substances 0.000 claims description 27
- 238000000034 method Methods 0.000 claims description 24
- 230000000737 periodic effect Effects 0.000 claims description 19
- 238000009826 distribution Methods 0.000 claims description 14
- 238000013507 mapping Methods 0.000 claims description 10
- 230000008030 elimination Effects 0.000 claims description 8
- 238000003379 elimination reaction Methods 0.000 claims description 8
- 239000000203 mixture Substances 0.000 claims description 4
- 238000001228 spectrum Methods 0.000 description 32
- 238000010586 diagram Methods 0.000 description 5
- 230000000694 effects Effects 0.000 description 5
- 229910052500 inorganic mineral Inorganic materials 0.000 description 5
- 239000011707 mineral Substances 0.000 description 5
- 238000002310 reflectometry Methods 0.000 description 5
- 238000010521 absorption reaction Methods 0.000 description 4
- 230000004069 differentiation Effects 0.000 description 3
- 238000003384 imaging method Methods 0.000 description 3
- 239000004615 ingredient Substances 0.000 description 3
- 238000005259 measurement Methods 0.000 description 3
- 230000003287 optical effect Effects 0.000 description 3
- 238000012545 processing Methods 0.000 description 3
- 239000002689 soil Substances 0.000 description 3
- 239000004973 liquid crystal related substance Substances 0.000 description 2
- 238000004088 simulation Methods 0.000 description 2
- 239000000126 substance Substances 0.000 description 2
- 230000009466 transformation Effects 0.000 description 2
- 230000005540 biological transmission Effects 0.000 description 1
- 238000007796 conventional method Methods 0.000 description 1
- 230000007613 environmental effect Effects 0.000 description 1
- 230000002068 genetic effect Effects 0.000 description 1
- 238000001727 in vivo Methods 0.000 description 1
- 230000001678 irradiating effect Effects 0.000 description 1
- 238000012544 monitoring process Methods 0.000 description 1
- 230000002265 prevention Effects 0.000 description 1
- 238000005070 sampling Methods 0.000 description 1
- 230000035945 sensitivity Effects 0.000 description 1
- 230000005236 sound signal Effects 0.000 description 1
- 230000003595 spectral effect Effects 0.000 description 1
- 230000002087 whitening effect Effects 0.000 description 1
Images
Landscapes
- Image Analysis (AREA)
- Investigating Or Analysing Materials By Optical Means (AREA)
- Geophysics And Detection Of Objects (AREA)
- Image Processing (AREA)
Description
【発明の属する技術分野】
本発明は、画像解析装置、画像解析プログラム、及び画像解析方法に関する。具体的には、上空から撮影された土地の画像データなどが有する混合スペクトルを複数のピュアスペクトルに分離して、その被覆パターンを解析する技術に関する。
【0002】
【従来の技術】
リモートセンシングは、航空機画像等を用いることで、立ち入り困難な地域の地理的状況を、調査に赴く事なくある程度把握できるという利点を有している。このため、リモートセンシングは、農業、環境、防災等の様々な分野においての利用が期待されている。
【0003】
通常、上空から撮影した衛星画像や航空機映像では、撮像系の解像度に限界があるため、対象とする構造物と背景とが、1画素内に複数混在することが多い。即ち、画像データにおける少なくとも一部の画素値(スペクトル)は、その画素内における各対象物の値(ピュアスペクトル)を、混合した値(混合スペクトル)となっている。このように混合スペクトルとなっている画素を、ミクセルという。
【0004】
ミクセルを含む画像データから被覆パターン(各対象物の占有率)を求めるためには、混合スペクトルを複数のピュアスペクトルに分離する作業を、全てのミクセルに対して行う必要がある。
ミクセルを含む画像データを解析する場合、従来の手法では、混合スペクトルは既知な構成要素のピュアスペクトルの線形和であると仮定している(例えば、非特許文献1参照)。しかし、このように構成要素を全て既知として扱うことは現実的でなく、各構成要素もそれらの混合比も未知である場合が殆どである。そして、各構成要素が推測不可能な場合に混合スペクトルをピュアスペクトルに分離することは、極めて困難であるため、極く一部の分野でしか利用可能な精度で実現できていない。
【0005】
混合スペクトルをピュアスペクトルに分離する方法が既に実用化されている分野の1つとして、ハイパースペクトル計測を用いた鉱物等の資源探査の分野が知られている(例えば、非特許文献2、及び非特許文献3参照)。なお、ハイパースペクトル計測は、数ナノメートルの波長分解能で、対象物の反射スペクトルを計測するものであり、これにより得られる観測データをハイパースペクトルデータという。また、ここでの反射スペクトルは、数十〜数百種類の波長に対して輝度値を有するハイパースペクトルとなる。
【0006】
上述した鉱物の資源探査方法は、ハイパースペクトル計測の高波長分解能という特徴と、鉱物毎に局所的な吸収特性が異なることとを利用したものである。具体的には、線形モデルを多次元空間に展開し、スペクトル空間における局所的なピークをピュアスペクトルとして抽出する。この後、鉱物毎の吸収特性との整合性をとることにより、高精度に特定の鉱物を検出する。
【0007】
【非特許文献1】
著者:小熊宏之、島崎彦人、山形与志樹著
題名:『森林樹冠率推定におけるミクセル分類手法の有効性に関する研究』
出典:写真測量とリモートセンシング、Vol.41、No2、p.4-14
出版社:写真測量学会 2002年4月出版
【非特許文献2】
著者:J.W.Boardman
題名:『AUTOMATING SPECTRAL UNMIXING OF AVIRIS DATA USING CONVEX
GEOMETRY CONCEPTS』
出典:Fourth JPL Airborne Geoscience Workshop, JPL
Publication 93-26, Vol.1, pp. 11-14 (1993)
http://makalu.jpl.nasa.gov/docs/workshops/93#docs/4.PDF
【非特許文献3】
著者:J.W.Boardman, F.A.Kruse and R.O.Green
題名:『MAPPING TARGET SIGNATURES VIA PARTIAL UNMIXING
OF AVIRIS DATA』
出典:Fifth JPL Airborne Geoscience Workshop
JPL Publication 95-1, Vol.1, pp. 23-26 (1995)
http://makalu.jpl.nasa.gov/docs/workshops/95#docs/7.PDF
【0008】
【発明が解決しようとする課題】
上述のように、異なる局所的な吸収特性を有する物質のピュアスペクトルは、混合スペクトルから分離可能である。しかし、局所的な吸収特性を有さない物質や、マルチスペクトル(数種類〜数十種類の波長に対して輝度値を有するスペクトル)のような低波長分解能のスペクトルには、同様の手法を適用できない。
【0009】
本発明の目的は、ミクセルを含む画像データから、周期的な空間分布を有する対象物の占有率分布を高精度に推定する技術を提供することである。
【0010】
【課題を解決するための手段】
本発明者は、各構成要素が未知である混合スペクトルからピュアスペクトルを抽出する手法として、独立成分分析による周期的ハイパースペクトルデータの分離方法を適用することを考えた。ここで、独立成分分析(Independent Component Analysis、以下、ICAと略記)は、複数の音源が発する音を複数の場所で観測して、複数の観測信号に基づいて、個々の音源が発した原信号をそれぞれ推定する技術である。具体的には、
複数の観測信号に対応する混合信号ベクトルをz=(z1,・・・,zK)T 、
原信号ベクトルをs=(s1,・・・,sU)T、
混合比が要素aij を持つK行U列の行列をA 、
無相関のガウス性ノイズ成分をz0
とする。このとき、次式の関係が成り立つ。
【数1】
ICAは、観測したz からガウス性ノイズ成分z0を分離して、混合比Aと原信号ベクトルsとを同時に推定する手法である。このICAを利用して、本発明は、以下のように構成される。
【0011】
<請求項1>
請求項1の画像解析装置は、周期的な空間分布を有する対象物の反射率及び背景の反射率の混合信号を含む画像データを解析するものであり、分析部と、基準値設定部と、任意性消去部とを備えていることを特徴とする。分析部は、独立成分分析における混合比を対象物の反射率と背景の反射率との差とすると共に、独立成分分析における独立成分を対象物の占有率として、独立成分分析により混合信号からノイズを分離して、混合比及び独立成分を任意性のある解として算出する。基準値設定部は、背景の反射率の波長微分値において、背景の種類に拘わらずに同じ値とみなせるものを、基準値として求めておく。任意性消去部は、背景の反射率が基準値をとると仮定して、任意性のある解として算出された混合比及び独立成分を、それぞれ1つの解に決定する。
【0012】
<請求項2>
請求項2の画像解析装置は、請求項1記載の発明において、以下の点を特徴とする。第1に、分析部は、画像データが撮像されたときの入射光強度をI0(λ)、ノイズ成分をN(λ,x)、対象物の占有率をηV(x)、背景の占有率をηS(x)、対象物の反射率をRV(λ)、背景の反射率をRS(λ)として(ここで、(x)は画像データにおける画素の位置xの関数であることを示し、(λ)は波長λの関数であることを示し、(λ,x)はx及びλの関数であることを示す。)、
ηS(x)=1−ηV(x) とみなす。
第2に、分析部は、次式
の関係に基づいて、RV(λ)とRS(λ)との差を混合比として独立成分分析を行う。
【0013】
<請求項3>
請求項3の画像解析装置は、請求項1または請求項2記載の発明において、『基準値設定部は、対象物の反射率の微分値がゼロとなる波長λi(i=1,2,・・・,p、pは2以上の自然数) における、背景の反射率値からなる行列をRefSJとして、RefSJの波長微分で与えられる行列RefSJ’と、RefSJ’のミラー行列とを判別分析により写像して、写像後の行列の要素から基準値を選択する』ことを特徴とする。
【0014】
<請求項4>
請求項4の画像解析プログラムは、請求項1〜請求項3のいずれか1項記載の分析部、基準値設定部、及び任意性消去部として、コンピュータを機能させるための画像解析プログラムである。
【0015】
<請求項5>
請求項5の画像解析方法は、周期的な空間分布を有する対象物の反射率及び背景の反射率の混合信号を含む画像データを解析するものであって、以下の手順を有することを特徴とする。
1つは、独立成分分析における混合比を『対象物の反射率』と『背景の反射率』との差とすると共に、独立成分分析における独立成分を『対象物の占有率』として、独立成分分析により混合信号からノイズを分離して、混合比及び独立成分を任意性のある解として算出する手順である。
1つは、背景の反射率の波長微分値において、背景の種類に拘わらずに同じ値とみなせるものを、基準値として求めておく手順である。
1つは、背景の反射率が基準値をとると仮定して、任意性のある解として算出された混合比及び独立成分を、それぞれ1つの解に決定する手順とである。
【0016】
【発明の実施の形態】
以下、図面を用いて本発明の実施の形態を説明する。
【0017】
<本実施形態の構成>
図1は、本実施形態の画像解析装置のブロック図である。本実施形態は、請求項1〜請求項3に対応する。画像解析装置12は、データ取得部16と、分析部20と、基準値設定部24と、任意性消去部28と、解析結果表示部32とで構成されている。
【0018】
データ取得部16は、画像データにおけるミクセル領域(画素がミクセルである画像領域)のデータ(即ち、複数の混合スペクトルで与えられるデータ)を、分析部20に入力する(後述するステップS1参照)。
分析部20は、入力された画像データにICAを施して、対象物の占有率を任意性のある解として算出する(後述するステップS2参照)。
【0019】
基準値設定部24は、背景の反射率の波長微分において安定した値を、基準値として求めておく(後述するステップS3参照)。
任意性消去部28は、背景の反射率が前記した基準値をとると仮定して、対象物の占有率を、1つの解に決定する(後述するステップS4、S5参照)。
解析結果表示部32は、任意性消去部28の解析結果を、不図示の記録媒体に記録すると共に液晶モニタに表示する(後述するステップS6参照)。
【0020】
<本実施形態の解析対象モデル>
次に、本実施形態の解析方法の説明に先立って、本実施形態の画像解析装置12の解析対象となる周期構造の観測モデルについて説明する。
図2は、周期構造を、円形の受容野を持つセンサで観測する際のモデルの一例を示している。図において、xは、円形窓(ピクセル)の中心位置であり、円形窓は、互いにその外縁で接するように(重ならないように)移動する。
【0021】
本実施形態では、観測対象領域には対象物と背景のみが存在し、同一観測対象領域内では対象物と背景のパターンは規則的であると仮定する。この場合、ピクセル位置xでの対象物の占有率をηV(x)、ピクセル位置xでの背景の占有率を
ηS(x)とすれば、任意のピクセル位置x1において、次式が成り立つ。
ηV(x1)+ηS(x1)=1 ・・・(2)
(2)式において、占有率ηV(x)とηS(x)との差が大きい場合(対象物、または背景が観測窓内の殆どを占有している場合)、ピュアスペクトルに近くなる。この場合、混合スペクトルを分離する必要がなくなる。従って、本実施形態では、どちらの占有率もある程度の大きさがある場合を想定する。
【0022】
図3は、ピクセル位置をx軸とし、対象物及び背景の占有率をy軸としたグラフの一例である。また、図4は、図3から得られる占有率のヒストグラムの一例である。
【0023】
図4が横軸方向に非対称であることから、周期的構造を有する対象領域上を低解像度のセンサでスキャンして得られる信号のヒストグラムは、非ガウス性であるといえる。このように構造物よりも大きな受容野を持つセンサでスキャンする代わりに、低解像度の撮像系を用いて空間的に連続したデータを取得しても、同様のヒストグラムが得られる。そして、このような非ガウス性の信号は、ICAを用いて、他のガウス性ノイズから高精度で分離できる。
【0024】
<本実施形態の解析方法>
図5は、本実施形態の画像解析装置12の解析手順を示す流れ図である。以下、図に示すステップ番号に従って、画像解析装置12の解析手順を説明する。
【0025】
[ステップS1]
ユーザは、画像解析装置12に、画像データを入力する。データ取得部16は、入力された画像データにおいて、ミクセルからなり、対象物及び背景の周期的な空間分布を示す画像領域があるか否かを判定する。あると判定された場合、このミクセル領域のデータを、分析部20に入力する。
なお、データ取得部16に機能に関しては、画像データにおけるユーザが選択した画像領域のデータを、分析部20に入力するようにしてもよい。或いは、入力された画像データをそのまま、分析部20に入力するようにしてもよい。
【0026】
[ステップS2]
分析部20は、以下に述べる理論に基づいて、入力されたミクセル領域のデータにICAを施し、対象物の占有率を任意性のある解として求める。
まず、ピクセル位置x、波長λでの観測信号値をI(λ,x)とすると、観測される反射光の強度は、次式で表される。
I(λ,x)=I0(λ)×{RV(λ)×ηV(x)+RS(λ)×ηS(x)} ・・・(3)
上式において、I0(λ)は入射光スペクトル強度、RV(λ)は対象物の反射率、
RS(λ)は背景の反射率である。
【0027】
十分な解像度を得られない画像観測分野において、ICAを適用するためには、ICAの新たなモデルを構成する必要がある。そこで、本実施形態では、以下の4つに示すように、従来のICAでの各パラメータを、ミクセルを含む画像データの各パラメータに置き換えて考える。
【0028】
第1に、従来のICAでの複数の観測点の位置を、波長(λ)とする。
第2に、従来のICAでのサンプリング点を、画像データ内での画素の位置、即ち、観測窓の位置(x)とする。
第3に、従来のICAでの音の原信号を、位置(x)の観測窓内における、対象物の占有率ηV(x)及び背景の占有率ηS(x)とする。
第4に、従来のICAでの混合比を、各波長における、対象物の反射率RV(λ)と背景の反射率RS(λ)との差とする。
【0029】
なお、ICAを適用するためには、以下の3つの条件を満たす必要がある。
第1に、混合信号は、原信号(占有率ηV(x)、ηS(x))の線形和である。
第2に、原信号の確率密度分布は、非ガウス性である。
第3に、各構成要素である原信号は、それぞれ独立である。
図4で説明したように、第2の条件は満たされている。また、本実施形態では、(2)式が満たされるモデルを想定しているため、第1及び第3の条件も満たされている。
【0030】
ここで、(2)式を(3)式に代入すると、次式が得られる。
I(λ,x)=I0(λ)×[{RV(λ)−RS(λ)}×ηV(x)+RS(λ)]・・・(4)
上式から明らかなように、本実施形態の解析対象モデルでは、1つの構成要素(この例では対象物)のみを推定することと等価になる。
また、実際の観測条件下では、撮像素子の熱雑音、撮像感度のばらつき、大気ノイズ等のノイズ成分を考慮する必要がある。従って、次式のように、(4)式にノイズ項を加算する必要がある。なお、ノイズ成分は通常、ガウス性である。
【0031】
I(λ,x)=I0(λ)×[{RV(λ)−RS(λ)}×ηV(x)+RS(λ)]+Rn(λ)×ηn(x)
・・・(5)
上式において、Rn(λ)はノイズの振幅であり、ηn(x)はノイズの周期パターンである。(5)式における、これらRn(λ)とηn(x)との積は、請求項記載のノイズ成分、N(λ,x)に対応する。また、実際には占有率(被覆率)にも各対象物サイズの個体差による揺らぎが存在するため、次式のように揺らぎ成分を考慮する必要がある。
【0032】
ηV(x)=ηVtrue(x)+ΔηV(x) ・・・(6)
また、揺らぎ成分を除いた対象物の占有率ηVtrue(x)を、交流成分ηV AC(x)と直流成分ηV DCに分け、次式のように表す。なお、直流成分ηV DCは、ピクセル位置xに依存しない。
ηV(x)=ηV AC(x)+ηV DC+ΔηV(x) ・・・(7)
従って、ノイズ成分や被覆率の揺らぎを考慮した場合、(7)式を(5)式に代入することで、観測される混合スペクトルは次式のように表される。
【0033】
次に、前記したICAの基本となる(1)式に基づいて、(8)式を次式のように定義する。なお、本発明では、ICAで算出される任意性を有する解を一意に定めるために、z0をピクセル位置xに対する定数項成分(直流成分)として扱う。
【0034】
ICAを適用する場合、周期性を示す信号(xに依存する交流成分)は独立成分(原信号)sとして見なされ、周期性を持たない信号(xに依存しない直流成分)は定数項成分z0と見なされる。即ち、(9)式は、以下の2式のように、周期性を示す独立成分sと定数項成分z0とに分離可能と考えられる。
【0035】
A×s=I0(λ)×{RV(λ)−RS(λ)}×{ηV AC(x)+ΔηV(x)}+Rn(λ)×ηn(x)
・・・(10)
z0=I0(λ)×[{RV(λ)−RS(λ)}×ηV DC+RS(λ)]・・・(11)
ここで、光源の状態により観測値が異なると評価が困難である。そこで、観測信号値を、入射光強度に対する比率として求めて、求めた比率をI(λ,x)とする。このようにすれば、I0(λ)=1とみなせる(但し、以下の式でも便宜上、I0(λ)を残して記載する)。次に、(10)式に対して、次式のようにICAを適用する。
【0036】
上式において、ICAにより求まるAV AC、An AC、sV AC、sn AC は、それぞれ以下の項に対応する。
【0037】
混合比の行列Aにおいて交流成分に対応するAV AC は、
I0(λ)×{RV(λ)−RS(λ)} に対応する。
独立成分の行列sにおいて交流成分に対応するsV AC は、
{ηV AC(x)+ΔηV(x)}に対応する。
混合比の行列Aにおいてノイズ成分に対応するAn AC は、Rn(λ)に対応する。
独立成分の行列sにおいてノイズ成分に対応するsn AC は、ηn(x)に対応する。
【0038】
ICAでは、推定される複数の原信号の絶対強度と、複数の原信号の混合比を表す混合行列の各要素の大きさのとり方に任意性が残されている。即ち、(12)式を満たすという条件において、交流成分に対応する混合比と独立成分とは、ICAにより算出された1つの解(AV AC、sV AC)以外に、幾つかの解を持つことができる。この任意性は、ICAでは白色化処理のために値を正規化(平均ゼロ、大きさ1に変換)することにより生じる。
【0039】
従って、反射率RV(λ)、RS(λ)、占有率ηV(x)、ηS(x)を特定するためには、混合比AV AC及び独立成分sV ACが持つ任意性を、一意に決定する必要がある。そこで、求めるべき真の混合比をnewAV AC、真の独立成分newsV ACを、次式ように定義する。
newAV AC=K×AV AC =I0(λ)×{RV(λ)−RS(λ)}・・・(13)
newsV AC=(1/K)×sV AC ・・・(14)
上式において、Kは較正係数である。そして、次のステップでは、較正係数Kを求めて、任意性を有する混合比AV AC及び独立成分sV ACを、1つの解に定める方法を説明する。
【0040】
[ステップS3]
較正係数Kを決定するためには、混合比AV ACが一定で安定した波長λで決定する方法と、波長微分を用いて混合比AV ACがより安定した波長λで決定する方法とが考えられる。本実施形態では、後者の手法を用いる。即ち、対象物の種類に依存しない波長λの微分値を特定することで任意性を消去する。
【0041】
対象物の反射率RV(λ)は様々な変化をするが、その微分値が0となる波長はほぼ決まっている。一方、背景の反射率RS(λ)は、スペクトルが滑らかに変化する特徴がある。このため、背景の反射率RS(λ)の波長微分値は、通常一定値となる。ここで、RS(λ)の波長微分値で最も安定したものを用いてもよいが、例えば背景が土の場合、土壌種類や水分含有量等によってばらつきがある。従って、特定の波長でのRS(λ)の波長微分値のみを用いると、微分値が変動するおそれがある。
【0042】
そこで、本実施形態では、まず、背景の反射率RS(λ)の波長微分値が一定値であり、且つ、対象物の反射率RV(λ)の波長微分値がゼロである条件が成立する波長を複数選択する。ここで複数の波長の値を利用するのは、様々な種類及び状態の背景に対して安定した値を得るためである。
次に、ミラー行列による判別分析を実施する。これは、選択した複数の波長での微分値に生じる分散を小さくし、且つ、背景の反射率RS(λ)の波長微分値の絶対値をできる限り大きくして誤差影響を軽減するためである。ここでの『誤差影響』とは、『ゼロと仮定した対象物の反射率RV(λ)の波長微分値が正確にゼロではないことによる誤差』である。また、ここでの『絶対値をできる限り大きく』する方法としては、元の行列に対して、座標の原点を中心として対称位置にくるミラー行列を用いることが有効と考えられる。
【0043】
なお、背景の反射率RS(λ)の波長微分値を求める手段としては、発光波長の僅かに異なる複数のレーザ光を照射し、それぞれの波長で得られる画像間の差分をとる方法が考えられる。或いは、波長選択性のない光源を用い、撮影レンズに装着する光学フィルタを切り替えて複数の画像を撮影し、その差分から波長微分値を算出してもよい。このとき、複数の光学フィルタは、透過波長域特性が僅かに異なるものを用いる。
【0044】
以上の理論に基づいて、基準値設定部24は、背景の反射率の波長微分において安定した値(後述するnewRefSJ’、及び請求項記載の基準値に対応)を、以下のように求めておく。
まず、対象物の反射率RV(λ)の微分値が0となる波長λi(i=1,2,・・・,p) における、特定の1つの背景の反射率をRefSJと定義する。ここで、Jは背景の種類であり、J=1、2、・・・mとする。RefSJは、次式で表されるように、p行1列である。
【0045】
RefSJ={RefSJ(λ1), RefSJ(λ2),・・・, RefSJ(λp)}T・・・(15)
また、特定の1つの背景の反射率RefSJの波長微分値をRefSJ’と定義し、これを次式で表す。
RefSJ’={RefSJ'(λ1), RRefSJ'(λ2),・・・, RefSJ'(λp)}T・・・(16)
次に、上式に基づいて、背景の反射率の(p×m)個の微分値分布をClass1と定義する。Class1は、次式で表されるように、p行m列である。
【0046】
Class1={RefS1’, RefS2’,・・・, RefSm’}・・・(17)
また、Class1のミラー行列を、Class2と定義する。Class2は、次式で表される。
これらClass1、Class2に対して、統計学で用いられている判別分析を適用する。なお、判別分析は、2つのグループX、Yに予め分けられたデータがあるときに、新たに入ってきたデータがどちらのグループに属するかを判定する方法である。判別分析では、グループX、Yは、それぞれ、軸変換の行列を乗じられて、別の座標系で表されたグループX’、Y’に写像される。判別分析の効果として、グループX’、Y’は、グループX、Yよりもグループ内の分散が小さくなる。また、グループX’、Y’間の距離は、グループX、Y間の距離よりも大きくなる。
【0047】
ここで、Class1、Class2に対して軸変換の行列を乗じて写像した後の座標系の軸を、Λとする。なお、写像前のClass1、Class2の座標系の軸は、波長λである。また、写像後のClass1をClass1’とし、写像後のClass2をClass2’とする。行列Class1’の行数及び列数は、写像前と同じである。同様に、行列Class2’の行数及び列数は、写像前と同じである。そして、Class1’の各列ベクトルを、次式のように定義する。
【0048】
Class1’={newRefS1’,newRefS2’,・・・,newRefSm’}
・・・(19)
上式右辺における1(または2、3、・・・)列目の全要素は、次式において、J=1(または2、3、・・・) としたものと定義する。
newRefSJ’={newRefSJ'(Λ1),・・・, newRefSJ'(Λp)}T・・・(20)
判別分析により得られるClass1’の各要素newRefSJ’ は、背景の反射率の微分値において、背景の種類に影響を受けにくく安定したものとなる(背景の種類に拘わらずほぼ同じ値となる)。
【0049】
[ステップS4] 任意性消去部28は、例えば、J=1のときのnewRefSJであるnewRefS1’を用いて、以下のように較正係数Kを求める。なお、J=1に限らず、他の値(2〜m)でもよい。
まず、(13)式において、RS(λ)がnewRefS1’に相当し、R(λ)が後述するnewRefV1’に相当すると見なせば、次式が成り立つ。
【0050】
前記したように、対象物の反射率RV(λ)の微分値がゼロであることを仮定しているので、newRefV1’はゼロとみなしてよい。なお、newRefV1’を具体的に求めるならば、まず、RefSJ’と同様に、対象物の反射率の微分値ベクトルを次式のように定義する(次式で、J=1、2、・・・m)。
【0051】
RefVJ’={RefVJ(λ1), ・・・, RefVJ(λp)}T・・・(22)
次に、RefVJ’を1列分の要素とするp行m列の行列をGroup1とし、そのミラー行列をGroup2とする。そして、これら Group1、Group2 に対して判別分析を施して同様に写像する。newRefV1’ は、写像後のGroup1における1列目の要素に対応する。
(21)式において、AV ACと、I0(λ)と、newRefV1’と、newRefS1’とは既に求めたから、較正係数Kを算出できる。そして、任意性消去部28は、求めた較正係数Kを(13)、(14)式に代入して、混合比newAV AC及び独立成分newsV ACを算出する。
【0052】
[ステップS5]
任意性消去部28は、(9)式において、スキャン方向xについてI(λ,x)の平均値を、mean{I(λ,x)}として求める。
これにより、次式のように、交流成分ηV AC(x)は消去され、揺らぎ成分ΔηV(x)及びガウス性ノイズ成分Rn(λ)×ηn(x)はゼロに収束する。
【0053】
mean{I(λ,x)}=I0(λ)×[{RV(λ)−RS(λ)}×ηV DC+RS(λ)]・・・(23)
そこで、(23)式に、(13)式を代入する。
任意性消去部28は、上式において、較正係数Kを決定したときと同様に、RS(λ)の1つの値として例えばnewRefS1’を代入して、直流成分ηV DCを求める。そして、任意性消去部28は、以下の5つを(9)式に代入して、背景の反射率RS(λ)を求める。
【0054】
第1に、本ステップで求めた直流成分ηV DC。
第2に、ステップS4で求めた混合比newAV AC を、
I0(λ)×{RV(λ)−RS(λ)} として代入する。
第3に、ステップS4で求めたsV ACを、{ηV AC(x)+ΔηV(x)}として代入する。
第4に、ステップS2において(12)式にICAを適用したときに求めた
An ACを、Rn(λ)に代入する。
第5に、ステップS2において(12)式にICAを適用したときに求めた
sn ACを、ηn(x)に代入する。
これにより、背景の反射率RS(λ)が求まる。従って、newAV ACとして求めた
I0(λ)×{RV(λ)−RS(λ)}を用いて、対象物の反射率RV(λ)も求まる。
【0055】
[ステップS6] 解析結果表示部32は、以上の解析により得られた較正係数K、混合比newAV AC、独立成分newsV AC、直流成分ηV DC、背景の反射率RS(λ)、対象物の反射率RV(λ)を、記録媒体に記録すると共に液晶モニタに表示する。
【0056】
図6及び図7は、信号強度よりも強いガウス性ノイズが混入された、図2のようにある程度のばらつきを有する周期構造の観測画像に、上述の解析処理を施したシミュレーション結果である。図6はノイズから分離された混合スペクトルを示しており、図7は分離されたガウス性ノイズ成分を示している。なお、対象物として植物、背景として土壌を想定している。
【0057】
<本実施形態の効果>
本実施形態では、対象物の反射率RV(λ)と背景の反射率RS(λ)との差を混合比AV ACとすると共に、対象物の占有率η(x)を独立成分sV ACとする。次に、ICAにより、混合信号I(λ,x)から、混合比AV AC及び独立成分sV ACを任意性のある解として算出する(ステップS2)。
【0058】
次に、背景の反射率RS(λ)の微分値が一定値であり、且つ、対象物の反射率RV(λ)の微分値がゼロである条件が成立する波長を複数選択する。そして、選択した波長での背景の反射率RS(λ)の微分値からなる行列と、この行列のミラー行列に判別分析を実施する。これにより、背景の種類に影響を受けにくく安定した、背景の反射率の微分値newRefSJ’を求める(ステップS3)。
【0059】
そして、newRefSJを用いて、混合比及び独立成分の任意性を一意に定める較正係数Kを求める(ステップS4)。従って、画像データにおける、ミクセルからなり、対象物の周期的な空間分布を示す画像領域のデータから、対象物の反射率RV(λ)、及び対象物の占有率ηV(x)を算出できる。
【0060】
また、ICAを用いるので、画像データに含まれるガウス性ノイズ成分が大きくても、ノイズ成分を分離して、周期的な空間分布を示す対象物の占有率を高精度で推定できる。なぜなら、空間に周期的に位置する対象物の占有率分布は、非ガウス性であるからである。
本実施形態の手法は、例えば以下に挙げるような、周期的分布を持つ対象領域に適用できる。内視鏡を用いた体内組織像の観測、光学顕微鏡の分解能限界を超えたナノ周期構造の観測、周期信号に変調を与える形で書き込まれた種々の微細記憶媒体からの信号の読み出し、遺伝子情報の解析、都市部の環境モニタリングなどである。
【0061】
<本実施形態の補足事項>
なお、上述したステップS1〜S6の処理をプログラムコード化して画像解析プログラムを作成し、コンピュータ上で上述した画像解析を実施してもよい。この場合、請求項4に対応し、上述した実施形態と同様の効果を得ることができる。また、ステップS1〜S6で説明した画像解析方法は、請求項5にも対応する。
【0062】
【発明の効果】
本発明では、周期的な空間分布を有する対象物の反射率及び背景の反射率の混合信号から、ICAにより、対象物の占有率を任意性のある解として求める。そして、背景の反射率が、背景の反射率の波長微分において安定した値をとると仮定する。これにより、前述の任意性を1つの解に定める。従って、前述の混合信号から、対象物の占有率を算出できる。
【図面の簡単な説明】
【図1】本発明の画像解析装置のブロック図である。
【図2】周期構造を、円形の受容野を持つセンサで観測する際のモデルの一例を示す説明図である。
【図3】ピクセル位置をx軸とし、対象物及び背景の占有率をy軸としたグラフの一例である。
【図4】図3から得られる占有率のヒストグラムの一例である。
【図5】本発明の画像解析装置の解析方法を示す流れ図である。
【図6】ガウス性ノイズを含む周期構造の観測画像に本発明の解析処理を施して得られた、ノイズから分離された混合スペクトルの波形図である。
【図7】図6のシミュレーション結果における、分離されたガウス性ノイズ成分の波形図である。
【符号の説明】
12 画像解析装置
16 データ取得部
20 分析部
24 基準値設定部
28 任意性消去部
32 解析結果表示部
Claims (5)
- 周期的な空間分布を有する対象物の反射率及び背景の反射率の混合信号を含む画像データを解析する画像解析装置であって、
独立成分分析における混合比を『前記対象物の反射率』と『前記背景の反射率』との差とすると共に、独立成分分析における独立成分を『前記対象物の占有率』として、独立成分分析により前記混合信号からノイズを分離して、前記混合比及び前記独立成分を任意性のある解として算出する分析部と、
前記背景の反射率の波長微分値において、前記背景の種類に拘わらずに同じ値とみなせるものを、基準値として求めておく基準値設定部と、
前記背景の反射率が前記基準値をとると仮定して、任意性のある解として算出された前記混合比及び前記独立成分を、それぞれ1つの解に決定する任意性消去部と
を備えていることを特徴とする画像解析装置。 - 請求項1または請求項2記載の画像解析装置において、
前記基準値設定部は、
前記対象物の反射率の微分値がゼロとなる波長λi(i=1,2,・・・,p、pは2以上の自然数) における、前記背景の反射率値からなる行列をRefSJとして、
前記RefSJの波長微分で与えられる行列RefSJ’と、前記RefSJ’のミラー行列とを判別分析により写像して、写像後の行列の要素から前記基準値を選択する
ことを特徴とする画像解析装置。 - 請求項1〜請求項3のいずれか1項記載の前記分析部、前記基準値設定部、及び前記任意性消去部として、コンピュータを機能させるための画像解析プログラム。
- 周期的な空間分布を有する対象物の反射率及び背景の反射率の混合信号を含む画像データを解析する画像解析方法であって、
独立成分分析における混合比を『前記対象物の反射率』と『前記背景の反射率』との差とすると共に、独立成分分析における独立成分を『前記対象物の占有率』として、独立成分分析により前記混合信号からノイズを分離して、前記混合比及び前記独立成分を任意性のある解として算出する手順と、
前記背景の反射率の波長微分値において、前記背景の種類に拘わらずに同じ値とみなせるものを、基準値として求めておく手順と、
前記背景の反射率が前記基準値をとると仮定して、任意性のある解として算出された前記混合比及び前記独立成分を、それぞれ1つの解に決定する手順と
を有することを特徴とする画像解析方法。
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
US43116902A | 2002-12-04 | 2002-12-04 |
Publications (2)
Publication Number | Publication Date |
---|---|
JP2004184391A JP2004184391A (ja) | 2004-07-02 |
JP4142966B2 true JP4142966B2 (ja) | 2008-09-03 |
Family
ID=32771726
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
JP2003076294A Expired - Fee Related JP4142966B2 (ja) | 2002-12-04 | 2003-03-19 | 画像解析装置、画像解析プログラム、及び画像解析方法 |
Country Status (1)
Country | Link |
---|---|
JP (1) | JP4142966B2 (ja) |
Families Citing this family (5)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
GB0822953D0 (en) | 2008-12-16 | 2009-01-21 | Stafforshire University | Image processing |
TWI493171B (zh) * | 2013-06-27 | 2015-07-21 | Univ China Medical | System and method for analyzing tissue cells by using hyperspectral image |
JP5668157B1 (ja) * | 2014-02-05 | 2015-02-12 | 佐鳥 新 | 対象物探索装置、対象物探索プログラムおよび人命救助用探索システム |
US10038706B2 (en) * | 2014-10-31 | 2018-07-31 | Verisign, Inc. | Systems, devices, and methods for separating malware and background events |
CN117115275B (zh) * | 2023-10-25 | 2024-03-12 | 深圳明锐理想科技股份有限公司 | 畸变参数的确定方法、装置和计算机设备 |
-
2003
- 2003-03-19 JP JP2003076294A patent/JP4142966B2/ja not_active Expired - Fee Related
Also Published As
Publication number | Publication date |
---|---|
JP2004184391A (ja) | 2004-07-02 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
Cabot et al. | On the robustness of analysis techniques for molecular detections using high-resolution exoplanet spectroscopy | |
Shiklomanov et al. | Enhancing global change experiments through integration of remote‐sensing techniques | |
Rivera et al. | On the semi-automatic retrieval of biophysical parameters based on spectral index optimization | |
Klein et al. | Quantitative hyperspectral reflectance imaging | |
Savorgnan et al. | Supermassive Black Holes and Their Host Spheroids. I. Disassembling Galaxies | |
Bothwell et al. | The star formation rate distribution function of the local Universe | |
Schmidt et al. | Long term data fusion for a dense time series analysis with MODIS and Landsat imagery in an Australian Savanna | |
Ali et al. | Retrieval of forest leaf functional traits from HySpex imagery using radiative transfer models and continuous wavelet analysis | |
Crucil et al. | Assessing the performance of UAS-compatible multispectral and hyperspectral sensors for soil organic carbon prediction | |
Peón et al. | Prediction of topsoil organic carbon using airborne and satellite hyperspectral imagery | |
Chapman et al. | Understanding the nature of optically faint radio sources and their connection to the submillimeter population | |
Ji et al. | Comparison of different multispectral sensors for photosynthetic and non-photosynthetic vegetation-fraction retrieval | |
Leitão et al. | From sample to pixel: multi‐scale remote sensing data for upscaling aboveground carbon data in heterogeneous landscapes | |
Kuhlmann et al. | An algorithm for in-flight spectral calibration of imaging spectrometers | |
Lahssini et al. | Influence of GEDI acquisition and processing parameters on canopy height estimates over tropical forests | |
JPWO2019077955A1 (ja) | スペクトル分析装置およびスペクトル分析方法 | |
Bi et al. | Estimating leaf chlorophyll and nitrogen contents using active hyperspectral LiDAR and partial least square regression method | |
Campbell et al. | Mapping individual tree and plot-level biomass using airborne and mobile lidar in piñon-juniper woodlands | |
JP4142966B2 (ja) | 画像解析装置、画像解析プログラム、及び画像解析方法 | |
De Grandi et al. | Wavelet based analysis of TanDEM-X and LiDAR DEMs across a tropical vegetation heterogeneity gradient driven by fire disturbance in Indonesia | |
Taubert et al. | Deriving tree size distributions of tropical forests from LiDAR | |
dos Santos et al. | Improving the generalization error and transparency of regression models to estimate soil organic carbon using soil reflectance data | |
Atkins et al. | Integrating forest structural diversity measurement into ecological research | |
Runnholm et al. | On the evolution of the size of Lyman alpha haloes across cosmic time: no change in the circumgalactic gas distribution when probed by line emission | |
Moncholi-Estornell et al. | Enhancing solar-induced fluorescence interpretation: Quantifying fractional sunlit vegetation cover using linear spectral unmixing |
Legal Events
Date | Code | Title | Description |
---|---|---|---|
A621 | Written request for application examination |
Free format text: JAPANESE INTERMEDIATE CODE: A621 Effective date: 20051024 |
|
A711 | Notification of change in applicant |
Free format text: JAPANESE INTERMEDIATE CODE: A711 Effective date: 20070720 |
|
A521 | Request for written amendment filed |
Free format text: JAPANESE INTERMEDIATE CODE: A821 Effective date: 20070720 |
|
A977 | Report on retrieval |
Free format text: JAPANESE INTERMEDIATE CODE: A971007 Effective date: 20071031 |
|
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: 20080610 |
|
A01 | Written decision to grant a patent or to grant a registration (utility model) |
Free format text: JAPANESE INTERMEDIATE CODE: A01 |
|
A61 | First payment of annual fees (during grant procedure) |
Free format text: JAPANESE INTERMEDIATE CODE: A61 Effective date: 20080613 |
|
FPAY | Renewal fee payment (event date is renewal date of database) |
Free format text: PAYMENT UNTIL: 20110620 Year of fee payment: 3 |
|
R150 | Certificate of patent or registration of utility model |
Free format text: JAPANESE INTERMEDIATE CODE: R150 |
|
FPAY | Renewal fee payment (event date is renewal date of database) |
Free format text: PAYMENT UNTIL: 20110620 Year of fee payment: 3 |
|
FPAY | Renewal fee payment (event date is renewal date of database) |
Free format text: PAYMENT UNTIL: 20120620 Year of fee payment: 4 |
|
FPAY | Renewal fee payment (event date is renewal date of database) |
Free format text: PAYMENT UNTIL: 20130620 Year of fee payment: 5 |
|
R250 | Receipt of annual fees |
Free format text: JAPANESE INTERMEDIATE CODE: R250 |
|
R250 | Receipt of annual fees |
Free format text: JAPANESE INTERMEDIATE CODE: R250 |
|
R250 | Receipt of annual fees |
Free format text: JAPANESE INTERMEDIATE CODE: R250 |
|
R250 | Receipt of annual fees |
Free format text: JAPANESE INTERMEDIATE CODE: R250 |
|
LAPS | Cancellation because of no payment of annual fees |