JP2020004033A - 画像処理方法、コンピュータプログラムおよび記録媒体 - Google Patents

Abstract

【課題】細胞を含む画像を、画像の特徴に応じて複数種の領域に適切に分割する。【解決手段】本発明に係る画像処理方法は、培養された細胞が含まれる原画像を取得する工程（ステップＳ１０１）と、原画像を所定数の画素からなるブロックに分割し、各ブロック内の画像が有する空間周波数成分を求める工程（ステップＳ１０４〜Ｓ１０６）と、ブロックの各々を、空間周波数成分のうちの直流成分の強度と、空間周波数が互いに異なる複数の交流成分それぞれの強度とを特徴量とする多次元特徴量空間における複数のクラスタに分類する工程（ステップＳ１０７、Ｓ１０８）と、原画像を、分類の結果に基づき複数の領域に分割する工程（ステップＳ１０９）とを備える。【選択図】図１

Description

この発明は、培養環境下で撮像された細胞とその背景とを含む画像から、細胞の領域とそれ以外の領域とを分割する画像処理技術に関するものである。

培養環境下で細胞を含んで撮像された画像においては、本来の撮像対象物である細胞の領域とその背景に当たる領域とが画像内で区別して扱われる場合がある。例えば細胞の大きさを求めるために、細胞と背景との境界を画定させる必要がある。従来、例えば細胞を所定の薬剤で染色したり蛍光試薬を投入したりすることで細胞の領域を特定することが行われてきた。しかしながら、例えば細胞の変化を経時的に観察する等の目的で、撮像対象である細胞を改変するような処理を用いずに細胞領域を特定することが求められる場合がある。

培地中の細胞は透明に近く屈折率も水とあまり違わないため、撮像された画像においては細胞と背景との境界が不明瞭であることも多い。この問題に対応するため、本願出願人はこれまでに、特に細胞の輪郭が明瞭に現れた画像を得ることのできる明視野撮像技術を提案している（例えば、特許文献１参照）。画像における細胞と背景との区別には、例えば細胞と背景との輝度の差を利用した２値化処理、細胞輪郭を抽出するためのエッジ検出などの画像処理が利用可能である。

特開２０１６−０７１５８８号公報

培養環境下の細胞を経時的に観察する場合、培養の初期段階では撮像された画像に含まれるのは細胞と概ね一様な背景たる培地とである。このため細胞と背景との識別は比較的容易である。しかしながら、培養が進むにつれて、細胞から排出される老廃物や死んだ細胞の残骸など（デブリ）が撮像対象である細胞の周りに現れ、画像に含まれるようになってくる。このため、２値化や輪郭抽出などの単純な処理では細胞とデブリとを適切に識別することが困難である。

また、例えばＥＳ細胞、ｉＰＳ細胞等の幹細胞を用いた培養においては、幹細胞がフィーダー細胞と共培養されることがある。また、培養に伴い幹細胞が体細胞に分化する過程で、幹細胞とは別の形態的特徴を有する細胞が現れてくる。このため、撮像される画像には複数種類の細胞が混在することとなり、画像内でそれらを識別する必要が生じる場合がある。

これらのことから、例えば画像に現れるテクスチャ等の画像特徴の違いによって、画像を複数種、特に３種以上の領域に分割することのできる技術が必要になってきている。

この発明は上記課題に鑑みなされたものであり、細胞を含む画像を、画像の特徴に応じて複数種の領域に適切に分割することのできる技術を提供することを目的とする。

この発明に係る画像処理方法の一の態様は、上記目的を達成するため、培養された細胞が含まれる原画像を取得する工程と、前記原画像を所定数の画素からなるブロックに分割し、各ブロック内の画像が有する空間周波数成分を求める工程と、前記ブロックの各々を、前記空間周波数成分のうちの直流成分の強度と、空間周波数が互いに異なる複数の交流成分それぞれの強度とを特徴量とする多次元特徴量空間における複数のクラスタに分類する工程と、前記原画像を、前記分類の結果に基づき複数の領域に分割する工程とを備えている。

このように構成された発明では、画像内における局所的な空間周波数成分の分布に基づいて画像が分割される。すなわち、画像がブロックに分割され、ブロックごとの空間周波数成分が求められる。空間周波数成分のうち直流成分は当該ブロックにおける平均的な輝度を表す。一方、交流成分はブロック内での輝度の変化を表し、比較的低い周波数成分は緩やかな輝度変化を、高い周波成分はより細かな輝度変化の状態をそれぞれ表す。したがって、このようにブロックを空間周波成分ごとの強度で表現した場合、直流成分は当該ブロックの明るさを、また交流成分はブロック内のテクスチャをそれぞれ反映したものとなる。特に周波数が高いほど細かい構造のテクスチャを表すものとなる。

これらの各周波数成分の強度を当該ブロックの特徴量と考えれば、これらの特徴量によってブロック内の画像内容の明るさやテクスチャの複雑さを定量的に表現することが可能である。画像中に例えば複数種の細胞が含まれている場合でも、その明るさやテクスチャの違いが特徴量の違いとなって現れる。

このことを利用して、原画像を画像内容に応じて複数の領域に分割することが可能である。具体的には、各ブロックが有する直流成分および複数の交流成分の強度を特徴量とする多次元特徴量空間において各ブロックが占める位置を考えると、細胞に該当するブロックとそれ以外のブロックとの間、あるいは互いに異なる種類の細胞に対応するブロックの間でその分布に有意な差が現れる。したがって、原画像の一部である各ブロックをその特徴量に基づき分類すれば、同じようなテクスチャを有するブロックが特徴量空間で１つのクラスタを形成する。特徴量空間で近接する複数のブロックは互いに同じ種類の細胞を含むと考えられる。一方、２つのブロックが特徴量空間で遠い位置にある場合には、それらブロックは異なる種類の細胞を含むと考えられる。このように、特徴量に基づく分類結果を利用して、原画像をその画像内容に応じて複数の領域に分割することが可能である。

本発明に係る画像処理方法は、コンピュータをその実行主体とすることが可能なものである。この意味において、本発明は、コンピュータに上記処理を実行させるためのコンピュータプログラムとして実現することが可能である。また、当該コンピュータプログラムを記録した記録媒体として実現することも可能である。

上記のように、本発明の画像処理方法によれば、原画像を分割したブロックの各々について、当該ブロック内の画像が有する空間周波数成分の特徴に基づき、原画像が複数の領域に分離される。直流成分の強度と、周波数の異なる複数の交流成分それぞれの強度とを特徴量とすることで、ブロック全体の明るさやテクスチャの違いに応じて適切な分類が可能である。

本実施形態の領域分割処理を示すフローチャートである。 この実施形態の領域分割処理を実行するコンピュータ装置の構成例である。 原画像の例を示す図である。 ＤＣＴ演算の概念を示す図である。 ＤＣＴ係数のグループ化の例を示す図である。 判定基準を策定するための処理を示すフローチャートである。 処理過程における指示入力態様を例示する図である。 特徴量空間におけるクラスタ分布の例を模式的に示す図である。

以下、本発明に係る画像処理方法の一実施形態について説明する。この実施形態における画像処理は、培地中で培養された細胞または細胞の集合体を撮像した画像から、細胞の領域とその他の領域とを分割し、また複数種の細胞が含まれる場合にはそれらに対応する領域を相互に分割する領域分割処理である。この画像処理は、細胞を含む試料を撮像する撮像機能を備えた撮像装置によって撮像された未加工の原画像に対してなされる。撮像装置が撮像後の画像データに対する後処理の１つとして本画像処理を実行してもよく、また撮像装置あるいは適宜の記憶手段から画像データを受け取ったコンピュータ装置が、本画像処理を実行してもよい。

ここでは、既に実行された撮像により生成された原画像データに対し、汎用のコンピュータ装置が本実施形態の画像処理を実行する態様を例示して説明する。撮像装置の構成については特に限定されず、培養された細胞を含む試料を培地とともに撮像し、デジタル画像データとして出力する機能を有するものであればよい。また画像は明視野撮像画像であることが好ましい。また、コンピュータ装置についても、例えばパーソナルコンピュータとして製品化されている一般的なハードウェア構成を有するものを利用可能である。以下ではこれらのハードウェアについての詳しい説明を省略する。

図１は本実施形態の領域分割処理を示すフローチャートである。まず、図１を参照してこの実施形態における領域分割処理の概略について説明し、その後で処理の具体的内容について詳しく説明する。コンピュータ装置は、培地中で培養された細胞を含む試料を適宜の撮像装置および撮像条件で明視野撮像することにより得られる原画像を取得する（ステップＳ１０１）。原画像は少なくとも１つの細胞の全体を含むものとする。したがって、原画像は細胞と細胞以外のデブリや背景等とを含み得る。

細胞としては形態的特徴の異なる複数種が含まれてもよい。例えば幹細胞の培養において幹細胞とフィーダー細胞とが共培養された環境では、これらのいずれもが原画像に含まれる可能性がある。また、幹細胞の一部が体細胞に分化した場合には、幹細胞と体細胞とが混在する状態があり得る。

取得された原画像については、解像度変換が必要であるか否かが判断される（ステップＳ１０２）。解像度変換は、撮像時の原画像の画素サイズを領域分割処理に適した画素サイズに変換する処理である。必要と判断された場合には解像度変換が実行される（ステップＳ１０３）。

以下の説明において特に断りなく単に「画素」というとき、ステップＳ１０４以降の処理に供される原画像の画素を指すものとする。解像度変換が実行されていれば変換後の画素が、また解像度変換が実行されていなければ撮像時の画素が、ここでいう「画素」である。なお特に区別する必要があるとき、変換前の（つまり撮像時の）画素を「原画素」と称することがある。すなわち、「原画素」は、解像度変換が実行されたか否かに関わらず撮像時の画素を意味するものとする。

原画像は複数のブロックに分割される（ステップＳ１０４）。この実施形態では、原画像における４画素×４画素の正方形が１つのブロックとされる。解像度変換が行われたときには、上記した通り、ここでいう画素は変換後の画素である。１ライン分の画素数が１ブロックの１辺の画素数４で割り切れない場合には、原画像の端部において処理に算入されない画素が生じてもよい。また、このような端数が出ないように原画像の端部を予め切除しておいてもよい。

次に、各画素の画素値から、原画像全体における平均画素値がオフセット値として減算される（ステップＳ１０５）。画素値が例えば８ビットデータとして表されている場合、各画素は０から２５５までの２５６段階の値を取り得るが、平均画素値を減算することにより、画素値として負の値を取る画素も出現することになる。この処理は、後の演算の便宜のためのレベルシフトのためのものである。

こうして分割されたブロックごとに、空間周波数成分が求められる。この目的のために、本実施形態では二次元の離散コサイン変換（Discrete Cosine Transform；以下「ＤＣＴ」と略称する）演算が利用される（ステップＳ１０６）。例えばＪＰＥＧ（Joint Photographic Experts Group）画像圧縮技術でも用いられているように、二次元ＤＣＴ演算により各ブロックの画像内容が二次元の空間周波数成分に分解される。

当該ブロックの画像内容は、各周波数成分を合成することによって再現することが可能である。画像の内容により、周波数成分ごとの強度が異なる。ＪＰＥＧ方式の画像圧縮理論において知られているように、各周波数成分の強度は、ＤＣＴ演算により得られるＤＣＴ係数として求められる。

空間周波数成分のうち低周波成分は、ブロック内の画像における緩やかな明るさの変化を表す。このうち直流成分は、当該ブロックの明るさの平均値に相当する。一方、高周波成分はより細かい明るさの変化を表しており、ブロック内の画像のテクスチャに関する情報を持つ。高い周波数ほど細かいテクスチャを表す。そこで、求められたＤＣＴ係数に基づき、各周波数成分の強度に対応する特徴量を算出する（ステップＳ１０７）。具体的には、直流成分に対応するＤＣＴ係数を１つの特徴量とするとともに、交流成分についてはＤＣＴ係数を対応する周波数の高さに応じて幾つかのグループに分け、グループ内でのＤＣＴ係数の合計を１つの特徴量とする。グループ化の具体的な方法については後述する。

こうすることで、１つのブロックを、明るさを示す直流成分特徴量とテクスチャを示す複数の交流成分特徴量とによって表現することができる。これらの特徴量で表されるブロックがいずれかの種類の細胞としての特徴を強く有しているか、あるいは細胞以外の非細胞の特徴を有しているかにより、当該ブロックがどの種類の細胞に対応するか、またそれ以外の非細胞に対応するかを判定することができる（ステップＳ１０８）。その結果に基づき、原画像が複数の領域、例えば細胞の領域と非細胞の領域、あるいは細胞の種類ごとの領域に分割される（ステップＳ１０９）。以上が本実施形態の領域分割処理の概要である。

図２はこの実施形態の領域分割処理を実行するコンピュータ装置の構成例である。コンピュータ装置１は例えばパーソナルコンピュータとして一般的な構成を有するものであり、ＣＰＵ（Central Processing Unit）１０、メモリ１４、ストレージ１５、入力デバイス１６、表示部１７、インターフェース１８およびディスクドライブ１９などを備えている。

ＣＰＵ１０は、予め用意された制御プログラムを実行することで、領域分割処理を実行するための各機能ブロック、すなわち画像加工部１１、ＤＣＴ演算部１２および分類部１３などの機能ブロックをソフトウェア的に実現する。なお、これらが専用ハードウェアにより構成されてもよい。メモリ１４はＣＰＵ１０の演算過程で生成される各種データを一時的に記憶する。ストレージ１５は、ＣＰＵ１０が実行すべき制御プログラムのほか、原画像の画像データや処理後の画像データ等を長期的に記憶する。

入力デバイス１６は、オペレータからの指示入力を受け付けるためのものであり、例えばマウス、キーボードなどを含む。また、表示部１７は画像を表示する機能を有する例えば液晶ディスプレイであり、原画像や処理後の画像、オペレータへのメッセージ等種々の情報を表示する。なお、入力デバイスと表示部とが一体化されたタッチパネルが設けられてもよい。

インターフェース１８は、電気通信回線を介して外部装置との間で各種データ交換を行う。ディスクドライブ１９は、画像データや制御プログラム等各種のデータを記録した外部の記録ディスク２を受け入れる。記録ディスク２に記憶された画像データや制御プログラム等は、ディスクドライブ１９により読み出され、ストレージ１６に記憶される。ディスクドライブ１９はコンピュータ装置１内で生成されたデータを記録ディスク２に書き込む機能を備えていてもよい。

本実施形態の領域分割処理をコンピュータ装置１に実行させるための制御プログラムについては、これを記録した記録ディスク２にディスクドライブ１９がアクセスして読み出される態様であってもよく、インターフェース１８を介して外部装置から与えられる態様であってもよい。原画像データについても同様である。

画像加工部１１は、原画像の解像度変換、ブロック分割および最終的な画像の領域分割等を実行する。また、ＤＣＴ演算部１２は上記したＤＣＴ演算をブロックごとに実行して空間周波数成分を求める。分類部１３は、ブロックごとに求められた特徴量に基づき、各ブロックを分類する。

図３は原画像の例を示す図である。図３（ａ）に示すように、ステップＳ１０１において取得された原画像Ｉａが、画像オブジェクトとして細胞Ｃ１〜Ｃ５と、培地に対応する背景Ｂとを含むものであったとする。ここで、細胞Ｃ１と細胞Ｃ４とは同一種類であり、その内部は同じようなテクスチャを持つものとする。同様に、細胞Ｃ２と細胞Ｃ３とは同一種類であるとする。細胞Ｃ１，Ｃ２，Ｃ５は互いに異なる種類である。図においては細胞内部の濃度の違いが種類の違いを表している。

図において画像の横方向座標をＸ座標、縦方向座標をＹ座標とする。これらの画像オブジェクトの形状および濃淡は単に一例を示したものであり、実際に得られる画像では種々の形状および濃淡を取り得る。なお、これらの細胞Ｃ１〜Ｃ５と背景Ｂとの他に、デブリが原画像に含まれる場合もあり得る。また、画像が細胞Ｃ１〜Ｃ５と、これを取り巻くデブリとで占められ、培地に対応する背景が含まれない場合もあり得る。また、培地の量や撮像時の照明条件等により、背景Ｂの明るさも常に同じとは限らない。

これらのことから、画像の中での細胞および背景の配置や濃淡等の関係は、原画像ごとに異なっていると考えるべきである。ここに説明する画像処理は、これらの関係によらず成立するものである。この実施形態における領域分割処理の基本的な考え方は、原画像Ｉａを小さなブロックに分割し、各ブロックが有する空間周波数成分によって、当該ブロックがいずれかの種類の細胞の画像特徴を強く有しているか、それ以外の特徴を強く有しているかを判定するというものである。

この領域分割処理では、原画像Ｉａが所定サイズを有する正方形のブロックに分割されて処理される。図３（ｂ）に示すように、１つのブロックＢＬはＸ方向およびＹ方向に４画素ずつ、計１６の画素Ｐからなる正方形の領域である。ブロックＢＬの大きさは、画像中の抽出対象となる細胞の大きさに応じて設定されることが望ましい。

ブロックサイズが小さすぎると、１つのブロックは細胞のごく一部のみを含むこととなり、画像において細胞が有する特徴的な性質を十分に処理に反映させることができず、領域分割の精度低下を招く。また、ブロックサイズが大きすぎると、１つのブロックに細胞の特徴とそれ以外のオブジェクトの特徴とが混在することになり、やはり精度低下の原因となる。本願発明者の知見では、ブロックＢＬの大きさは、画像中の抽出対象となる細胞の内部に概ね１つのブロックが入る程度に選ばれるのが好ましい。

細胞が例えばヒトの体細胞である場合、その直径は概ね１０μｍ〜２０μｍ程度であるが、幹細胞はこれよりもう少し小さい。したがって、ブロックＢＬの大きさも、試料サイズに換算したときの１辺の長さがこれと同程度か少し小さいことが好ましい。このことから１つの画素Ｐのサイズは、１辺が２．５μｍ〜５μｍまたはこれより少し小さいことが好ましい。

撮像における解像度を例えば６４００ｄｐｉ（dots per inch）とすると、画素サイズは約４μｍである。また、解像度を１６０００ｄｐｉとすると、画素サイズは約１．６μｍである。このことから、上記した画素Ｐのサイズを得るために必要な解像度としては６４００ｄｐｉないし１６０００ｄｐｉが好適である。ただし、これらは本実施形態の領域分割処理を行うのに適した値であって、撮像時にはこれに制約されずできるだけ高い解像度で撮像を行うことが基本である。

例えば３２０００ｄｐｉ（画素サイズ０．８μｍ）の解像度で原画像の撮像が可能であるとき、領域分割処理の都合のために敢えて解像度を低下させる必然性はない。この場合、図３（ｃ）に示すように、撮像時の原画像Ｉａにおける４画素を擬似的に１つの画素Ｐとみなし、それらの平均画素値を当該画素Ｐの画素値とする。こうすることで、撮像時の解像度を低下させることなく、かつ領域分割処理に適した画素サイズ（１．６μｍ）を得ることができる。本実施形態ではこのような処理を「解像度変換」と称している。

解像度変換が必要であるか否かは、撮像時の画素サイズと処理対象である細胞のサイズとの関係によって決まる。これらの関係が予めわかっている場合には、その既知の関係から解像度変換が必要であるか否かを判断することができる。また、解像度変換を行うかどうかをユーザからの指示入力に基づいて決定するようにしてもよい。必要と判断された場合には解像度変換が実行される（ステップＳ１０３）。

解像度変換が行われたとき、ステップＳ１０４では、変換後の画像における４画素×４画素の正方形が１つのブロックＢＬとされる。解像度変換が行われていなければ、原画像Ｉａの撮像時の画像における４画素×４画素の正方形が１つのブロックＢＬとされる。Ｘ方向およびＹ方向の少なくとも一方において、１ライン分の画素数が１ブロックの１辺の画素数４で割り切れない場合には、原画像の端部において処理に算入されない画素が生じてもよい。また、このような端数が出ないように原画像の端部を予め切除しておいてもよい。

次に、各画素の画素値から、原画像Ｉａ全体における平均画素値がオフセット値として減算される（ステップＳ１０５）。画素値が例えば８ビットデータとして表されている場合、各画素は０から２５５までの２５６段階の値を取り得るが、平均画素値を減算することにより、画素値として負の値を取る画素も出現することになる。本願発明者の知見では、画像における最大画素値と最小画素値との差が比較的小さい細胞画像では、各画素から画像全体の平均画素値を差し引くことで、各画素値の正負のバランスが良好なものとなる。各画素Ｐの画素値をｐ（Ｘ，Ｙ）、平均画素値をｐ_avとするとき、減算後の画素Ｐの画素値ｇ（Ｘ，Ｙ）は次式（式１）で表される。

こうして分割されたブロックごとに、ＤＣＴ演算が実行される（ステップＳ１０６）。ＤＣＴ演算により、各ブロックの画像内容が、Ｘ方向およびＹ方向それぞれについて２次元の空間周波数成分に分解される。ＤＣＴ演算処理はＪＰＥＧ方式の画像圧縮処理でも利用されているように広く知られたものであるため、以下ではＪＰＥＧ画像圧縮処理と同じ考え方を適用することができる事項については説明を省略する。

図４はＤＣＴ演算の概念を示す図である。図左は、原画像の一部であって（４×４）画素からなる、１つのブロックＢＬを表している。原画像での座標位置と区別してブロック内での各画素の位置を表すために小文字を用い、横方向座標をｘ座標、縦方向座標をｙ座標とする。ＤＣＴ演算は、ブロックＢＬ内の各画素の画素値ｇ（ｘ，ｙ）から、図右に示す４行×４列のＤＣＴ係数テーブルＴを求める演算と考えることができる。

ＤＣＴ係数テーブルＴにおける行位置ｕ（ｕ＝０，１，２，３）はｘ方向における空間周波数成分を表し、具体的には、ｕ＝０は直流成分を、ｕ＝１は基本周波数成分を、ｕ＝２，３はそれぞれ基本周波数の２倍、３倍の成分を表す。同様に、列位置ｖ（ｖ＝０，１，２，３）はｙ方向における空間周波数成分を示しており、ｖ＝０は直流成分を、ｖ＝１は基本周波数成分を、ｖ＝２，３はそれぞれ基本周波数の２倍、３倍の成分を表す。

各周波数成分の強度を表すＤＣＴ係数Ｃ（ｕ，ｖ）は、次式（式２）により表される。

ここで、定数Ｎはブロックサイズであり、ＪＰＥＧ方式では、Ｎ＝８、すなわち（８×８）画素が１ブロックとされるが、前述したように、この実施形態では（４×４）画素を１ブロックとしているので、Ｎ＝４である。本願発明者の知見では、画像の圧縮・復元を目的とせずＤＣＴ係数を特徴量として利用する本処理方法では、ブロックサイズを（４×４）画素とした方が良好な結果が得られている。またこうすることで演算量も削減することができ高速化にも適している。

ＤＣＴ係数テーブルＴを構成するＤＣＴ係数Ｃ（ｕ，ｖ）のうち、直流成分を表す係数Ｃ（０，０）が直流成分特徴量とされる。また、交流成分を表す各係数は複数のグループに分けられ、グループ内で合計される。この合計値、または合計値を合計された係数の数で除して正規化した値が、当該グループに対応する１つの特徴量とされる。なお、ＤＣＴ係数は正負いずれの値も取り得る。合計時にそれらの値が相殺されることを回避するためには、例えばＤＣＴ係数の合計をそれぞれの絶対値の和とすることができる。

図５はＤＣＴ係数のグループ化の例を示す図である。基本的な考え方は、図５（ａ）に示すように、最大周波数が基本周波数（ｕ＝１、ｖ＝１）である３つのＤＣＴ係数、符号で表せばＣ（１，０）、Ｃ（０，１）、Ｃ（１，１）を第１のグループＧ１とし、最大周波数が基本周波数の２倍（ｕ＝２、ｖ＝２）である５つのＤＣＴ係数Ｃ（２，０）、Ｃ（２，１）、Ｃ（０，２）、Ｃ（１，２）、Ｃ（２，２）を第２のグループＧ２とし、上記および直流成分Ｃ（０，０）を除く７つのＤＣＴ係数Ｃ（３，０）、Ｃ（３，１）、Ｃ（３，２）、Ｃ（０，３）、Ｃ（１，３）、Ｃ（２，３）、Ｃ（３，３）を第３のグループＧ３とする、というものである。これは、画像において同程度のテクスチャの細かさを表す周波数成分の強度を合計する、との考え方に基づくものである。

各グループＧ１，Ｇ２，Ｇ３に含まれるＤＣＴ係数をそれぞれ合計した値、または合計値を当該グループのＤＣＴ係数の個数（例えばグループＧ１では３）で除して正規化した値を、当該周波数に対応する交流成分特徴量とすることができる。各グループを代表する周波数は、例えば当該グループにおける最大周波数により表すことができる。

正規化を行う理由は、後の分類工程において各周波数成分が同じ重みで評価されるようにするためである。これにより、特定の周波数成分に対応するテクスチャのみが分類結果に大きな影響を与えることが回避される。なお直流成分特徴量はテクスチャに関わるものでなくブロック全体の輝度を表すものであるため、テクスチャを表す交流成分特徴量との間での正規化は必ずしも必要ない。

また、図５（ｂ）に示すように、例えば直流成分の係数Ｃ（０，０）および第１グループＧ１を除く１２個の係数を全て第２グループＧ２に含める、というように、上記原則に限定されず、適宜にグループ化を実行することが可能である。

ただし、図５（ａ）に破線で示すように、ＤＣＴ係数テーブルＴにおいて変数ｕ、ｖを座標と見たときの直線ｖ＝ｕに対して互いに対称な位置のＤＣＴ係数は、同じグループに含まれることが望ましい。というのは、変数ｕは画像のｘ座標、変数ｖは画像のｙ座標に対応しており、撮像対象物である培養細胞等の画像においては原理的にｘ方向とｙ方向との間で特別な異方性はないことから、ｕ方向とｖ方向とで扱いを異ならせるべきではないからである。

また、図５（ｃ）に示すグループ化の例では、係数Ｃ（１，０）、Ｃ（２，０）、Ｃ（０，１）、Ｃ（１，１）、Ｃ（０，２）の５つを第１グループＧ１、これらに隣接する係数Ｃ（３，０）、Ｃ（２，１）、Ｃ（１，２）、Ｃ（２，２）、Ｃ（０，３）の５つを第２グループＧ２、他の５つの係数Ｃ（３，１）、Ｃ（３，２）、Ｃ（１，３）、Ｃ（２，３）、Ｃ（３，３）を第３グループＧ３に含ませている。このようにすると、各グループにおいて値が合計されるＤＣＴ係数の個数が全て同じとなり、各グループで求められる交流成分特徴量の間で重みが同じなることから、正規化を行う必要はなくなる。

構成されるグループの数が多いほど画像の特徴を表す情報は多く保存されるが、特徴量空間の次元が高くなるため後の処理における演算量が増大する。周波数成分が比較的近い係数同士をグループ化することで、画像の情報を維持しつつ、演算負荷を低減することが可能である。

上記のようにして求められた各周波数成分の特徴量によって、各ブロックＢＬの輝度およびテクスチャを表す特徴量ベクトルが構成される。すなわち、当該ブロックＢＬの明るさを表す直流成分特徴量と、ブロック内のテクスチャを表す複数の交流成分特徴量とが特徴量ベクトルをなす。

この特徴量ベクトルに基づき、各ブロックが複数種の「細胞」、「非細胞」のいずれのカテゴリに相当するかが判定される（ステップＳ１０８）。この判定は、各周波数成分の特徴量の組み合わせにより構成される多次元特徴量空間において当該ブロックに対応する特徴量ベクトルが示す点の位置に基づいてなされる。判定を行うためには、その判定基準、具体的には多次元特徴量空間において各カテゴリの特徴が占める範囲を特定しておく必要がある。

複数の撮像の間で、撮像される細胞の種類、培地の種類および量、撮像条件（撮像倍率や照明強度など）等が固定されている場合には、同じ条件で撮像された典型的な細胞の画像を事例として多数収集し、適宜の機械学習アルゴリズムを用いることで、多次元特徴量空間において各カテゴリに対応するクラスタが占める範囲を予め特定しておくことができる。この場合、未知の原画像から選出されたブロックの特徴量がどのクラスタに属するかによって、当該ブロックがどの種類の細胞あるいは非細胞に対応するかを判定することが可能である。

ただし、一般的には撮像条件等は原画像ごとに異なる。このため、より精度の高い判定を行うためには、処理対象である原画像から読み取れる情報から判定基準を策定することができると便宜である。以下、その方法の一例について説明する。

図６は判定基準を策定するための処理を示すフローチャートである。また、図７は処理過程における指示入力態様を例示する図である。この処理は、取得された原画像Ｉａのうち複数種の細胞それぞれの特徴を強く有する典型例についてオペレータから教示を受け、その情報から多次元特徴量空間におけるカテゴリ間の境界を特定するものである。この処理は、コンピュータ装置１が原画像Ｉａを取得した時（ステップＳ１０１）以降、ブロックの分類（ステップＳ１０８）を開始するまでの任意のタイミングで実行することができる。処理は、ＣＰＵ１０がストレージ１５に記憶されている制御プログラムを実行することにより実現される。

最初に、取得された原画像Ｉａが表示部１７に表示されオペレータに提示される（ステップＳ２０１）。そして、表示された原画像Ｉａの中で、細胞の領域に属する代表点と、非細胞の領域に対応する代表点とを指定するオペレータからの教示入力を受け付ける（ステップＳ２０２）。細胞が複数種ありこれらを識別する必要がある場合、種類ごとに教示入力を受け付ける。なお、細胞の領域と非細胞の領域とがはっきりした外見上の差異を有している場合には、細胞の領域に属する代表点のみが指定されてもよい。

精度向上のために、代表点は複数指定されることが望ましい。そこで、図７（ａ）に示すように、例えば、マウス操作あるいはタッチペン入力等により細胞Ｃ２の領域内で連続的にポインタを動かす動作、または何点かを選択する動作をオペレータに実行させることで、その軌跡上の複数の点を細胞Ｃ２の代表点として取得収集することができる。図における白抜き矢印は、ユーザによるマウス等の操作に連動して画面上で移動するポインタを表している。種類の異なる細胞Ｃ４，Ｃ５や背景Ｂについても、同様にして代表点の指定を受け付ける。

また、図７（ｂ）に点線で示すように、細胞内に属する矩形領域、細胞外の矩形領域等を二次元の領域として指定させてもよい。画像内にデブリが存在する場合、細胞外の領域の指定においてはデブリの領域と背景の領域との両方を含ませることが望ましい。こうすることにより、細胞外の領域から抽出される代表点にはデブリの特徴を有するものと背景の特徴を有するものとがともに含まれることとなる。そのため、特徴量空間における細胞の領域と非細胞の領域との境界を決めるに当たり、非細胞に関する情報をより強く反映させることができる。もちろんデブリと背景とで個別に代表点を指定するようにしてもよい。

こうしていくつかの代表点が指定されると、上記した領域分割処理と同様に、各画素の画素値から原画像Ｉａ全体における平均画素値が減算され（ステップＳ２０３）、代表点を含む（４×４）画素サイズの仮想的なブロックにつきＤＣＴ演算が実行され（ステップＳ２０４）、ＤＣＴ係数から直流成分特徴量および複数の交流成分特徴量が算出される（ステップＳ２０５）。これにより、代表点を含むブロックを典型例としたときの多次元特徴量空間における分布が明らかになり、その結果から、多次元特徴量空間における各種の「細胞」のクラスタが占める領域と「非細胞」のクラスタが占める領域とが特定される（ステップＳ２０６）。

図８は特徴量空間におけるクラスタ分布の例を模式的に示す図である。特徴量空間の次元は種々の値を取り得るが、ここでは３次元の例を図示している。図７（ａ）または図７（ｂ）に示した方法を用いて、ユーザにより複数種の細胞領域および非細胞領域が指定されると、それらの領域内の代表点を含むいくつかのブロックにつき直流成分特徴量および交流成分特徴量が求められる。各ブロックを多次元特徴量空間にプロットすると、図８（ａ）に示すように、各種の細胞および非細胞それぞれの画像特徴の差異に起因して、細胞Ｃ２に対応するクラスタＣＬ１と、細胞Ｃ４に対応するクラスタＣＬ２と、細胞Ｃ５に対応するクラスタＣＬ３と、それ以外の非細胞、すなわち背景Ｂおよびデブリに対応するクラスタＣＬ４とが形成されると考えられる。

そして、原画像Ｉａを構成する各ブロックＢＬの分類に際しては、当該ブロックを表す特徴量ベクトルが多次元特徴量空間において占める位置が、どのクラスタに属するべきものかにより、当該ブロックがいずれかのクラスタに分類される。例えば、特徴量ベクトルが指し示す位置がいずれかのクラスタ領域の内部にあれば、そのブロックは当該クラスタに属するものと判定することができる。

また、図８（ｂ）に示すように、分類すべきブロックの特徴量ベクトルにより位置を示される点Ｐａと各クラスタとの距離をそれぞれ求めて、その距離が最も近いクラスタに分類するようにしてもよい。点Ｐａとクラスタとの距離については、例えば次のようにして求めることができる。各クラスタにつき、それぞれ距離基準点Ｐ１〜Ｐ４を予め定めておく。ここで、距離基準点Ｐ１〜Ｐ４としては、例えば各クラスタが特徴量空間において占める空間領域の重心点とすることができる。また、ユーザから教示された代表点の平均位置に対応する点を距離基準点としてもよい。そして、点Ｐａと距離基準点との距離によって、点Ｐａと各クラスタとの距離を表すことができる。

以上のように、この実施形態の領域分割処理では、細胞を含む原画像を複数の小ブロックに分割し、ブロックごとの空間周波数成分を求める。そして、空間周波数成分のうち直流成分の強度に基づく直流成分特徴量と、交流成分のうちの複数の周波数成分それぞれの強度に基づく複数の交流成分特徴量とを含む特徴量ベクトルが多次元特徴量空間において示す位置が、どのクラスタに属するものであるかによって、当該ブロックの画像が細胞に対応するものかそれ以外の非細胞に対応するものか、あるいはどの種類の細胞に対応するかが判定される。

直流成分特徴量はブロックの明るさを表す。また、交流成分特徴量のうち比較的周波数の低い成分に対応するものは比較的単純な構造を、また高周波成分に対応するものはブロック内の細かいテクスチャをそれぞれ表す特徴量である。このように明るさとテクスチャとを総合的に評価することで、単なる２値化処理や輪郭抽出処理による領域分割よりも優れた精度で、画像内の細胞と非細胞とを区分することができ、複数種類の細胞が画像に含まれる場合でも、それらを区別することができる。

なお、本発明は上記した実施形態に限定されるものではなく、その趣旨を逸脱しない限りにおいて上述したもの以外に種々の変更を行うことが可能である。例えば、上記実施形態における原画像の領域分割（ステップＳ１０９）では、原画像がブロック単位で細胞や非細胞に分割されることになる。しかしながら、ブロック境界は自動的に策定されたものであり、原画像における細胞とその周囲との境界に一致するとは限らない。このことに鑑み、例えば公知の輪郭抽出処理と組み合わせて原画像の分割が行われてもよい。すなわち、上記のようにしてブロック単位で特定された領域の境界の位置を、輪郭抽出処理の結果によって修正するようにしてもよい。

また、上記実施形態における多次元特徴量空間でのクラスタを特定する処理は単なる例示であり、これ以外にも各種のクラスタリング方法を適用することが可能である。また、細胞および非細胞の典型例を収集する方法についても、上記に限定されるものではない。

また、上記実施形態では、離散コサイン変換処理によって求められるＤＣＴ係数によって各周波数成分の強度を表しているが、ＤＣＴ演算以外にも適宜の空間周波数成分を求める方法があれば、それが利用されてもよい。

また、上記実施形態では、汎用のコンピュータ装置１によって領域分割処理が実行されているが、前述したように例えば撮像装置にこの処理機能が組み込まれていてもよい。また、既存の撮像装置に本実施形態の領域分割処理を実行するためのプログラムを追加的に実装して機能拡張を行うことも可能である。

また、上記説明では画素のサイズやブロックサイズの具体的数値を例示しているが、本実施形態の領域分割処理は、これらの数値範囲外の条件であっても実行することが可能である。また、画像に含まれる細胞の種類が単一であり、細胞とその他の領域とを分割する際にも、本発明の処理方法を適用可能である。

以上、具体的な実施形態を例示して説明してきたように、本発明においては、ブロックごとに、当該ブロックに含まれる画素の画素値に基づく二次元離散コサイン変換を実行し、各周波数成分の離散コサイン変換係数に基づき特徴量が求められてもよい。このような処理はＪＰＥＧ方式の画像圧縮処理にも利用されており、演算アルゴリズムが既に確立されているため、本発明にもこれを適用することができる。

また、二次元離散コサイン変換により求められる複数の離散コサイン変換係数を、空間周波数の高さに応じてグループ化した複数のグループそれぞれの離散コサイン変換係数の合計、または合計された離散コサイン変換係数の数に応じて合計を正規化した値を特徴量としてもよい。このような構成によれば、例えば比較的周波数成分の近い離散コサイン変換係数同士を１つにまとめることで、画像の特徴に関する情報を失うことなく特徴量空間の次元を低くし、後工程での演算負荷を低減することができる。

この場合、画素それぞれの画素値から原画像全体の平均画素値を減じた値に基づき二次元離散コサイン変換が実行されてもよい。透明に近い細胞を撮像した画像では画素値の分布幅が小さいため、離散コサイン変換演算のためのレベルシフトは、原画像全体の平均画素値を用いて行われるのが好ましい。

また、複数の周波数成分の強度の合計は、各周波数成分の離散コサイン変換係数の絶対値の和として求められてもよい。離散コサイン変換係数は正負いずれの値をも取り得るため、単純な合計では各周波数成分が相殺されて良好な結果を得られないおそれがある。絶対値の和を取ることにより、この問題を解消することができる。

また、各ブロックは、１辺が４画素に相当するブロックサイズの正方形であってもよい。ＪＰＥＧ方式の画像圧縮ではブロックサイズは８画素とされているが、本願発明者の知見では、画像の圧縮・復元を目的としない本発明においてはブロックサイズを４画素とした方が良好な結果が得られる。

また、原画像をブロックに分割するのに先立って、原画像の画素サイズを、取得したときの原画素サイズを、原画素サイズと異なるブロック分割用の画素サイズに変換する工程がさらに設けられてもよい。細胞とそれ以外の非細胞との間における明るさやテクスチャの違いを利用して領域分割を行う本発明では、演算時の好ましい画素サイズが例えば細胞等、画像中のオブジェクトのサイズに依存する。このための画素サイズは撮像時の好ましい画素サイズとは必ずしも一致しないが、画素サイズの変換を行うことにより、それぞれ個別に画素サイズの最適化を図ることができる。

また、分類に先立って、原画像中で複数のクラスタ各々に対応する代表点に関するユーザからの教示入力を受け付ける工程と、教示入力の結果に基づき、代表点を含むブロックについて特徴量を求め、その結果に基づき多次元特徴量空間におけるクラスタ各々の範囲を特定する工程とがさらに設けられてもよい。このような構成によれば、細胞の種類や撮像条件が一定でない場合でも、細胞やそれ以外の領域の分割を安定的に行うことが可能である。

また、原画像は細胞を明視野撮像することで得られた画像であってもよい。一般に細胞は透明に近く、また培地との屈折率の差も小さいため、明視野撮像画像では目視による細胞と非細胞との区別がつきにくい。このような画像に本発明を適用することで、画像中の細胞の領域とそれ以外の領域とを良好にかつ安定的に分割することが可能となる。

本発明は、細胞を撮像し画像を評価する生化学や医療の分野に適用可能であり、特に培養されている細胞にダメージを与えるような加工を行うことのない撮像が必要とされる技術分野に好適なものである。

１ コンピュータ装置
２ 記録ディスク（記録媒体）
Ｂ 背景
ＢＬ ブロック
Ｃ１〜Ｃ５ 細胞
Ｉａ 原画像

Claims (11)

1. 培養された細胞が含まれる原画像を取得する工程と、
   前記原画像を所定数の画素からなるブロックに分割し、各ブロック内の画像が有する空間周波数成分を求める工程と、
   前記ブロックの各々を、前記空間周波数成分のうちの直流成分の強度と、空間周波数が互いに異なる複数の交流成分それぞれの強度とを特徴量とする多次元特徴量空間における複数のクラスタに分類する工程と、
   前記原画像を、前記分類の結果に基づき複数の領域に分割する工程と
   を備える画像処理方法。
2. 前記ブロックごとに、当該ブロックに含まれる前記画素の画素値に基づく二次元離散コサイン変換を実行し、各周波数成分の離散コサイン変換係数に基づき前記特徴量を求める請求項１に記載の画像処理方法。
3. 前記二次元離散コサイン変換により求められる複数の前記離散コサイン変換係数を、空間周波数の高さに応じてグループ化した複数のグループそれぞれの前記離散コサイン変換係数の合計、または合計された前記離散コサイン変換係数の数に応じて前記合計を正規化した値を前記特徴量とする請求項２に記載の画像処理方法。
4. 前記画素それぞれの画素値から前記原画像全体の平均画素値を減じた値に基づき前記二次元離散コサイン変換を実行する請求項２または３に記載の画像処理方法。
5. 前記合計は、前記離散コサイン変換係数の絶対値の和として求められる請求項３に記載の画像処理方法。
6. 前記ブロックは、１辺が４画素に相当するブロックサイズの正方形である請求項２ないし５のいずれかに記載の画像処理方法。
7. 前記原画像を前記ブロックに分割するのに先立って、前記原画像の画素サイズを、取得したときの原画素サイズから、前記原画素サイズと異なるブロック分割用の画素サイズに変換する工程を備える請求項１ないし６のいずれかに記載の画像処理方法。
8. 前記分類に先立って、
   前記原画像中で前記複数のクラスタ各々に対応する代表点に関するユーザからの教示入力を受け付ける工程と、
   前記教示入力の結果に基づき、前記代表点を含む前記ブロックについて前記特徴量を求め、その結果に基づき前記多次元特徴量空間における前記クラスタ各々の範囲を特定する工程と
   を備える請求項１ないし７のいずれかに記載の画像処理方法。
9. 前記原画像は、前記細胞を明視野撮像することで得られた画像である請求項１ないし８のいずれかに記載の画像処理方法。
10. コンピュータに、請求項１ないし９のいずれかに記載の画像処理方法の各工程を実行させるためのコンピュータプログラム。
11. 請求項１０に記載のコンピュータプログラムを記録した、コンピュータ読み取り可能な記録媒体。

Images