JP2009014354A - 画像処理装置および画像処理プログラム - Google Patents
画像処理装置および画像処理プログラム Download PDFInfo
- Publication number
- JP2009014354A JP2009014354A JP2007173245A JP2007173245A JP2009014354A JP 2009014354 A JP2009014354 A JP 2009014354A JP 2007173245 A JP2007173245 A JP 2007173245A JP 2007173245 A JP2007173245 A JP 2007173245A JP 2009014354 A JP2009014354 A JP 2009014354A
- Authority
- JP
- Japan
- Prior art keywords
- representative
- pixel
- image
- extraction
- image processing
- 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.)
- Withdrawn
Links
- 230000003595 spectral effect Effects 0.000 claims abstract description 189
- 239000000523 sample Substances 0.000 claims abstract description 93
- 239000013074 reference sample Substances 0.000 claims abstract description 48
- 238000000605 extraction Methods 0.000 claims description 112
- 238000012545 processing Methods 0.000 claims description 112
- 238000002834 transmittance Methods 0.000 claims description 71
- 238000000034 method Methods 0.000 claims description 67
- 238000004364 calculation method Methods 0.000 claims description 24
- 238000005457 optimization Methods 0.000 claims description 24
- 238000012360 testing method Methods 0.000 claims description 23
- 238000003384 imaging method Methods 0.000 claims description 19
- 239000000284 extract Substances 0.000 claims description 13
- 230000002596 correlated effect Effects 0.000 abstract 2
- 230000000875 corresponding effect Effects 0.000 abstract 1
- 239000006185 dispersion Substances 0.000 abstract 1
- 239000000975 dye Substances 0.000 description 66
- 239000011159 matrix material Substances 0.000 description 52
- 230000008569 process Effects 0.000 description 41
- 239000000049 pigment Substances 0.000 description 25
- 238000005315 distribution function Methods 0.000 description 24
- 230000001575 pathological effect Effects 0.000 description 16
- 238000010186 staining Methods 0.000 description 16
- WZUVPPKBWHMQCE-UHFFFAOYSA-N Haematoxylin Chemical compound C12=CC(O)=C(O)C=C2CC2(O)C1C1=CC=C(O)C(O)=C1OC2 WZUVPPKBWHMQCE-UHFFFAOYSA-N 0.000 description 12
- 230000014509 gene expression Effects 0.000 description 12
- 238000011156 evaluation Methods 0.000 description 10
- 210000003855 cell nucleus Anatomy 0.000 description 9
- 230000003287 optical effect Effects 0.000 description 9
- 238000005259 measurement Methods 0.000 description 8
- 230000001419 dependent effect Effects 0.000 description 7
- YQGOJNYOYNNSMM-UHFFFAOYSA-N eosin Chemical compound [Na+].OC(=O)C1=CC=CC=C1C1=C2C=C(Br)C(=O)C(Br)=C2OC2=C(Br)C(O)=C(Br)C=C21 YQGOJNYOYNNSMM-UHFFFAOYSA-N 0.000 description 7
- 210000003743 erythrocyte Anatomy 0.000 description 7
- 230000035945 sensitivity Effects 0.000 description 7
- 210000000805 cytoplasm Anatomy 0.000 description 6
- 238000004043 dyeing Methods 0.000 description 5
- 238000005286 illumination Methods 0.000 description 5
- 230000005855 radiation Effects 0.000 description 5
- 239000000126 substance Substances 0.000 description 5
- 238000007476 Maximum Likelihood Methods 0.000 description 4
- 239000000243 solution Substances 0.000 description 4
- 210000001519 tissue Anatomy 0.000 description 4
- BMUDPLZKKRQECS-UHFFFAOYSA-K 3-[18-(2-carboxyethyl)-8,13-bis(ethenyl)-3,7,12,17-tetramethylporphyrin-21,24-diid-2-yl]propanoic acid iron(3+) hydroxide Chemical compound [OH-].[Fe+3].[N-]1C2=C(C)C(CCC(O)=O)=C1C=C([N-]1)C(CCC(O)=O)=C(C)C1=CC(C(C)=C1C=C)=NC1=CC(C(C)=C1C=C)=NC1=C2 BMUDPLZKKRQECS-UHFFFAOYSA-K 0.000 description 3
- QTBSBXVTEAMEQO-UHFFFAOYSA-N Acetic acid Chemical compound CC(O)=O QTBSBXVTEAMEQO-UHFFFAOYSA-N 0.000 description 3
- 230000005540 biological transmission Effects 0.000 description 3
- 238000006243 chemical reaction Methods 0.000 description 3
- 230000007423 decrease Effects 0.000 description 3
- 238000010586 diagram Methods 0.000 description 3
- 229940109738 hematin Drugs 0.000 description 3
- 230000003902 lesion Effects 0.000 description 3
- 238000003672 processing method Methods 0.000 description 3
- 102000053602 DNA Human genes 0.000 description 2
- 108020004414 DNA Proteins 0.000 description 2
- 235000005811 Viola adunca Nutrition 0.000 description 2
- 240000009038 Viola odorata Species 0.000 description 2
- 235000013487 Viola odorata Nutrition 0.000 description 2
- 235000002254 Viola papilionacea Nutrition 0.000 description 2
- 230000008859 change Effects 0.000 description 2
- 102000004169 proteins and genes Human genes 0.000 description 2
- 108090000623 proteins and genes Proteins 0.000 description 2
- 241000196324 Embryophyta Species 0.000 description 1
- 230000002378 acidificating effect Effects 0.000 description 1
- 230000002411 adverse Effects 0.000 description 1
- 150000001413 amino acids Chemical class 0.000 description 1
- 210000000988 bone and bone Anatomy 0.000 description 1
- 244000309464 bull Species 0.000 description 1
- 238000004891 communication Methods 0.000 description 1
- 210000002808 connective tissue Anatomy 0.000 description 1
- 239000000470 constituent Substances 0.000 description 1
- 238000007796 conventional method Methods 0.000 description 1
- 238000012937 correction Methods 0.000 description 1
- 230000003247 decreasing effect Effects 0.000 description 1
- 230000002950 deficient Effects 0.000 description 1
- 238000003745 diagnosis Methods 0.000 description 1
- 230000001747 exhibiting effect Effects 0.000 description 1
- 230000006870 function Effects 0.000 description 1
- 230000010354 integration Effects 0.000 description 1
- 239000004973 liquid crystal related substance Substances 0.000 description 1
- 238000013507 mapping Methods 0.000 description 1
- 229940050561 matrix product Drugs 0.000 description 1
- 230000007246 mechanism Effects 0.000 description 1
- 238000012986 modification Methods 0.000 description 1
- 230000004048 modification Effects 0.000 description 1
- 239000005445 natural material Substances 0.000 description 1
- 238000013188 needle biopsy Methods 0.000 description 1
- 210000000056 organ Anatomy 0.000 description 1
- 238000010827 pathological analysis Methods 0.000 description 1
- 125000002467 phosphate group Chemical group [H]OP(=O)(O[H])O[*] 0.000 description 1
- 238000011002 quantification Methods 0.000 description 1
- 238000000611 regression analysis Methods 0.000 description 1
- 230000003252 repetitive effect Effects 0.000 description 1
- 238000012552 review Methods 0.000 description 1
- 239000012192 staining solution Substances 0.000 description 1
Images
Landscapes
- Investigating Or Analysing Materials By Optical Means (AREA)
- Image Processing (AREA)
- Color Television Image Signal Generators (AREA)
- Image Analysis (AREA)
Abstract
【課題】標本上の各標本点における分光特性を、意図しないばらつきを生じさせることなく高精度に推定できること。
【解決手段】画像処理装置1は、基準標本が記録された画像であって画素ごとに基準標本上の対応する標本点における基準分光特性が対応付けられた教師画像の各画素を、少なくとも主要要素を特徴付ける主要特徴量に基づいて分類する主要要素分類部7bと、主要要素分類部7bが分類した画素群の中から複数の代表画素を抽出する代表点抽出部7cと、代表点抽出部7cが抽出した複数の代表画素に対応付けられた基準分光特性をもとに、標本の分光特性の推定に用いる基準分光特性の統計量を算出する統計量算出部7dと、を備える。
【選択図】図1
【解決手段】画像処理装置1は、基準標本が記録された画像であって画素ごとに基準標本上の対応する標本点における基準分光特性が対応付けられた教師画像の各画素を、少なくとも主要要素を特徴付ける主要特徴量に基づいて分類する主要要素分類部7bと、主要要素分類部7bが分類した画素群の中から複数の代表画素を抽出する代表点抽出部7cと、代表点抽出部7cが抽出した複数の代表画素に対応付けられた基準分光特性をもとに、標本の分光特性の推定に用いる基準分光特性の統計量を算出する統計量算出部7dと、を備える。
【選択図】図1
Description
本発明は、1以上の主要要素を含む標本を撮像した観察画像を処理する画像処理装置および画像処理プログラムに関するものである。
生体組織標本、特に病理標本では、臓器摘出によって得たブロック標本や針生検によって得た標本を厚さ数ミクロン程度に薄切した後、様々な所見を得る為に顕微鏡を用いて拡大観察することが広く行われている。中でも光学顕微鏡を用いた透過観察は、機材が比較的安価で取り扱いが容易である上、歴史的に古くから行われてきたこともあって、最も普及している観察方法の一つである。この場合、薄切された生体標本は光を殆ど吸収及び散乱せず無色透明に近い為、観察に先立って色素による染色を施すのが一般的である。
染色手法としては種々のものが提案されており、その総数は100種類以上にも達するが、特に病理標本に関しては、青紫色のヘマトキシリンと赤色のエオジンとの2つの色素を用いるヘマトキシリン−エオジン染色(以下、H&E染色と呼ぶ。)が標準的に用いられている。
ヘマトキシリンは植物から採取された天然の物質であり、それ自身には染色性は無い。しかし、その酸化物であるヘマチンは好塩基性の色素であり、負に帯電した物質と結合する。細胞核に含まれるデオキシリボ核酸(DNA)は、構成要素として含むリン酸基によって負に帯電している為、ヘマチンと結合して青紫色に染色される。尚、前述の通り、染色性を有するのはヘマトキシリンでは無く、その酸化物であるヘマチンであるが、色素の名称としてはヘマトキシリンを用いるのが一般的である為、以下それに従う。
エオジンは好酸性の色素であり、正に帯電した物質と結合する。アミノ酸やタンパク質が正負どちらに帯電するかはpH環境に影響を受け、酸性下では正に帯電する傾向が強くなる。その為、エオジン溶液には酢酸を加えて用いることがある。細胞質に含まれるタンパク質は、エオジンと結合して赤から薄赤に染色される。
H&E染色後の標本では細胞核、骨組織等が青紫色に、細胞質、結合織、赤血球等が赤色に染色され、容易に視認できるようになる。その結果、観察者は、細胞核等の組織を構成する要素の大きさや位置関係等を把握することができ、標本の状態を形態学的に判断することが可能となる。
生体組織標本の染色は、元々個体差を有する生体組織に対し、化学反応を用いて色素を固定する作業である為、常に均一な結果を得ることが難しい。具体的には、同一濃度の染色液に同一時間標本を反応させた場合でも、固定される色素の量が同程度であるとは限らない。標本によっては比較的多くの色素が固定される場合や、比較的少ない色素しか固定されない場合がある。前者では通常より濃く染色された標本となり、後者は薄く染色された標本となる。このような標本間での染色のばらつきを抑える為、施設によっては専門の技能を有した染色技師を配置している。しかし、染色技師の職人的な調整作業によって同一施設内での染色ばらつきはある程度軽減できるが、異なる施設間での染色ばらつきは依然として生じたままである。
染色ばらつきは、2つの点で問題がある。1つは、染色された標本を観察者が目視観察する場合、標本の染色状態が不揃いであることが観察者のストレスに繋がる可能性がある。重度の染色ばらつきが生じている場合、決定的な所見が見落とされる可能性がある。2つ目は、染色された標本をカメラで撮像して画像処理する場合、染色ばらつきが画像処理精度に悪影響を及ぼす可能性がある。例えば、ある病変が特定の色を呈することが判っていたとしても、標本を撮像した観察画像から自動的にその病変に対応する画像領域を抽出することが難しくなる。病変の特性による色変化を染色ばらつきが撹乱してしまうからである。
このような染色ばらつきの問題を、マルチバンド画像を用いた画像処理によって解決する手法が非特許文献1に開示されている。非特許文献1に記載された画像処理方法では、観察画像としてマルチバンド画像を取得し、このマルチバンド画像の画像データをもとに、まず標本上の各標本点における分光特性をウィナー(Wiener)推定によって推定する。つぎに、この推定した分光特性を用い、標本上の各標本点に固定された色素量を物理モデルに基づいて推定する。そして、この推定した色素量を適宜増減させて補正するとともに、この補正した色素量をもとに画像データを生成することで、染色状態を補正した観察画像を取得している。これによると、推定した色素量の補正を適切に行うことで、濃く染色された標本あるいは薄く染色された標本の観察画像を、適切に染色された標本と同等の色を呈する観察画像に補正することができる。
ところで、非特許文献1に記載の画像処理方法では、標本上の各標本点における分光特性をウィナー推定するために、その標本に含まれる主要要素を含んだ基準標本の分光特性の統計量が必要となる。ウィナー推定の場合、その統計量は、分光特性としての分光透過率の自己相関行列もしくは共分散行列に相当する。この自己相関行列もしくは共分散行列は、基準標本上の複数の標本点における分光透過率を事前に分光計によって計測した計測結果をもとに算出される。以下、本発明の実施の形態では、上記統計量に自己相関行列を用いる場合を説明し、この自己相関行列の算出に用いる基準標本上の標本点を代表点と呼ぶ。
基準標本が記録された基準標本画像をもとに基準標本の全域から代表点を均一に抽出した場合、基準標本画像上で基準標本中の主要要素を示す主要領域のうち面積が大きい主要領域に対応する主要要素からより多くの代表点が抽出される。ウィナー推定では、平均二乗誤差を最小化する手法であるため、代表点が多く集まる主要要素に対する推定精度が相対的に高くなる。すなわち、ウィナー推定では、基準標本画像上で面積が大きい主要領域に対応する主要要素に対して推定精度が高く、面積が小さい主要領域に対応する主要要素に対して推定精度が低くなる傾向を示す。
したがって、非特許文献1に記載の画像処理方法では、標本が病理標本である場合、診断で重要とされる細胞核が小さいため、細胞核の分光透過率を精度良く推定することが困難であった。これに対して、特許文献1に開示された技術では、典型的な病理標本を撮像した基準標本画像をもとに、ユーザーが各主要要素からほぼ同数の代表点を抽出するようにしているため、細胞核に対する推定精度を他の主要要素に対する推定精度と同等にすることができる。ここで、病理標本における主要要素には、細胞核のほか、細胞質、赤血球および背景等がある。
しかしながら、ユーザーが代表点の抽出を行うユーザー抽出では、抽出される代表点がユーザーによって異なるため、ウィナー推定による分光透過率の推定結果もユーザーによってばらつくこととなる。つまり、ユーザー抽出を行った場合、抽出される代表点の最適性および再現性が保証されておらず、分光透過率の推定精度にばらつきが生じるという問題があった。また、ユーザー抽出では、ユーザーが代表点を抽出するために試行錯誤せざるを得ず、その作業に多大な労力と時間を要するという問題があった。
本発明は、上記に鑑みてなされたものであって、標本上の各標本点における分光特性を、意図しないばらつきを生じさせることなく高精度に推定することができる画像処理装置および画像処理プログラムを提供することを目的とする。
上記の目的を達成するために、本発明にかかる画像処理装置は、1以上の主要要素を含む標本を撮像した観察画像の画像データをもとに、前記1以上の主要要素を含む基準標本上の複数の標本点における基準分光特性の統計量を用いて前記標本の分光特性を推定する画像処理装置において、前記基準標本が記録された画像であって画素ごとに前記基準標本上の対応する標本点における前記基準分光特性が対応付けられた教師画像の各画素を、少なくとも前記主要要素を特徴付ける主要特徴量に基づいて分類する画素分類手段と、前記画素分類手段が分類した画素群の中から複数の代表画素を抽出する代表抽出手段と、複数の前記代表画素に対応付けられた前記基準分光特性をもとに前記統計量を算出する統計量算出手段と、を備えたことを特徴とする。
また、本発明にかかる画像処理プログラムは、1以上の主要要素を含む標本を撮像した観察画像の画像データをもとに、前記1以上の主要要素を含む基準標本上の複数の標本点における基準分光特性の統計量を用いて前記標本の分光特性を推定する画像処理装置に、前記基準標本が記録された画像であって画素ごとに前記基準標本上の対応する標本点における前記基準分光特性が対応付けられた教師画像の各画素を、少なくとも前記主要要素を特徴付ける主要特徴量に基づいて分類する画素分類手順と、前記画素分類手順によって分類された画素群の中から複数の代表画素を抽出する代表抽出手順と、複数の前記代表画素に対応付けられた前記基準分光特性をもとに前記統計量を算出する統計量算出手順と、を実行させることを特徴とする。
本発明にかかる画像処理装置および画像処理プログラムによれば、標本上の各標本点における分光特性を、意図しないばらつきを生じさせることなく高精度に推定することができる。
以下、添付図面を参照して、本発明にかかる画像処理装置および画像処理プログラムの好適な実施の形態を詳細に説明する。なお、この実施の形態により本発明が限定されるものではない。また、図面の記載において、同一部分には同一符号を付して示している。
(実施の形態)
図1は、本実施の形態にかかる画像処理装置1の要部構成を示すブロック図である。この図に示すように、画像処理装置1は、色素によって染色された標本を撮像した観察画像を取得する観察画像取得部2と、基準標本の分光特性を取得する基準分光特性取得部3と、各種情報の入力、表示および記憶をそれぞれ行う入力部4、表示部5および記憶部6と、観察画像取得部2が取得した観察画像および基準分光特性取得部3が取得した分光特性に対して各種演算処理を行う画像処理部7と、これらの各部に電気的に接続され、その接続された各部の処理および動作を制御する制御部8とを備える。
図1は、本実施の形態にかかる画像処理装置1の要部構成を示すブロック図である。この図に示すように、画像処理装置1は、色素によって染色された標本を撮像した観察画像を取得する観察画像取得部2と、基準標本の分光特性を取得する基準分光特性取得部3と、各種情報の入力、表示および記憶をそれぞれ行う入力部4、表示部5および記憶部6と、観察画像取得部2が取得した観察画像および基準分光特性取得部3が取得した分光特性に対して各種演算処理を行う画像処理部7と、これらの各部に電気的に接続され、その接続された各部の処理および動作を制御する制御部8とを備える。
観察画像取得部2は、例えば図2に示すように、複数の色素によって染色された標本10が載置されたステージ11と、ステージ11を介して標本10を透過照明する光源12と、標本10からの光を集光して標本10の観察像を結像させる結像光学系13と、観察像を結像する光の波長帯域を所定範囲に制限するバンドパスフィルタ14と、観察像を撮像して観察画像を生成するカメラ16とを備える。
カメラ16は、CCD等を用いたモノクロの撮像素子上に、図3に示すようなモザイク状のRGBフィルタを配置して構成されており、R,G,Bの各フィルタ領域において、図4に示すような分光感度特性を有している。カメラ16は、各画素においてR,G,Bいずれかの成分しか撮像することができないが、例えば特許文献2に開示されているように、近傍の画素値を利用して各々不足するR,G,B成分を補間することができる。これによってカメラ16は、各画素においてR,G,B成分をすべて取得している。なお、カメラ16は、RGBフィルタを用いた構成に限らず、例えば3CCDタイプのカメラを用いることができる。その場合、各画素において補間を行うことなくR,G,B成分をすべて取得することができる。
バンドパスフィルタ14は、例えば図5−1および図5−2に示すような異なる分光透過率特性を有したバンドパスフィルタ14a,14bを含み、それぞれ回転式のターレット15に保持されている。これによって、バンドパスフィルタ14a,14bは、結像光学系13とカメラ16との間に切換自在に配置可能とされている。観察画像取得部2は、観察像を結像する光の波長帯域をバンドパスフィルタ14a,14bによって順次変化させながら撮像を行うことで、標本10の観察画像として6バンドのマルチバンド画像(6バンド画像)を取得することができる。観察画像取得部2が取得した観察画像は、制御部8に転送され、記憶部6に記憶される。
なお、ここでは観察画像取得部2は、2枚のバンドパスフィルタ14a,14bを用いて6バンドの観察画像を取得するものとしたが、2枚に限定されず3枚以上用いてもよい。それによって、より多くのバンド数のマルチバンド画像を観察画像として取得することができる。例えば、16枚のバンドパスフィルタをフィルタホイールで回転させて切り換えながら面順次方式でマルチバンド画像を撮像する技術が特許文献1に開示されている。これによると、16バンドの観察画像を得ることができる。
また、バンドパスフィルタ14は、結像光学系13とカメラ16との間に配置されるものとして説明したが、光源12からカメラ16に至る光路中のいずれの位置に配置させてもよい。さらに、観察像を結像する光の波長帯域を制限するものとして、バンドパスフィルタ14の代わりに透過帯域可変型のフィルタを用いてもよく、例えばバンドパスフィルタ14およびターレット15に替えて液晶チューナブルフィルタを用いることができる。
なお、本実施の形態では、複数の色素によって染色された標本10として、ヘマトキシリンとエオジンとの2つの色素によって染色された病理標本を用いるものとし、画像処理装置1は、その病理標本におけるヘマトキシリンの色素量と、エオジンの色素量と、赤血球の色素量との3種類の色素量を推定するものとして説明する。また、以下の説明では、これら3種類の色素をそれぞれ色素H、色素Eおよび色素Rと略称する。ただし、画像処理装置1によって推定可能な色素は、色素H、色素Eおよび色素Rに限定されるものではなく、画像処理装置1によって処理可能な標本は、病理標本に限定されるものではない。なお、色素を施さない状態であっても赤血球はそれ自身特有の色を有しており、H&E染色後は、赤血球自身の色と染色過程において変化したエオジンの色が重畳して観察される。よって、正確には両者を併せたものが色素Rとなる。
基準分光特性取得部3は、図2に示した観察画像取得部2の構成をもとに、カメラ16に替えて図示しない分光計を備えるとともに、バンドパスフィルタ14が取り除かれた構成とされる。この基準分光特性取得部3は、観察画像取得部2の一部を兼用し、観察画像取得部2からバンドパスフィルタ14を取り除くとともに、カメラ16を分光計に置き換えることで適宜実現される。その際、ステージ11上には、標本10に替えて図示しない基準標本が載置される。基準標本には、標本10に含まれる主要要素としての細胞核、細胞質、赤血球等を含むとともにH&E染色された典型的な病理標本が用いられる。なお、基準分光特性取得部3は、観察画像取得部2と個別に構成することもできる。
基準分光特性取得部3は、分光計によって、基準標本上の各標本点における分光特性としての分光透過率(以下、基準分光透過率と呼ぶ。)を取得する。基準分光特性取得部3が取得した標本点ごとの基準分光透過率は、制御部8に転送され、各々対応する標本点位置(標本点座標)と対応付けられて記憶部6に記憶される。
入力部4は、キーボード、マウス等の各種入力デバイスを備え、制御部8に対して観察画像の処理に用いる処理パラメータ等、各種情報の入力を行う。特に、入力部4は、後述する代表点の抽出または最適化に用いるパラメータ値の入力を行う。表示部5は、FPD(Flat Panel Display)等の各種表示器を用いて構成され、観察画像、後述する教師画像、分光特性の推定結果および色素量の推定結果等、各種情報の表示をする。
記憶部6は、ハードディスク、ROMおよびRAM等を用いて構成され、制御部8が実行させる処理プログラム、観察画像取得部2が取得した観察画像、基準分光特性取得部3が取得した基準分光透過率、画像処理部7が用いる処理パラメータ、画像処理部7の処理結果等、各種情報を記憶する。また、記憶部6は、後述の色素量推定に用いる色素H、色素Eおよび色素Rの標準的な分光特性をあらかじめ記憶する。
画像処理部7は、例えばCPUによって実現され、制御部8が実行させる所定の画像処理プログラムに基づき、記憶部6に記憶された観察画像および基準分光透過率に対して種々の演算処理を行う。特に、画像処理部7は、基準分光特性取得部3が取得した基準分光透過率をもとに基準標本を記録した基準標本画像を生成する基準標本画像生成部7aと、基準標本画像の各画素を、少なくとも標本10中の主要要素を特徴付ける主要特徴量に基づいて分類することで主要要素を分類する主要要素分類部7bと、主要要素分類部7bが分類した画素群の中から複数の代表画素を抽出することで代表点を抽出する代表点抽出部7cと、を備える。
また、画像処理部7は、代表点抽出部7cが抽出した複数の代表画素に対応付けられた基準分光透過率をもとにその統計量を算出する統計量算出部7dと、統計量算出部7dが算出した統計量をもとに代表画素の更新要否を判断し、更新が必要であると判定した場合、代表画素を抽出するための抽出条件を更新することで代表点を最適化する代表点最適化部7eと、観察画像取得部2が取得した観察画像の画像データをもとに、統計量算出部7dが算出した統計量を用いて標本10上の各標本点における分光特性を推定する分光特性推定部7fと、分光特性推定部7fが推定した分光特性をもとに各標本点における色素量を推定する色素量推定部7gと、を備える。
制御部8は、CPUによって実現され、記憶部6に記憶された所定の処理プログラムを実行させることで、画像処理装置1が備える各部の処理および動作を制御する。例えば、制御部8は、記憶部6に記憶された所定の画像処理プログラムを実行させることで、記憶部6に記憶された観察画像と基準分光透過率を画像処理部7に処理させるとともに、画像処理部7による処理結果を記憶部6に記憶させ、表示部5に表示させる制御をする。
つづいて、画像処理装置1が行う画像処理手順として、観察画像をもとに標本10の分光特性を推定するとともに、この推定した分光特性をもとに標本10中の色素H、色素Eおよび色素Rの各色素量を推定する処理手順について説明する。図6は、その処理手順の概要を模式的に示す図である。この図に示すように、画像処理装置1は、まず基準標本をもとに基準分光透過率を取得するとともに、取得した基準分光透過率をもとに基準標本画像を生成する。そして、この基準標本画像をもとに基準標本から複数の代表点を抽出するとともに代表点の最適化を行い、ウィナー推定に用いる統計量としての分光透過率の自己相関行列を最適化する。
つぎに、画像処理装置1は、観察画像として取得したマルチバンド画像をもとに、標本10上の各標本点における分光特性をウィナー推定によって推定し、この推定した推定分光特性と、あらかじめ取得した色素H、色素Eおよび色素Rの標準的な分光特性とをもとに、標本10上の各標本点における各色素量(H色素量、E色素量およびR色素量)を推定する。そして、画像処理装置1は、その推定結果として推定色素量と推定分光特性とを表示部5に表示させる。なお、画像処理装置1は、この一連の処理過程の中で、観察画像、基準標本画像、基準分光透過率、推定分光特性および推定色素量等を記憶部6に適宜記憶させる。
図7は、図6に対応する処理手順を示すフローチャートである。この処理手順は、制御部8が記憶部6に記憶された所定の画像処理プログラムを実行させることで行われる。図7に示すように、まず観察画像取得部2は、標本10の観察画像としてのマルチバンド画像を取得する観察画像取得処理を行い(ステップS101)、基準分光特性取得部3は、基準標本の基準分光透過率を取得する基準分光特性取得処理を行う(ステップS102)。
つぎに、基準標本画像生成部7a、主要要素分類部7b、代表点抽出部7c、統計量算出部7dおよび代表点最適化部7eは、基準分光透過率をもとに、ウィナー推定行列に用いる分光透過率の自己相関行列を算出する統計量算出処理を協働して行い(ステップS103)、分光特性推定部7fは、観察画像の画像データをもとに、ステップS103によって算出された自己相関行列を用い、ウィナー推定によって標本10上の各標本点における分光特性を推定する分光特性推定処理を行う(ステップS104)。
その後、色素量推定部7gは、ステップS104によって推定された分光特性と、記憶部6に記憶された色素H、色素Eおよび色素Rの標準分光特性とをもとに、この色素ごとに標本10上の各標本点における各色素量を推定する色素量推定処理を行い(ステップS105)、制御部8は、ステップS104によって推定された推定分光特性およびステップS105によって推定された推定色素量等を表示部5に表示させる推定結果表示処理を行い(ステップS106)、一連の処理を終了させる。
なお、ステップS106では、制御部8は、推定結果に限らず、観察画像、基準標本画像、基準分光透過率等を表示部5に適宜表示させることができる。また、制御部8は、ステップS101〜S106の各ステップにおいて、ステップごとの処理結果および処理パラメータ等を適宜記憶部6に記憶させる。
つづいて、図7に示した各ステップについて詳細に説明する。まずステップS101の観察画像取得処理では、観察画像取得部2は、バンドパスフィルタ14a,14bを順次切り換えながら標本10を撮像することで、例えば図8に示すような観察画像としてのマルチバンド画像を取得する。ただし、図8ではモノクロ表示している。また、観察画像取得部2は、観察画像を取得した後、標本10をステージ11上から取り除いた状態で、背景光のマルチバンド画像として光源12が発した照明光のみのマルチバンド画像を、観察画像の場合と同様の撮像手順によって取得する。制御部8は、このステップS101によって取得された観察画像と背景光のマルチバンド画像とを記憶部6に記憶させる。
つぎに、ステップS104の分光特性推定処理について説明する。この分光特性推定処理では、分光特性推定部7fは、標本10上の各標本点の分光特性をウィナー推定によって推定する。ここで、標本10上の標本点とは、観察画像上の各画素もしくは各画素群に対応する標本10上の被観測点である。
通常、色素は標本10内に3次元的に分布しているが、観察画像では標本10を3次元像として捉えることはできず、標本10内をその厚さ方向に透過した照明光をカメラ16の撮像素子上に結像させた2次元像として捉えている。このため、観察画像では、標本10を平面的に捉えたその平面内での各点を標本点として観察している。この場合、各標本点における色素量とは、その標本点において標本10の厚さ方向に分布した色素の全色素量(積分量)に相当する。
なお、標本10上の標本点と、それに対応する観察画像上の観察点(画素もしくは画素群)との各位置(座標)は、観察画像取得部2における各種の装置パラメータ等に基づく換算式によってあらかじめ対応付けられており、画像処理部7は、その換算式を用いて標本点と観察点との座標を適宜変換することができる。
一般に、観察画像上の観察点位置uおよびバンドbにおける画素値g(u,b)と、その観察点位置uに対応する標本10上の標本点における分光透過率t(u,λ)との間には、次式(1)によって示される関係が成り立つ。次式(1)は、波長λと、観察画像の撮像に用いたバンドパスフィルタのバンドbに対応する分光透過率f(b,λ)と、観察画像の撮像に用いたカメラの分光感度特性s(λ)と、観察画像の撮像に用いた照明光の分光放射特性e(λ)と、バンドbにおける撮像ノイズn(b)とを用いて表されている。なお、バンドbは、観察画像としてのマルチバンド画像における各バンドを識別する通し番号である。
式(1)を波長について離散化すると、画素値g(u,b)、分光透過率t(u,λ)、分光透過率f(b,λ)、分光感度特性s(λ)、分光放射特性e(λ)および撮像ノイズn(b)のそれぞれに対応する行列G(u),T(u),F,S,E,Nを用いて次式(2)によって示される。
G(u)=FSET(u)+N ・・・(2)
G(u)=FSET(u)+N ・・・(2)
ここで、波長のサンプル点数Dとバンド数Bとをもとに、行列G(u)はB行1列、行列T(u)はD行1列、行列FはB行D列、行列SはD行D列、行列EはD行D列、行列NはB行1列とされている。また、行列Sは、分光感度特性s(λ)を対角要素とする対角行列であり、行列Eは、分光放射特性e(λ)を対角要素とする対角行列である。なお、式(2)では、行列を用いて複数のバンドに関する式を集約しているため、変数としてバンドbが陽に記述されていない。また、波長λに関する積分は、行列の積に置き換えられている。
ウィナー推定によると、観察点位置uに対応する標本10上の標本点における分光透過率の推定値を示す行列T^(u)は、式(2)をもとに、行列RSS,RNNを用いて次式(3)によって算出される。ここで、行列RSSは、D行D列の行列であり、推定対象とする標本の分光透過率の自己相関行列である。また、行列RNNは、B行B列の行列であり、観察画像の撮像に用いたカメラのノイズの自己相関行列である。なお、記号()tは転置行列を示し、記号()−1は逆行列を示している。
T^(u)=RSS(FSE)t((FSE)RSS(FSE)t+RNN)−1G(u) ・・・(3)
T^(u)=RSS(FSE)t((FSE)RSS(FSE)t+RNN)−1G(u) ・・・(3)
ステップS104の分光特性推定処理では、分光特性推定部7fは、式(3)を用いて標本10上の各標本点における推定分光透過率を算出する。算出された各推定分光透過率は、記憶部6に記憶される。なお、行列F,S,Eは、観察画像取得部2におけるバンドパスフィルタ14、カメラ16および光源12が選定された後、分光計等を用いてあらかじめ測定され、記憶部6に記憶される。ここでは、結像光学系13の分光透過率は「1.0」と近似されているが、「1.0」からの乖離が許容できない場合に、この分光透過率もあらかじめ測定し、行列Eに乗じておくとよい。
また、行列RSSは、ステップS102によって取得された基準分光透過率をもとに、ステップS103によって後述のようにあらかじめ算出され、記憶部6に記憶される。また、行列RNNは、ステップS101によって撮像された背景光のマルチバンド画像の各バンドについて画素値の自己相関行列を生成することであらかじめ求められ、記憶部6に記憶される。ただし、この場合、バンド間でノイズの相関はないものと仮定される。
つぎに、ステップS105の色素量推定処理について説明する。一般に、光を透過する物質では、波長λごとに入射光の強度I0(λ)と射出光の強度I(λ)との間で、物質に固有の係数k(λ)と物質の厚さdとを用い、次式(4)によって示されるLambert-Beer則が成り立つことが知られている。
H&E染色された標本10では、式(4)をもとに、色素H、色素Eおよび色素Rごとの係数kH(λ),kE(λ),kR(λ)と、観察点位置uに対応する標本点位置Uの厚さdH(U),dE(U),dR(U)とを用いて次式(5)が成り立つ。
式(5)の左辺は、分光透過率であり、係数kH(λ),kE(λ),kR(λ)は、それぞれ色素H、色素Eおよび色素Rに固有の標準分光特性である。係数kH(λ),kE(λ),kR(λ)は、各々対応する単一の色素で染色した標準的な標本を作成し、その分光透過率を分光計によって測定することであらかじめ求められ、記憶部6に記憶される。
また、厚さdH(U),dE(U),dR(U)は、それぞれ色素H、色素Eおよび色素Rに対応した仮想的な厚さであって、標本10が色素H、色素Eおよび色素Rごとに単一の色素で染色されていると仮定した場合と比較して、相対的にどの程度の量の色素が存在しているか、という相対的な色素量を表す指標となる。すなわち、厚さdH(U),dE(U),dR(U)は、それぞれ標本10上の標本点位置Uにおける色素H、色素Eおよび色素Rの色素量を表していると言える。
式(6)の左辺における分光透過率I(λ)/I0(λ)として、式(3)を用いて算出される推定分光特性t^(U,λ)を代入すると、次式(7)を得る。ここで、推定分光特性t^(U,λ)は、式(3)によって算出される行列T^(u)の波長λに対応する要素t^(u,λ)を標本点位置Uについて換算することで求められる。
式(7)において未知変数は厚さdH(U),dE(U),dR(U)の3つであるから、少なくとも3つの異なる波長λに対して式(7)を連立させることで、これを解くことができる。あるいは、より算出精度を高めるために、4つ以上の異なる波長λに対して式(7)を連立させて重回帰分析を行うこともできる。
ステップS105の色素量推定処理では、色素量推定部7gは、式(9)を用いて、標本10上の各標本点位置Uにおける推定色素量として厚さdH(U),dE(U),dR(U)を算出する。算出された各推定色素量は、記憶部6に記憶される。なお、以下の説明では、厚さdH(U),dE(U),dR(U)をそれぞれ推定色素量dH(U),dE(U),dR(U)と呼ぶ。
つぎに、ステップS102の基準分光特性取得処理について説明する。この基準分光特性取得処理では、基準分光特性取得部3は、基準標本上の標本点ごとに、その標本点位置に対応付けて基準分光透過率を取得する。その際、基準分光特性取得部3は、基準標本上の2つの標本領域から基準分光透過率を取得し、制御部8は、その一方の標本領域に対応する基準分光透過率を、後述のように代表点を抽出するための教師データとする。また、他方の標本領域に対応する基準分光透過率を、後述のように代表点の最適化を行うためのテストデータとして、それぞれ記憶部6に記憶させる。なお、テストデータに対応する標本領域を以下、テスト領域と呼ぶ。
つぎに、ステップS103の統計量算出処理について詳細に説明する。図9は、統計量算出処理の処理手順を示すフローチャートである。この図に示すように、まず基準標本画像生成部7aは、ステップS102によって取得された教師データおよびテストデータをもとに、基準標本画像としての教師画像およびテスト画像を生成する(ステップS111)。つづいて主要要素分類部7bは、少なくとも標本10中の主要要素を特徴付ける主要特徴量に基づいて教師画像の各画素を分類し、これによって各主要要素を分類する主要要素分類処理を行う(ステップS112)。
そして、代表点抽出部7cは、ステップS112によって分類された画素群ごとに代表画素を抽出するための抽出条件を、主要要素ごとに代表点を抽出するための抽出条件として設定し(ステップS113)、この設定した抽出条件に基づいて画素群から代表画素を抽出し、これによって代表点を抽出する(ステップS114)。
その後、統計量算出部7dは、ステップS113によって抽出された各代表画素に対応付けられた基準分光透過率をもとに、その統計量としての自己相関行列RSS(式(3)参照)を算出し(ステップS115)、代表点最適化部7eは、ステップS115によって算出された自己相関行列RSSを用いてウィナー推定処理行列を算出する(ステップS116)。
つづいて、代表点最適化部7eは、ステップS111によって生成されたテスト画像の画像データをもとに、ステップS116によって算出したウィナー推定処理行列を用い、基準標本のテスト領域における分光特性としての分光透過率を推定し(ステップS117)、この推定した推定分光透過率の基準分光透過率に対する分光透過率誤差を算出する(ステップS118)。
そして、代表点最適化部7eは、算出した分光透過率誤差が所定の判定条件を満足するか否かを判断し(ステップS119)、判定条件を満足していないと判定した場合(ステップS119:No)、代表画素の抽出条件つまり代表点の抽出条件を更新する(ステップS120)。このステップS120の後、制御部8は、代表点抽出部7c、統計量算出部7dおよび代表点最適化部7eにステップS114からの処理を繰り返させる。一方、分光透過率誤差が所定の判定条件を満足している場合には(ステップS119:Yes)、制御部8は、統計量算出処理を終了し、ステップS103へリターンする。
ステップS111では、基準標本画像生成部7aは、式(1)に基づく次式(10)をもとに、教師データおよびテストデータからそれぞれ教師画像およびテスト画像の画像データを算出する。次式(10)は、基準標本上の各標本点位置Vに対応する教師画像上またはテスト画像上の画素位置vにおける画素値g(v,b)が、教師データまたはテストデータとして取得された基準分光透過率t(v,λ)を用いて算出されることを示している。
ここで、基準標本画像生成部7aは、分光透過率f(b,λ)、分光感度特性s(λ)および分光放射特性e(λ)に、それぞれバンドパスフィルタ14の分光透過率、カメラ16の分光感度特性および光源12の分光放射特性を用いて画素値g(v,b)を算出することで、観察画像と同様の6バンド画像であるテスト画像を生成する。さらに、バンドパスフィルタの分光透過率f(b,λ)を「1」として各画素値g(v,b)を算出することで、RGB画像である教師画像を生成する。生成された教師画像およびテスト画像の画像データは、画素値g(v,b)ごとに、各々対応する基準分光透過率が対応付けられて記憶部6に記憶される。
つぎに、ステップS112の主要要素分類処理について説明する。図10は、その主要要素分類処理の処理手順を示すフローチャートである。この図に示すように、主要要素分類部7bは、まず教師画像の各画素を色空間としてのa*b*空間へ写像する(ステップS131)。このa*b*空間は、L*a*b*表色系におけるa*値およびb*値を各軸にもつ2次元の特徴空間である。主要要素分類部7bは、教師画像上の各画素についてRGB値からa*b*値を算出し、この算出したa*b*値をもとに各画素をa*b*空間へ写像する。これによって、教師画像上の各画素はa*b*空間において、例えば図11に黒点で示すように、各主要要素を特徴付ける主要特徴量として各主要要素を特徴付ける色を示すa*b*値ごとに画素群を形成する。各画素のa*b*値は、教師画像上の画素位置v(xy座標)に対応付けて記憶部6に記憶される。
つづいて、主要要素分類部7bは、a*b*空間におけるヒストグラムを算出するとともに、このヒストグラムを教師画像の総画素数で除算することで、a*b*空間における確率分布を算出する(ステップS132)。このステップS132では、主要要素分類部7bは、例えば図12に示すようなヒストグラムを得る。図12に示すヒストグラムは、図11に示したXII−XII断面におけるヒストグラムに相当する。ここで算出された確率分布は、記憶部6に記憶される。
つづいて、主要要素分類部7bは、教師画像上の画素位置vにおけるa*b*値と、a*b*空間における確率分布とを用いて、a*b*空間における分布関数を画素群ごとに算出する。具体的には、主要要素分類部7bは、EMアルゴリズム(Expectation Maximization Algorithm)を用いて各分布関数を算出する。そのため、主要要素分類部7bは、まず画素群ごとに確率分布に対する分布関数を混合ガウス分布(GMD)によってモデル化する(ステップS133)。
GMDは、一般的な測定値の分布に対して幅広く近似できる関数として知られている。GMDのパラメータθ(以下、GMDパラメータθと呼ぶ。)は、重み係数α、平均μおよび分散σ2を用いて表される。ステップS133では、主要要素分類部7bは、画素群ごとにGMDパラメータθに適当な初期値を設定する。この初期値は、例えば画像処理装置1の操作者によって入力部4から適宜入力される。
つづいて、主要要素分類部7bは、画素群ごとにEMアルゴリズムによって、a*b*空間における分布関数を近似するGMDパラメータθを算出する。EMアルゴリズムは、非特許文献2に開示されているように、GMDを用いた繰り返し演算を用いる最尤推定アルゴリズムの一種として広く知られている。主要要素分類部7bは、そのEMアルゴリズムにおけるEステップ(Expectation Step)(ステップS134)、Mステップ(Maximization Step)(ステップS135)および判定処理(ステップS136)の繰り返し計算によって画素群ごとにGMDパラメータθの最尤解を推定し、この最尤解としてのGMDパラメータθを用いて画素群ごとに最適な分布関数を算出する。
具体的には、ステップS135では、主要要素分類部7bは、ステップS131によって算出した各画素のa*b*値をGMDにおける測定値とし、ステップS132によって算出したヒストグラムの総数(総画素数)NをGMDにおける測定値の個数として、GMDパラメータθを構成する重み係数α、平均μおよび分散σ2をそれぞれ次式(11)〜(13)によって順次更新する。ここで、次式(11)〜(13)は、EMアルゴリズムの繰り返し数t、画素群の識別番号jおよび注目画素の識別番号zを用いて表されている。
さらに、主要要素分類部7bは、正規分布関数pj、条件付き期待値pおよび評価値l(θj)をそれぞれ次式(14)〜(16)によって順次更新する。ここで、次式(14)〜(16)は、自己相関行列Σj、自己相関行列の行列式|Σj|、次元数dを用いて表されている。本実施の形態では、a*b*空間であるため次元数d=2となる。なお、符号Tは、転置行列であることを示している。
そして、ステップS136の判定処理では、主要要素分類部7bは、評価値l(θj)が所定の収束条件を満足しているか否かを判断し、収束条件が満足されていない場合(ステップS136:No)、ステップS134からの処理を繰り返し、収束条件が満足されている場合には(ステップS136:Yes)、EMアルゴリズムを終了させるとともに主要要素分類処理を終了させ、ステップS112へリターンする。
ここで、所定の収束条件とは、例えば評価値l(θj)の収束値が所定の閾値以内であること、もしくは評価値l(θj)が所定の閾値以内であることとされる。また、評価値l(θj)の収束値は、対数尤度として算出される評価値l(θj)について、現評価時点までの一連のEMアルゴリズムの繰り返し演算の中で最も高い評価を得た対数尤度から現評価時点の対数尤度を差し引いた差分として算出される。
以上の主要要素分類処理によって、主要要素分類部7bは、最適なGMDパラメータθと、a*b*空間における確率分布を最も近似する分布関数p(a*,b*:θ)とを求めることができる。分布関数p(a*,b*:θ)のうちの画素群ごとの分布関数pj(a*,b*:θj)は、式(15)における右辺の分子に相当するものであって、次式(17)によって示される。
画素群ごとの分布関数pj(a*,b*:θj)は、例えば図12に示すように求められる。図12に示す分布関数pj(a*,b*:θj)は、図11に示したXII−XII断面における分布関数pj(a*,b*:θj)に相当する。このようにして求められた画素群ごとの分布関数pj(a*,b*:θj)によって、教師画像上の各画素が分類されるとともに各主要要素が分類される。ここで算出された分布関数pj(a*,b*:θj)は、教師画像上の画素位置vにおけるa*b*値と、a*b*空間における確率分布とに対応付けて記憶部6に記憶される。
つぎに、ステップS113による代表点の抽出条件の設定処理について説明する。このステップS113では、代表点抽出部7cは、ステップS112によって分類された画素群ごとに代表画素を抽出するための抽出条件として、画素群ごとの代表画素の抽出範囲と、その抽出範囲内の抽出個数とを設定する。具体的には、代表点抽出部7cは、ステップS112によって算出したGMDパラメータと分布関数p(a*,b*:θ)とを用い、代表点の抽出範囲として分類特徴量範囲C(a*,b*)を算出する。そして、ステップS112によって分類された画素群ごと、つまり主要要素ごとの分類特徴量範囲Cj(a*,b*)を、画素群ごとの分布関数pj(a*,b*:θj)における平均μを中心位置、この中心位置からの広がりを分散σ2とする等高線の内側の範囲とする。これによって、画素群ごとの分類特徴量範囲Cj(a*,b*)は、例えば図11に示すように求められる。ここで算出された分類特徴量範囲C(a*,b*)は、分布関数p(a*,b*:θ)に対応付けて記憶部6に記憶される。
また、代表点抽出部7cは、分類特徴量範囲C(a*,b*)から代表画素を抽出する抽出方法として、例えばランダム抽出を設定する。ランダム抽出では、分類特徴量範囲C(a*,b*)におけるa*値およびb*値の複数の組み合わせがランダムに決定され、決定されたa*値およびb*値の複数の組み合わせによって示される各画素が代表画素として抽出される。また、a*値およびb*値の1つの組み合わせに対して複数の画素が対応付けられている場合には、その複数の画素の中からもランダムに抽出が行われる。代表画素の抽出は、各画素群の分類特徴量範囲Cj(a*,b*)において、その範囲ごとに設定された抽出個数mjまで繰り返し行われる。ただし、一度代表画素として抽出された画素は、それ以降、その一連の抽出処理における選択対象から除外される。
本実施の形態では、標本が病理標本であることを鑑みて、代表点抽出部7cは、抽出範囲内の抽出個数つまり主要要素ごとの抽出個数mjを、代表画素の総抽出個数Mを主要要素の種類数qで除算した数M/qとする。また、代表点抽出部7cは、複数の主要要素における特定の主要要素の推定精度を向上させる場合、重み係数αを用いて抽出個数の重み付けをすることができる。具体的には、画素群ごとつまり主要要素ごとの重み係数αjと総抽出個数Mとの積αjMによって抽出個数mjを定め、推定精度を向上させたい主要要素の重み係数αjを他の主要要素の重み係数αjに比して大きい値とする。ただし、重み係数αjは、次式(18)を満足する。
また、代表点抽出部7cは、分類特徴量範囲C(a*,b*)から代表画素を抽出する抽出方法として、ランダム抽出以外にも確率依存型抽出あるいは測定値依存型抽出を用いることができる。ここで、確率依存型抽出とは、ステップS112によって算出された分布関数p(a*,b*:θ)をもとに重み付けを行って代表画素を抽出する方法である。つまり、分布関数p(a*,b*:θ)の確率をもとに、分類特徴量範囲C(a*,b*)における代表画素のa*値およびb*値を決定する。確率依存型抽出は、主要要素内の特徴量の分布に依存した代表点の抽出を行いたい場合に用いることができ、例えば、病理標本中の各主要要素に対応する教師画像上の平均的なRGB値を持つ画素を代表画素として優先的に抽出したい場合に用いることができる。
一方、測定値依存型抽出は、主要要素内の特徴量に依存した代表点の抽出を行いたい場合に用いることができ、例えば、病理標本中の主要要素に対応する教師画像上の画素のうち、次式(19)によって示されるa*b*空間上の原点からのユークリッド距離Lyが大きい画素を代表画素として優先的に抽出したい場合、もしくはその逆の場合に用いることができる。
このとき、代表画素の抽出測定強度範囲を限定し、測定値強度範囲を抽出条件に付加する。測定値強度範囲は、分類特徴量範囲Cj(a*,b*)における最大のユークリッド距離Lyを100%として分類特徴量範囲Cj(a*,b*)の対強度比率で表せばよい。例えば、主要要素の分類特徴量範囲Cj(a*,b*)のうち、a*b*空間上の原点からのユークリッド距離Lyが最小の場合は0%となり、測定値強度を0%から30%の範囲で代表画素を抽出したい場合には、測定値強度範囲を0%〜30%の範囲に指定すればよい。ここで、強度とはa*b*空間におけるa*値およびb*値に相当し、強度範囲とはa*値およびb*値の範囲に相当する。
なお、代表画素の抽出範囲および抽出個数は、以降の代表点の最適化処理によって最適解が導出されるため、このステップS113によって設定する初期値に制限はない。しかしながら、本実施の形態では、標本が病理標本であることを鑑みて、各主要要素に対する抽出個数mjは、総抽出個数Mを主要要素の種類数q=4で除算した数M/4とすることで最適解の収束性が向上する。同様の目的から分類特徴量範囲Cj(a*,b*)は、各分布関数pj(a*,b*:θj)の確率分布に対して相当の割合、例えば90%となる範囲を算出する平均μおよび分散σ2を初期値にするとよい。
つぎに、ステップS116〜S120による代表点の最適化処理について説明する。まずステップS116では、代表点最適化部7eは、システム行列Hを次式(20)によって算出し、この算出したシステム行列Hをもとにウィナー推定処理行列としての推定オペレータWを次式(21)によって算出する。なお、式(20)および式(21)に用いる行列F,S,E,RSS,RNNは、式(3)に用いた各行列に対応する。
H=(FSE)t ・・・(20)
W=RSSH(HRSSHt+RNN)−1 ・・・(21)
H=(FSE)t ・・・(20)
W=RSSH(HRSSHt+RNN)−1 ・・・(21)
ここで、式(3)は、テスト画像上の画素位置vに対し、式(21)を用いて次式(22)のように書き直すことができる。ステップS117では、代表点最適化部7eは、テスト画像における複数の画素に対してそれぞれ次式(22)を用いて行列T^(v)を算出し、これによって基準標本のテスト領域における分光透過率を推定する。
T^(v)=WG(v) ・・・(22)
T^(v)=WG(v) ・・・(22)
その後、ステップS118では、代表点最適化部7eは、ステップS117によって算出した行列T^(v)と、ステップS102によってテストデータとして取得した基準分光透過率t^(v)とをもとに、次式(23)によって分光透過率誤差ΔEを算出する。分光透過率誤差ΔEは、式(23)に示すように、個数Nの画素における推定分光透過率と基準分光透過率との差分の平均値として求められる。
ステップS119では、代表点最適化部7eは、ステップS118によって算出した分光透過率誤差ΔEを自己相関行列の評価値として用い、分光透過率誤差ΔEが所定の判定条件を満足するか否かによって、ステップS114で抽出された代表画素つまり代表点の更新要否を判定する。
ここで、所定の判定条件とは、分光透過率誤差ΔEが所定の閾値内であること、もしくはステップS114〜S120の処理を繰り返すことによって分光透過率誤差ΔEが収束していることとされる。分光透過率誤差ΔEの収束条件としては、例えばステップS114〜S120を複数回繰り返した際の分光透過率誤差ΔEの最大値と最小値との差分が所定の閾値内であることとされる。
ステップS120では、代表点最適化部7eは、例えば代表画素の抽出範囲と、抽出範囲内の抽出個数との少なくとも一方を変更することで抽出条件の更新を行う。抽出範囲の更新では、代表点最適化部7eは、分類特徴量範囲C(a*,b*)を定める分散σ2を変更することで分類特徴量範囲C(a*,b*)を更新する。また、代表点最適化部7eは、抽出範囲および抽出個数をそれぞれ変化させた方向を更新方向として記憶部6に随時記憶させ、この記憶させた更新方向を次回の更新に反映させる。例えば、一回目の更新方向に対して分光透過率誤差ΔEの値が減小した場合、つぎの更新時には抽出範囲および抽出個数を同じ更新方向へ変化させる。更新方向は、具体的には、抽出確率幅としての分散σ2の増減方向と、総抽出個数Mもしくは主要要素ごとの抽出個数mjの増減方向とに相当する。
なお、制御部8は、ステップS120において代表点最適化部7eに抽出条件を更新させずに、ステップS114からの処理を繰り返させることもできる。この場合、抽出条件をそのままにして抽出範囲内から代表画素をランダムに抽出し直すこととなる。これによって、抽出条件を更新することなく代表画素つまり代表点の組み合わせを変更することができる。また、代表点最適化部7eは、ステップS120によって代表画素の抽出方法をランダム抽出から確率依存型抽出方法もしくは測定値依存型抽出へ変更し、これによって代表点の組み合わせを変更することもできる。
以上説明したように、本実施の形態にかかる画像処理装置1は、1以上の主要要素を含む標本10を撮像した観察画像の画像データをもとに、この1以上の主要要素を含む基準標本上の複数の標本点における基準分光特性の統計量を用いて標本10の分光特性を推定する画像処理装置であり、基準標本が記録された画像であって画素ごとに基準標本上の対応する標本点における基準分光特性が対応付けられた教師画像の各画素を、少なくとも主要要素を特徴付ける主要特徴量に基づいて分類する主要要素分類部7bと、主要要素分類部7bが分類した画素群の中から複数の代表画素を抽出する代表点抽出部7cと、代表点抽出部7cが抽出した複数の代表画素に対応付けられた基準分光特性をもとに、標本10の分光特性の推定に用いる基準分光特性の統計量を算出する統計量算出部7dと、を備えている。
このため、画像処理装置1では、基準標本中の各主要要素から均等もしくは所望の割合で複数の代表点を抽出することができ、この複数の代表点における基準分光特性をもとに、標本10の分光特性の推定に用いる基準分光特性の統計量を算出することができるとともに分光特性の推定に用いる処理行列としてのウィナー推定処理行列を算出することができる。これによって、画像処理装置1では、標本10上の各標本点における分光特性を、意図しないばらつきを生じさせることなく高精度に推定することができる。
また、画像処理装置1は、統計量算出部7dが算出した統計量をもとに代表画素の更新要否を判断し、更新が必要であると判定した場合、代表画素を抽出するための抽出条件を更新する代表点最適化部7eを備え、代表点抽出部7cは、代表点最適化部7eが更新した抽出条件に基づいて複数の代表画素を再抽出し、統計量算出部7dは、代表点抽出部7cが再抽出した複数の代表画素に対応付けられた基準分光特性をもとに統計量を算出する。このため、画像処理装置1では、1以上の主要要素を含む標本10の分光特性の推定に対して最適な代表点を基準標本から抽出することができるとともに、この抽出した代表点における基準分光特性をもとに最適な統計量およびウィナー推定処理行列を算出することができる。これによって、画像処理装置1では、従来技術に比して標本10の分光特性の推定精度を向上させることができる。
なお、画像処理装置1では、以上のようにして推定した標本10の推定分光特性をもとに、標本10に含まれた複数の色素の色素量を推定しているため、高精度に各色素量を推定することもできる。
ここまで、本発明を実施する最良の形態を実施の形態として説明したが、本発明は、上述した実施の形態に限定されず、本発明の趣旨を逸脱しない範囲であれば、種々の変形が可能である。
例えば、上述した実施の形態では、主要要素分類部7bは、標本10中の主要要素を特徴付ける主要特徴量として主要要素を特徴付ける色を示すa*b*値をもとに基準標本画像上の各画素を分類するものとしたが、a*b*値に限定されず、RGB値もしくはXYZ値等、他の表色系によって示される値を用いて基準標本画像上の各画素を分類することもできる。ここで、XYZ値は、XYZ表色系における三刺激値を意味する。さらに、主要特徴量として色情報に限定されず、濃淡情報、形状情報、テクスチャ情報等を用いることができる。
ただし、形状情報やテクスチャ情報は、単一の画素では特徴量を示すことができないため、基準標本画像の領域分割、例えばブロック分割等の処理が必要となる。また、病理標本を推定対象とする場合、テクスチャを用いた特徴量としてエントロピーを用いてもよい。ここで、主要要素分類処理が各主要要素に対応する画素の分類であることを鑑みると、単一の画素を単位にして分類することにこだわる必要はなく、複数の画素からなる画素群を単位として分類を行っても構わない。
なお、a*b*値以外の特徴量を用いる場合にも、上述のように教師画像上の各画素をa*b*空間に写像して分類を行った場合と同様に、各特徴量を示す特徴空間へ教師画像上の各画素を写像して分類を行うことで、それ以降の処理を、上述した実施の形態と同様の処理によって行うことができる。
また、画像処理装置1では、上述したステップS112の主要要素分類処理におけるEMアルゴリズムの近似精度を向上させるために、さらにつぎのような処理を行うことができる。すなわち、まず基準分光特性取得部3は、ステップS102の基準分光特性取得処理によって、背景領域の分光透過率を背景データとして取得する。そして、基準標本画像生成部7aは、ステップS111によって、背景データに対応する背景画像を、教師画像と同様に生成する。その後、主要要素分類部7bは、ステップS112の主要要素分類処理において、背景画像の各画素をa*b*空間に写像し、この写像した背景画像のa*b*空間上で2次元正規分布の最尤推定法から最も近似する単一正規分布関数のパラメータを算出する。この処理は、上述したEMアルゴリズムを単一正規分布関数として用いることで処理を行うことができる。つぎに、主要要素分類部7bは、算出した単一正規分布関数を用いて背景領域に相当するa*b*値範囲を決定し、教師画像を用いて写像したa*b*空間から背景領域に相当する範囲を除去する。教師データの色空間から事前に背景領域を除去することで、EMアルゴリズムの近似対象を細胞核、細胞質および赤血球の3領域に限定することができる。一般的な病理標本では、細胞質と背景の境界は色空間上で接していることが多く、事前に背景をEMアルゴリズムの処理対象外とすることで、より高精度にその他の領域の分布関数の近似式を得ることができる。
また、上述した実施の形態では、画像処理装置1は、推定した分光透過率をもとに色素量を推定するものとして説明したが、色素量の推定に限定されず、推定した分光透過率および/または推定した色素量をもとに、病理診断支援につながる種々の処理を行うようにすることができる。
また、上述した実施の形態では、観察画像取得部2は、観察画像を撮像して取得するものとして説明したが、あらかじめ外部装置によって撮像された観察画像を取得する構成とすることもできる。この場合、観察画像取得部2は、外部装置等から観察画像の画像データを入力可能なデータ通信インターフェースを備えるとよい。
なお、上述した実施の形態では、画像処理装置1は、標本10および基準標本の分光特性として分光透過率を扱うものとして説明したが、分光透過率に限定されず、分光反射率を扱うこともできる。その場合、図2に示した観察画像取得部2および基準分光特性取得部3は、光源12の代わりに、標本10および基準標本に対して落射照明を行う照明機構を備えるとよい。
1 画像処理装置
2 観察画像取得部
3 基準分光特性取得部
4 入力部
5 表示部
6 記憶部
7 画像処理部
7a 基準標本画像生成部
7b 主要要素分類部
7c 代表点抽出部
7d 統計量算出部
7e 代表点最適化部
7f 分光特性推定部
7g 色素量推定部
8 制御部
10 標本
11 ステージ
12 光源
13 結像光学系
14,14a,14b バンドパスフィルタ
15 ターレット
16 カメラ
2 観察画像取得部
3 基準分光特性取得部
4 入力部
5 表示部
6 記憶部
7 画像処理部
7a 基準標本画像生成部
7b 主要要素分類部
7c 代表点抽出部
7d 統計量算出部
7e 代表点最適化部
7f 分光特性推定部
7g 色素量推定部
8 制御部
10 標本
11 ステージ
12 光源
13 結像光学系
14,14a,14b バンドパスフィルタ
15 ターレット
16 カメラ
Claims (13)
- 1以上の主要要素を含む標本を撮像した観察画像の画像データをもとに、前記1以上の主要要素を含む基準標本上の複数の標本点における基準分光特性の統計量を用いて前記標本の分光特性を推定する画像処理装置において、
前記基準標本が記録された画像であって画素ごとに前記基準標本上の対応する標本点における前記基準分光特性が対応付けられた教師画像の各画素を、少なくとも前記主要要素を特徴付ける主要特徴量に基づいて分類する画素分類手段と、
前記画素分類手段が分類した画素群の中から複数の代表画素を抽出する代表抽出手段と、
複数の前記代表画素に対応付けられた前記基準分光特性をもとに前記統計量を算出する統計量算出手段と、
を備えたことを特徴とする画像処理装置。 - 前記統計量算出手段が算出した前記統計量をもとに前記代表画素の更新要否を判断し、更新が必要であると判定した場合、前記代表画素を抽出するための抽出条件を更新する代表最適化手段を備え、
前記代表抽出手段は、前記代表最適化手段が更新した前記抽出条件に基づいて複数の前記代表画素を再抽出し、
前記統計量算出手段は、前記代表抽出手段が再抽出した複数の前記代表画素に対応付けられた前記基準分光特性をもとに前記統計量を算出することを特徴とする請求項1に記載の画像処理装置。 - 前記代表最適化手段は、前記基準標本の領域のうち前記教師画像に記録されていないテスト領域が記録された画像であって画素ごとに該テスト領域上の対応する標本点における前記基準分光特性が対応付けられたテスト画像の画像データをもとに、前記統計量算出手段が算出した前記統計量を用いて前記テスト領域の分光特性を推定するとともに、この推定した分光特性と前記テスト領域の前記基準分光特性とを比較し、この比較結果に基づいて前記代表画素の更新要否を判断することを特徴とする請求項2に記載の画像処理装置。
- 前記代表抽出手段は、前記画素分類手段の分類結果をもとに前記画素群における前記代表画素の抽出範囲と該抽出範囲内の抽出個数との少なくとも一方を前記代表画素を抽出するための抽出条件として設定することを特徴とする請求項1〜3のいずれか一つに記載の画像処理装置。
- 前記代表抽出手段は、前記分類結果をもとに前記画素群における前記代表画素の抽出範囲と該抽出範囲内の抽出個数との少なくとも一方を前記代表画素を抽出するための抽出条件として設定し、
前記代表最適化手段は、前記代表画素の更新が必要であると判定した場合、前記分類結果をもとに前記抽出範囲と該抽出範囲内の抽出個数との少なくとも一方を更新することを特徴とする請求項2または3に記載の画像処理装置。 - 前記画素分類手段は、前記主要特徴量をもとに前記教師画像の各画素を特徴空間に写像し、該特徴空間に写像された各画素の分布に基づいて前記教師画像の各画素を分類することを特徴とする請求項1〜5のいずれか一つに記載の画像処理装置。
- 前記主要特徴量は、色情報、濃淡情報およびテクスチャ情報の少なくとも1つであることを特徴とする請求項6に記載の画像処理装置。
- 前記代表抽出手段は、前記特徴空間に写像された各画素の分布に基づいて、前記画素群における前記代表画素の抽出範囲と該抽出範囲内の抽出個数との少なくとも一方を、前記代表画素を抽出するための抽出条件として設定することを特徴とする請求項6または7に記載の画像処理装置。
- 前記代表抽出手段は、前記特徴空間に写像された各画素の分布と前記主要特徴量とに基づいて、前記画素群における前記代表画素の抽出範囲と該抽出範囲内の抽出個数との少なくとも一方を、前記代表画素を抽出するための抽出条件として設定することを特徴とする請求項6または7に記載の画像処理装置。
- 前記代表抽出手段は、前記代表画素を抽出するための抽出条件に基づいて前記画素群の中から複数の前記代表画素をランダムに抽出することを特徴とする請求項1〜9のいずれか一つに記載の画像処理装置。
- 前記画素分類手段は、前記特徴空間に写像された各画素の分布確率を算出し、
前記代表抽出手段は、前記代表画素を抽出するための抽出条件と前記分布確率とに基づいて、前記画素群の中から複数の前記代表画素を抽出することを特徴とする請求項6〜9のいずれか一つに記載の画像処理装置。 - 前記分光特性および前記基準分光特性は、分光透過率特性または分光反射率特性であることを特徴とする請求項1〜11のいずれか一つに記載の画像処理装置。
- 1以上の主要要素を含む標本を撮像した観察画像の画像データをもとに、前記1以上の主要要素を含む基準標本上の複数の標本点における基準分光特性の統計量を用いて前記標本の分光特性を推定する画像処理装置に、
前記基準標本が記録された画像であって画素ごとに前記基準標本上の対応する標本点における前記基準分光特性が対応付けられた教師画像の各画素を、少なくとも前記主要要素を特徴付ける主要特徴量に基づいて分類する画素分類手順と、
前記画素分類手順によって分類された画素群の中から複数の代表画素を抽出する代表抽出手順と、
複数の前記代表画素に対応付けられた前記基準分光特性をもとに前記統計量を算出する統計量算出手順と、
を実行させることを特徴とする画像処理プログラム。
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| JP2007173245A JP2009014354A (ja) | 2007-06-29 | 2007-06-29 | 画像処理装置および画像処理プログラム |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| JP2007173245A JP2009014354A (ja) | 2007-06-29 | 2007-06-29 | 画像処理装置および画像処理プログラム |
Publications (1)
| Publication Number | Publication Date |
|---|---|
| JP2009014354A true JP2009014354A (ja) | 2009-01-22 |
Family
ID=40355455
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| JP2007173245A Withdrawn JP2009014354A (ja) | 2007-06-29 | 2007-06-29 | 画像処理装置および画像処理プログラム |
Country Status (1)
| Country | Link |
|---|---|
| JP (1) | JP2009014354A (ja) |
Cited By (3)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| JP2011099823A (ja) * | 2009-11-09 | 2011-05-19 | Olympus Corp | バーチャル顕微鏡システム |
| JP2012078156A (ja) * | 2010-09-30 | 2012-04-19 | Olympus Corp | 検査装置 |
| JP2014055793A (ja) * | 2012-09-11 | 2014-03-27 | Tokai Kogaku Kk | クラス判別方法 |
-
2007
- 2007-06-29 JP JP2007173245A patent/JP2009014354A/ja not_active Withdrawn
Cited By (5)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| JP2011099823A (ja) * | 2009-11-09 | 2011-05-19 | Olympus Corp | バーチャル顕微鏡システム |
| US8780191B2 (en) | 2009-11-09 | 2014-07-15 | Olympus Corporation | Virtual microscope system |
| JP2012078156A (ja) * | 2010-09-30 | 2012-04-19 | Olympus Corp | 検査装置 |
| US8836779B2 (en) | 2010-09-30 | 2014-09-16 | Olympus Corporation | Inspection device |
| JP2014055793A (ja) * | 2012-09-11 | 2014-03-27 | Tokai Kogaku Kk | クラス判別方法 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| EP1470411B1 (en) | Method for quantitative video-microscopy and associated system and computer software program product | |
| AU2014230469B2 (en) | Spectral unmixing | |
| EP3053138B1 (en) | Systems and methods for adaptive histopathology image unmixing | |
| JP5508792B2 (ja) | 画像処理方法、画像処理装置および画像処理プログラム | |
| US8977017B2 (en) | System and method for support of medical diagnosis | |
| JP2010169467A (ja) | 画像処理装置、データセット生成装置、画像処理プログラムおよびデータセット生成プログラム | |
| JP5154844B2 (ja) | 画像処理装置および画像処理プログラム | |
| JP5305618B2 (ja) | 画像処理装置および画像処理プログラム | |
| JP5117274B2 (ja) | 画像処理装置および画像処理プログラム | |
| US9406118B2 (en) | Stain image color correcting apparatus, method, and system | |
| JP5178226B2 (ja) | 画像処理装置および画像処理プログラム | |
| JP4740068B2 (ja) | 画像処理装置、画像処理方法、および画像処理プログラム | |
| JP5137481B2 (ja) | 画像処理装置、画像処理プログラム、および画像処理方法 | |
| JP4920507B2 (ja) | 画像処理装置および画像処理プログラム | |
| CN113570538A (zh) | 一种叶片rgb图像偏态分布参数信息采集及分析方法 | |
| CN113610965A (zh) | 基于点云的植物叶绿素含量三维空间立体分布可视化方法 | |
| JP2009014354A (ja) | 画像処理装置および画像処理プログラム | |
| JP5210571B2 (ja) | 画像処理装置、画像処理プログラムおよび画像処理方法 | |
| JP2008304205A (ja) | 分光特性推定装置および分光特性推定プログラム | |
| US11037294B2 (en) | Image processing device, image processing method, and computer-readable recording medium | |
| JP2018169390A (ja) | 染色体異常判定装置 | |
| JP2010169592A (ja) | 画像処理装置および画像処理プログラム | |
| US9122904B2 (en) | Method for optimization of quantitative video-microscopy and associated system | |
| WO2020037255A1 (en) | Automatic identification and analysis of a tissue sample | |
| JP2009025147A (ja) | 画像処理装置および画像処理プログラム |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| A300 | Application deemed to be withdrawn because no request for examination was validly filed |
Free format text: JAPANESE INTERMEDIATE CODE: A300 Effective date: 20100907 |