JP2015052663A - 画像処理方法、画像処理装置、撮像装置およびプログラム - Google Patents

Abstract

【課題】明視野顕微鏡の観察で得られた１枚の画像から標本の透過光の振幅および位相分布を高速かつ高精度に再構成すること。【解決手段】部分コヒーレントまたはコヒーレントな結像光学系が形成した未知の標本の評価画像を、複数の既知の標本の画像から生成された基底の線形結合により近似するための結合係数を算出するステップ（Ｓ２１５）と、複数の既知の標本の各々に関する複素量のデータと前記結合係数とに基づき中間データを算出するステップ（Ｓ２１７）と、中間データから前記未知の標本に関する複素量のデータを算出するステップ（Ｓ２１８）と、を有する画像処理方法を提供する。【選択図】図３

Description

本発明は、部分コヒーレントまたはコヒーレントな結像光学系を介して撮像した標本の画像処理方法、画像処理装置、撮像装置およびプログラムに関する。

生物学や病理診断に用いられる標本の観察に使用される顕微鏡として位相コントラスト方式、微分干渉方式などの方式は、標本の内部構造の詳細を把握するために透過光の強度および位相変化に関する定量的な情報を得ることは困難であった。この理由の１つとして、観測過程において標本の透過光の振幅および位相分布と像の間に非線形な関係があるために、逆問題を解く計算が困難になることが挙げられる。なお以降では、逆問題とは観測データ（像）から原信号（標本の情報）を数値的に推定する問題を意味し、逆問題を解くことを再構成と呼ぶ。

そこで、近年では標本の位相情報を定量化する方式が提案されている。例えば、デジタルホログラフィ顕微鏡を利用した干渉法（特許文献１）、標本をレーザービームで照明し空間光変調器（ＳＬＭ；Ｓｐａｔｉａｌ Ｌｉｇｈｔ Ｍｏｄｕｌａｔｏｒ）を挿入した光学系を介して複数枚の画像を取得する手法（非特許文献１）がある。一般に干渉法は、標本を透過していない参照光と標本を透過した物体光との干渉の結果像面に形成される光強度分布を観測し、取得した画像データに演算を行い標本の情報を再構成する手法である。一方、レーザービームで照明し物体光を直接撮像する手法は参照光を用いる必要はないが、焦点ずれを変えて取得した複数の画像データに演算を行い、標本の情報を再構成する。

また、非線形なモデルに対する逆問題を解く新しい計算方法として、Ｋｅｒｎｅｌ Ｃｏｍｐｒｅｓｓｉｖｅ Ｓｅｎｓｉｎｇ（以下、「ＫＣＳ」と略す）と呼ばれる手法（非特許文献２）が挙げられる。Ｃｏｍｐｒｅｓｓｉｖｅ Ｓｅｎｓｉｎｇ（以下、「ＣＳ」と略す）は、線形な観測過程を前提として観測データ量を大幅に減らす手法である。一方、ＫＣＳは、データに適切な非線形写像を行うことで疎表現（Ｓｐａｒｓｅ Ｒｅｐｒｅｓｅｎｔａｔｉｏｎ）が可能になるような場合において、観測データ量を更に減らす手法である。しかしながらＫＣＳは非線形な観測過程を対象とはしていないため、前述した顕微鏡の逆問題に直接適用することはできない。一方で、Ｃｏｍｐｒｅｓｓｉｖｅ Ｓｅｎｓｉｎｇの概念を用いて、通常の明視野顕微鏡（部分コヒーレント結像）の逆問題を解く方法が、非特許文献３に開示されている。また、非特許文献４にはカーネル法の一般的な定義が記載されている。

特許第４７７２９６１号

Ｖｌａｄｉｍｉｒ Ｋａｔｋｏｖｎｉｋ ａｎｄ Ｊａａｋｋｏ Ａｓｔｏｌａ著， 「Ｐｈａｓｅ ｒｅｔｒｉｅｖａｌ ｖｉａ ｓｐａｔｉａｌ ｌｉｇｈｔ ｍｏｄｕｌａｔｏｒ ｐｈａｓｅ ｍｏｄｕｌａｔｉｏｎ ｉｎ ４ｆ ｏｐｔｉｃａｌ ｓｅｔｕｐ： ｎｕｍｅｒｉｃａｌ ｉｎｖｅｒｓｅ ｉｍａｇｉｎｇ ｗｉｔｈ ｓｐａｒｓｅ ｒｅｇｕｌａｒｉｚａｔｉｏｎ ｆｏｒ ｐｈａｓｅ ａｎｄ ａｍｐｌｉｔｕｄｅ」， Ｊｏｕｒｎａｌ ｏｆ ｔｈｅ Ｏｐｔｉｃａｌ Ｓｏｃｉｅｔｙ ｏｆ Ａｍｅｒｉｃａ Ａ， アメリカ合衆国， Ｏｐｔｉｃａｌ Ｓｏｃｉｅｔｙ ｏｆ Ａｍｅｒｉｃａ， ２０１２年， ２９巻， １号， ｐ．１０５−１１６ Ｋａｒｔｈｉｋｅｙａｎ Ｎａｔｅｓａｎ Ｒａｍａｍｕｒｔｈｙ ａｎｄ Ａｎｄｒｅａｓ Ｓｐａｎｉａｓ著， 「Ｏｐｔｉｍｉｚｅｄ Ｍｅａｓｕｒｅｍｅｎｔｓ ｆｏｒ Ｋｅｒｎｅｌ Ｃｏｍｐｒｅｓｓｉｖｅ Ｓｅｎｓｉｎｇ」， Ｐｒｏｃｅｅｄｉｎｇｓ ｏｆ Ａｓｉｌｏｍａｒ Ｃｏｎｆｅｒｅｎｃｅ ｏｎ Ｓｉｇｎａｌｓ， Ｓｙｓｔｅｍｓ， ａｎｄ Ｃｏｍｐｕｔｅｒｓ， アメリカ合衆国， ２０１１年， ｐ．１４４３―１４４６ Ｙ． Ｓｈｅｃｈｔｍａｎ， Ｙ． Ｃ． Ｅｌｄａｒ， Ａ． Ｓｚａｍｅｉｔ， ａｎｄ Ｍ． Ｓｅｇｅｖ著， 「Ｓｐａｒｓｉｔｙ ｂａｓｅｄ ｓｕｂ−ｗａｖｅｌｅｎｇｔｈ ｉｍａｇｉｎｇ ｗｉｔｈ ｐａｒｔｉａｌｌｙ ｉｎｃｏｈｅｒｅｎｔ ｌｉｇｈｔ ｖｉａ ｑｕａｄｒａｔｉｃ ｃｏｍｐｒｅｓｓｅｄ ｓｅｎｓｉｎｇ」， Ｏｐｔｉｃｓ Ｅｘｐｒｅｓｓ， アメリカ合衆国， Ｏｐｔｉｃａｌ Ｓｏｃｉｅｔｙ ｏｆ Ａｍｅｒｉｃａ， ２０１１年， １９巻， ｐ．１４８０７−１４８２２ Ｊｏｈｎ Ｓｈａｗｅ−Ｔａｙｌｏｒ， Ｎｅｌｌｏ Ｃｒｉｓｔｉａｎｉｎｉ著， 「Ｋｅｒｎｅｌ Ｍｅｔｈｏｄｓ ｆｏｒ Ｐａｔｔｅｒｎ Ａｎａｌｙｓｉｓ」， Ｃａｍｂｒｉｄｇｅ Ｕｎｉｖｅｒｓｉｔｙ Ｐｒｅｓｓ， アメリカ合衆国， ２００４年

干渉法は光の干渉が環境の変化に敏感であるために高精度の調整が必要であり、観測データにノイズが加わりやすい。また多くの場合、参照光と物体光の２つの光路が必要であるために装置の複雑化や高コスト化につながる。また、レーザービームで照明し撮像する手法は、焦点ずれを正確に制御して複数枚の画像を取得する必要があるため、機械的制御手段やデータ量の多さが装置コストやデータ処理時間を増大させる。さらには、標本にレーザービームを照射するコヒーレント照明であるために、スペックルノイズが発生する。スペックルを抑制するために拡散板を光路中に挿入すると解像力が低下する。また、いずれの手法にも共通の課題として、以下の２点が挙げられる。１つは標本の情報を再構成するための演算に反復処理が必要であるため計算コストが大きいこと、もう１つは同じ装置構成で通常の顕微鏡観察はできないことである。非特許文献３の方法の課題として、以下の２点が挙げられる。１つは位相が２値しかとらず振幅分布が十分に疎である標本にしか適用できないこと、もう１つは前述と同じで計算コストが大きいことである。

本発明は、明視野顕微鏡から得られた１枚の画像に基づき標本の透過光の振幅および位相分布を高速かつ高精度に再構成することが可能な画像処理方法、画像処理装置、撮像装置およびプログラムを提供することを例示的な目的とする。

本発明の画像処理装置は、部分コヒーレントまたはコヒーレントな結像光学系が形成した既知の複数の第１の標本の光学像を光電変換して得られる複数の第１の画像から生成された基底の線形結合により、第２の画像を近似するための結合係数を算出する第１の算出部と、前記複数の第１の標本のデータの非線形写像によって得られる複数の複素量のデータと、前記第１の算出部によって算出された前記結合係数と、に基づいて中間データを算出する第２の算出部と、前記第２の算出部によって算出された前記中間データに基づいて未知の第２の標本に関する複素量のデータを算出する第３の算出部と、を有することを特徴とする。

本発明によれば、明視野顕微鏡から得られた１枚の画像に基づき標本の透過光の振幅および位相分布を高速かつ高精度に再構成することが可能な画像処理方法、画像処理装置、撮像装置およびプログラムを提供することができる。

本実施形態の撮像装置のブロック図である。 本実施形態の画像処理の概念を説明する図である。 本実施形態の画像処理を説明するフローチャートである。 実施例のシミュレーションに用いる照明および光学素子の図である。 本発明の訓練データを示す図である。（実施例１） 本発明の真の振幅、位相分布および評価画像を示す図である。（実施例１） 本発明の再構成された振幅および位相分布を示す図である。（実施例１） 本発明の評価画像を示す図である。（実施例２） 本発明の再構成された振幅および位相分布を示す図である。（実施例２） 本発明の真の振幅、位相分布および評価画像を示す図である。（実施例３） 本発明の再構成された振幅および位相分布を示す図である。（実施例３） 図３（ｃ）のフローを用いた場合の再構成結果を示す図である。（実施例３） 本発明の再構成された振幅および位相分布を示す図である。（実施例３）

本発明は、撮像装置に関するものであり、特に、明視野顕微鏡から取得されたデジタル画像から標本の情報を再構成するシステムに関する。本実施形態の画像処理部（画像処理装置）は、記憶手段と結合係数演算部（第１の算出部）と中間データ生成部（第２の算出部）と変換部（第３の算出部）と判断部とを有する。特に、結合係数演算部が評価画像を基底の線形結合で近似する点、中間データ生成部および変換部が複素量のデータを出力する点に特徴を有している。このような撮像装置または撮像システムはデジタル顕微鏡として好適であり、例えば、医学・生物学の研究や病理診断の用途に有用である。

以下、添付図面を参照して、本発明の一実施形態を説明する。

図１は、本発明の実施形態に係る撮像システムの概略構成を示している。図１に示すように、本実施形態の撮像システムは、撮像装置１０、コンピュータ（ＰＣ）５０１、画像処理部（画像処理装置）５０２、記憶装置５０３から構成され、ＰＣ５０１には入力装置５０４および表示装置５０５が接続されている。ただし図１のシステム構成は一例に過ぎない。

撮像装置１０は、照明光学系１００、光源１０１、標本ステージ２０１、ステージ駆動部２０２、結像光学系３００、光学素子３０１、撮像素子４０１から構成される。

光源１０１として、例えば、ハロゲンランプやＬＥＤ（Ｌｉｇｈｔ Ｅｍｉｔｔｉｎｇ Ｄｉｏｄｅ）を用いることができる。

照明光学系１００は、絞りなどの照明光変調装置を備えていてもよい。明視野顕微鏡の光学系において、照明光学系１００の絞りは解像力や焦点深度を変化させる。そのため、標本の種類やユーザのニーズに応じて観察画像を調整する目的で、照明光変調装置は有用である。また、後述の再構成精度を向上させる目的で照明光変調装置を使用してもよい。例えば、照明光変調装置として開口数の小さい絞りや複雑な透過率分布を有するマスクを用いると観察画像における標本の解像力は低下するが、再構成精度が照明光変調装置を用いない場合より向上するのであれば利用する価値がある。

照明光学系１００からの射出光は、標本ステージ２０１に載置された標本２０３に導かれる。ステージ駆動部２０２は、標本ステージ２０１を結像光学系３００の光軸方向と光軸に直交する方向とに移動する機能を備えており、さらには標本を交換するための機能を備えていてもよい。標本の交換は、別の機構（不図示）が自動的に行う形態でもよいし、ユーザが手動で行う形態でもよい。

標本２０３を同定する情報としては、透過光の振幅および位相分布が対応付けられる。この振幅および位相分布は、（通常は標本表面に垂直な方向からの）平行光を照明された標本２０３の直後における透過光の振幅と位相の２次元分布を意味し、以降では単に、標本２０３の「振幅」および「位相分布」と呼ぶ。

またこれを一般化して、標本２０３の構造を反映した複素量の２次元分布としてもよい。例えば、標本２０３において染色の濃度が高い部位や光の散乱が大きい部位は透過率が低いため、透過光の振幅は入射光に対して減衰する。また、標本２０３において屈折率が相対的に高い部位は光路長が相対的に長くなるために、入射光に対する位相変化量が相対的に大きい。なお「光路長」とは、光が透過する媒質の屈折率と厚みの積で与えられる量であり、透過光の位相変化量に比例する。このように標本２０３の透過光の振幅および位相分布は標本２０３の構造を反映しており、これらの分布から標本２０３の３次元的な構造を推定することも可能である。

標本２０３からの透過光は、結像光学系３００を介して撮像素子４０１上に標本２０３の光学像を形成する。ここで、結像光学系３００の光路中に配置された光学素子３０１は、結像光学系３００の瞳面近傍における投影光の強度または位相の少なくとも一方の分布を変調する手段である。光学素子３０１の目的は、再構成の対象である標本の振幅または位相分布の情報を効率的に観察像に埋め込むことである。言い換えると、取得する画像の枚数または解像度を最小限にして標本の振幅または位相分布の高精度な再構成を可能にすることである。光学素子３０１はＳＬＭのような可変型の素子であってもよいし、光学特性が固定された位相板のような素子であってもよい。また、光路中に挿脱可能な可動型の構造であってもよい。

なお、光学素子３０１の光学特性は製造誤差や制御誤差の影響を受け、再構成結果に影響することが懸念される。しかし、Ｓ２１３の説明で述べるように、あらかじめこの光学特性を計測するか、もしくは訓練標本の標本データが完全に既知であることにより、この問題は解決される。また、光学系の配置は必ずしも標本の透過光を結像する透過型である必要はなく、落射型の配置であってもよい。

撮像素子４０１は、結像光学系３００が投影する標本２０３の光学像を光電変換し、コンピュータ（ＰＣ）５０１、画像処理部５０２、記憶装置５０３のいずれかに画像データとして伝送する。

標本２０３に関する再構成を画像取得後に直ちに行わない場合には、撮像素子４０１からＰＣ５０１あるいは記憶装置５０３に画像データが伝送され、蓄積される。再構成を画像取得後に直ちに行う場合には、画像処理部５０２に画像データが伝送され、再構成するための演算処理を行う。

画像処理部５０２は、記憶手段、結合係数演算部（第１の算出部）、中間データ生成部（第２の算出部）、変換部（第３の算出部）、判断部を有する。これらは、別個のモジュールとしてハードウェアまたはソフトウェアで構成される。画像処理部５０２は、ＰＣ５０１によって制御されるが、画像処理装置のようにマイクロコンピュータ（プロセッサ）を含んで構成されてもよい。

結合係数演算部は、結像光学系３００が形成した複数の既知の標本の光学像を光電変換して得られる複数の画像から生成された基底の線形結合で、第２の画像を近似するための結合係数を算出する。中間データ生成部は、複数の既知の標本のデータの非線形写像によって得られる複数の複素量のデータと算出された結合係数と、に基づいて中間データを算出する。変換部は、算出された中間データから未知の標本に関する複素量のデータを算出する。判断部は、算出された結合係数に基づき、再構成に用いる訓練データを交換し、再構成をやり直すかどうかを判断する。

生成されたデータは、表示装置５０５に表示されるか、ＰＣ５０１あるいは記憶装置５０３に伝送され蓄積されるか、もしくは両方が実行される。この処理内容は、入力装置５０４を介したユーザの指示か、ＰＣ５０１あるいは記憶装置５０３に格納されている情報に基づき決定される。

図１における撮像装置１０以外の全ての装置は、必ずしも撮像装置１０と物理的に直接接続されている必要はない。例えば、撮像装置１０はＬＡＮ（Ｌｏｃａｌ Ａｒｅａ Ｎｅｔｗｏｒｋ）やクラウドサービスを介して外部に接続されていてもよい。このメリットとして、撮像装置１０と画像処理部５０２が一体化されていないために個々の撮像装置は低コストでコンパクトにできること、複数のユーザ間で実時間でデータを共有できることが挙げられる。

以下、本発明の様々な実施例における標本２０３の情報の再構成方法について説明する。

本発明は、未知の標本２０３の振幅および位相分布を、撮像装置を介して取得した評価画像から訓練データを用いて再構成するための手段を開示する。再構成方法の説明に先立ち、まず本実施形態の画像処理の概念について、図２を参照しながら説明する。

図２に示す画像処理は、独立の画像処理装置で行われてもよいし、撮像装置１０と一体の画像処理部５０２で行われてもよい。また、図２に示す再構成方法は画像処理方法として機能する。

まず、訓練データについて説明する。振幅および位相分布が既知の複数の標本（第１の標本）２０３を「訓練標本」と呼び、訓練標本の振幅および位相分布を再構成に用いる。本実施例では訓練標本の個数をＴと記す。Ｔ個の訓練標本の各々の振幅および位相分布を「標本データ」、未知の標本（第２の標本）の振幅および位相分布を「再構成データ」と呼ぶ。このＴ個の訓練標本の各々から撮像装置を用いて得られるＴ個の観察像（第１の画像）を「訓練画像」、未知の標本から同様に得られる１つの観察像（第２の画像）を「評価画像」と呼ぶ。また、この訓練標本と訓練画像を合わせて「訓練データ」と呼ぶ。

従って、第１の標本と第１の画像とは、部分コヒーレントまたはコヒーレントな結像光学系が形成した複数の既知の第１の標本の光学像を光電変換することによって得られる。具体的には、複数の第１の画像は、前記結像光学系または前記結像光学系と同一の光学特性を有する結像光学系が形成した複数の既知の第１の標本の光学像を光電変換することによって得られる。あるいは、複数の第１の画像は、複数の第１の標本の各々に対応する複素量のデータから計算によって生成されてもよい。

同様に、前記第２の画像は、前記結像光学系または前記結像光学系と同一の光学特性を有する結像光学系が形成した未知の第２の標本の光学像を光電変換することによって得られる。

本実施例では明視野顕微鏡を前提とするので、コヒーレント結像または部分的コヒーレント結像によって標本２０３の観察像が形成される。例えば、非特許文献３に開示されているように、この場合に標本２０３の振幅および位相分布と観察像の関係は非線形である。具体的には、観察像を表すベクトルＩと標本の振幅および位相分布を表すベクトルｘとの間には式（１）の関係が成り立つ。

ここで、ｘは標本２０３の振幅および位相分布を複素数で表す列ベクトル、Ｇはコヒーレント結像または部分的コヒーレント結像を表現する複素行列、ｖｅｃは行列を列ベクトルに分解して縦に連結する演算、Ｈはベクトルの転置共役である。また、ｋｒｏｎはクロネッカー積、右肩の＊は複素共役を表す。行列Ｇは結像光学系３００の回折限界に起因する光学的なぼけの情報に加えて、結像光学系３００の収差および焦点ずれ、光学素子３０１の情報、撮像装置自体または環境に起因する振動や温度変化などの画像劣化要因の全ての情報を含む。

なお、式（１）には明示していないが、現実には観測過程において撮像素子４０１などに起因して観測ノイズが加わるので、式（１）の右辺にノイズを表現する実定数ベクトルが加わると考えてよい。

次に、図２に示す入力空間、特徴空間、像空間の３つを定義する。

「入力空間」とは、標本２０３の振幅および位相分布に関するデータが形成する空間であり、個々のデータはＮ次元の複素ベクトルであるとする。標本データおよび再構成データは入力空間に含まれる。これらのデータは、標本の表面に略平行な平面内のＮ点の既定の座標でサンプリングされた振幅および位相の値を複素数として表現したものである。

「特徴空間」とは、入力空間のデータに対する非線形写像φによって得られるデータが形成する空間である。ここで、非線形写像φは、式（１）を踏まえ式（２）で定義する。

これは入力空間のＮ次元の複素ベクトルを特徴空間のＮ^２次元の複素ベクトルに変換する写像である。非線形写像φの導入により、式（１）で表されるコヒーレント結像または部分的コヒーレント結像は、特徴空間のデータに対する線形写像と解釈できる。この線形写像は式（１）に示すように行列Ｇを乗じる変換であるので、以降ではこの線形写像をＧと表す。

本実施形態の特徴の１つは、コヒーレント結像または部分的コヒーレント結像において逆問題を解く困難の原因となっていた非線形性を、非線形写像と線形写像とに分解する点にある。以降では、図２に示すように、標本データを特徴空間に写像したデータを「変換データ」、再構成データを特徴空間に写像したデータを「中間データ」と呼ぶ。

なお、図２に示すφ^−１はφの逆写像であり、φで写像されたデータをφ^−１で写像すると元に戻るものとする。しかし、式（２）のφに対してはφ^−１を陽に表現することはできず、φ^−１の結果を得るには実際には数値的に推定するしかない。そのための具体的な方法として、行列の特異値分解（ＳＶＤ；Ｓｉｎｇｕｌａｒ Ｖａｌｕｅ Ｄｅｃｏｍｐｏｓｉｔｉｏｎ）が挙げられる。本実施例において特徴空間の全てのデータは階数１の正方行列であるため、そのような行列に特異値分解を行うと、特異ベクトルとしてＮ次元のベクトルｘが一意に求められる。言い換えると、特異値分解を行い特異ベクトルを出力する操作としてφ^−１を定義することが可能である。

「像空間」とは、特徴空間のデータを線形写像Ｇで写像したデータ、すなわち観測像が形成する空間である。図２に示すように、変換データのＧによる写像が訓練画像であり、中間データのＧによる写像が評価画像という関係がある。像空間の個々のデータは現実に観測される像強度分布を既定のＭ点でサンプリングしたデータであり、Ｍ次元の実数ベクトルとする。

次に、カーネル行列を定義する。近年、機械学習の分野で非線形なモデルに基づく学習などの目的でカーネル法と呼ばれる手法が使われている。カーネル法の一般的な定義は、非特許文献４などで解説されている。一般的にカーネル法においては、カーネル行列Ｋを適当な非線形写像φ’を用いて式（３）で定義する。

ここで、Ｋ_ｉｊはカーネル行列Ｋのｉ行ｊ列成分、Ｘ_ｉはサンプルの母集団の中のｉ番目のサンプルに対応するデータ、＜・，・＞は内積をそれぞれ表す。内積を２つのデータの類似度と解釈するならば、カーネル行列Ｋとは特徴空間におけるデータの全ての組合せに対する類似度を表現した行列であると言うことができる。式（３）の非線形写像φ’を式（２）の非線形写像φに置き換えると、式（４）に示すように、カーネル行列Ｋはφを陽に含まない形で表すことができる。

ここで、ｘ_ｉはｉ番目の訓練標本の標本データを表す複素ベクトルである。カーネル行列Ｋを算出するために、式（３）では入力空間のデータを特徴空間に写像した後にＮ^２次元の２ベクトルの内積を計算するのに対し、式（４）では入力空間のＮ次元の２ベクトルの内積を計算するだけで済むため、式（３）に比べ計算量が大幅に減る。このようにカーネル行列の計算量を低減する手法は、一般的にカーネルトリックと呼ばれる。

本発明の特徴の１つは、カーネル法を機械学習ではなくコヒーレント結像または部分的コヒーレント結像の逆問題に応用する点にある。さらには、式（４）のカーネルを用いること、およびその結果として計算量を抑えられる点に特徴を有する。

式（３）においてφ’（Ｘ）に対して平均値を除去する処理（「センタリング」又は「中心化」と呼ばれる）を行う場合があるが、非特許文献４に開示されているように、その場合にもカーネル行列Ｋから直接センタリング後のカーネル行列Ｋを算出できる。

次に、カーネル行列Ｋに基づき訓練画像と評価画像の関係を抽出する。なおこの関係は、特徴空間における変換データと中間データの関係と同一である。その理由は、特徴空間と像空間が線形写像Ｇで対応しているためである。

関係を抽出するために、複数の訓練画像の線形結合で新たな基底を生成する。この新たな基底を固有像と呼ぶ。基底を生成する具体的な方法として、例えばカーネル行列Ｋに対し主成分分析を行い、得られた複数の固有ベクトルの各々を線形結合係数として用いて前記複数の訓練画像を線形結合することで、複数の固有像を生成する。これは式（５）で表すことができる。

ここで、Ｅは各列が固有像Ｅ_１、Ｅ_２、…Ｅ_ＬであるＭ行Ｌ列の行列、Ｉは各列が訓練画像Ｉ_１、Ｉ_２、…Ｉ_ＴであるＭ行Ｔ列の行列をそれぞれ表す。Ｌはカーネル行列Ｋの階数以下の自然数であり、ユーザが任意に指定してもよいし、カーネル行列Ｋの固有値に基づいて決定してもよい。後者の場合、例えば既定の閾値以上である固有値の数を自動的にＬとする決定方法が挙げられるが、これに限定されることはない。αはカーネル行列Ｋの複数の固有ベクトルからなるＴ行Ｌ列の行列である。行列αを求める具体的な方法として、例えばカーネル行列Ｋに対する特異値分解を用いることができる。その場合には、対応する特異値が大きい順、または他の基準によってＬ個の特異ベクトルを選択し、それらの列ベクトルを連結することで行列αを生成する。

次に、Ｌ個の固有像の線形結合で評価画像を近似することで、訓練画像と評価画像の関係を完全に決定する。これは、固有像を基底とするＬ次元空間において評価画像に最も近い解を探索すると言い換えることもできる。近似は、例えば、固有像の線形結合と評価画像との差分のノルムを最小化する結合係数を求める問題、すなわち最小二乗問題として定式化することができる。結合係数をＬ次元の実数ベクトルγとして表すことにすると、最小二乗問題は式（６）で表すことができる。

ここで、γ上のハット（＾）は推定解であることを表すために付けている。また、式（６）の右辺の第２項は解が異常値をとることを防止する目的で付加される正則化項の一種である。特に式（６）のように解のＬ２ノルム（最小二乗解）で定義される場合には、Ｔｉｋｈｏｎｏｖ正則化またはリッジ回帰と呼ばれる。なお、正則化項の係数λは任意の実数である。また、第１の算出部は評価画像にふくまれる観測ノイズの大きさを推定して推定値に基づき正則化係数λを決定してもよい。なお、正則化を行わないために正則化係数λを０としても良い。その場合、第１の算出部は、基底からなる行列に関するＭｏｏｒｅ―Ｐｅｎｒｏｓｅの疑似逆行列を用いて最小二乗解を求め、結合係数を算出できる。
正則化係数λが０でない場合、式（６）の解は、非特許文献４にも開示されているように式（７）で与えられ、反復処理なしに比較的高速に算出可能である。なお、式（７）においては、行列Ｅが縦長か横長かによって条件分岐される。

ここで、右肩のＴは行列の転置、右肩の−１は逆行列、行列ＬはＬ次元の単位行列、ベクトルＩ’は評価画像をそれぞれ表す。式（６）の解の算出方法は式（７）に限られることはなく、代わりに、例えば、式（８）を用いることができる。

ここで、Ｕ_ＬとＵ_ＲはＥの特異ベクトルからなる行列、ΣはＥの特異値を対角要素にもつ対角行列、ｔｈｒｅｓｈｏｌｄは引数の行列要素のうち閾値θを超える要素を０などの定数に置き換える関数である。このように式（５）および（６）に基づき、評価画像Ｉ’はＴ個の訓練画像に対応する行列Ｉに行列αとベクトルγを乗じることで近似される。これが訓練画像と評価画像の関係である。この関係を特徴空間にそのまま適用すると、中間データφ（ｚ）はＴ個の変換データを列成分とする行列Φに行列αとベクトルγを乗じることで得られる。この関係は式（９）で表すことができる。

ここで、行列Φは各列がφ（ｘ_１）、φ（ｘ_２）…φ（ｘ_Ｔ）であるＮ^２行Ｔ列の行列、Ｖは各列が訓練基底ｖ_１、ｖ_２、…ｖ_ＬであるＮ^２行Ｌ列の行列である。図２に示すように、訓練基底とは特徴空間において固有像に対応するＬ個のベクトルの集合であり、中間データを形成するための基底である。訓練基底は説明の都合上形式的に導入した概念であり、実施においては陽に計算する必要はない。

このように、線形写像Ｇの逆写像を直接計算することを回避していることが、本実施形態の特徴の１つである。像空間の次元Ｍよりも特徴空間の次元Ｎ^２の方が圧倒的に大きいために、線形写像Ｇの逆写像は一意に定まらず、この条件下での逆問題はいわゆる不良設定（ｉｌｌ−ｐｏｓｅｄ）問題に属する。

これに対し、本発明の方法ではこのような非適切性を含む計算を一切含まないため、安定して正しい再構成結果が得られる。また、前述したとおり、特徴空間のデータを入力空間に逆写像するには特異値分解などを用いればよく、その結果として中間データφ（ｚ）から再構成データｚが一意に求められる。以上のことを総合すると、式（５）、式（６）および式（９）の計算を順に行うことにより、既知のデータである訓練データと評価画像とから再構成データｚを再構成できる。

本発明の特徴の１つは、コヒーレント結像または部分的コヒーレント結像における再構成を、反復処理なしに行列積の演算だけで高速に実行する点にある。さらには、再構成の処理が図２において、評価画像Ｉ’から直接中間データφ（ｚ）ないしは再構成データｚに直接進むのではなく、線形写像Ｇの逆写像を回避して訓練データを経由する迂回したルートで進むという点に特徴を有する。

また、一般的なＣＳは再構成対象の情報が疎（ｓｐａｒｓｅ）であることを仮定するが、本発明は前述の説明から分かるように疎性の仮定を用いないため、原理的には疎でない情報を再構成できる点にも特徴を有する。

この点、非特許文献２は、物理的な観測量（観測像）どうしの関係をそのまま特徴空間に適用するという概念はなく、本発明とは本質的に異なる。

次に図３を参照して、図１に示す撮像システムにおいて標本２０３の情報を再構成するための画像処理方法について説明する。「Ｓ」はステップを表し、図３に示すフローチャートはコンピュータに各ステップの機能を実現させるためのプログラムとして具現化が可能である。まず全体の処理の手順について、図３（ａ）に示すフローチャートを参照して説明する。

Ｓ２０１では、標本ステージ２０１上に標本２０３を設置する。例えば、標本ステージ２０１自身、または連携する自動送り込み機構が、カセットなどの標本保持手段から標本２０３を取り出して標本ステージ２０１上に載置する。なお、この処理を装置が自動で行う代わりにユーザ自身が手動で行ってもよい。ステージ駆動部２０２による標本ステージ２０１の駆動を制御することによって、標本２０３の観察対象領域を所定の合焦状態で撮像するために適切な位置に標本を移動する。この移動は、Ｓ２０３以前ならばどのタイミングで行ってもよい。

Ｓ２０２では、標本２０３の情報の再構成を行う場合には、必要に応じて光学素子３０１が投影光の強度または位相の少なくとも一方の分布を変調する状態にする。一方、通常の顕微鏡像を取得する場合には、光学素子３０１を光路の外に移動する、またはＳＬＭを制御するなどの方法で投影光の強度と位相に有意な変調を与えない状態にする。

Ｓ２０３では、振幅および位相分布が未知である標本２０３の画像を撮像する。

Ｓ２０４では、Ｓ２０３において取得された画像データを撮像素子４０１からＰＣ５０１、画像処理部５０２、記憶装置５０３のいずれかに伝送する。

標本２０３の情報の再構成を行う場合には、Ｓ２０６の処理を実行する。

Ｓ２０６では、画像処理部５０２において前記画像データに基づき標本の振幅または位相分布に関する再構成データを出力する。この処理内容は、後に図３（ｂ）を用いて詳細に説明する。

Ｓ２０７では、ユーザからの指示または既定の設定に従い、前記再構成データを記憶装置５０３またはＰＣ５０１に格納するか、表示装置５０５に表示する。

次に、Ｓ２０６において画像処理部５０２が行う画像処理について、図３（ｂ）に示すフローチャートを参照しながら説明する。

Ｓ２１１では、ＰＣ５０１または記憶装置５０３のいずれかに蓄積されている訓練データと未知の標本に対する評価画像とを読み出す。この標本と評価画像の組は１つでもよいし、複数でもよい。

訓練データは、図３（ａ）に示した一連の処理に先立ち予め取得し蓄積されている。訓練データは、撮像装置１０における標本２０３と観察像の関係に影響する要因、例えば、結像光学系３００の収差情報や焦点ずれ量、光学素子３０１が変調する強度及び位相分布についての情報が既知である場合には全て計算によって生成してもよい。即ち、既定のルールに従いＴ個の標本データを計算上で生成し、標本データと撮像装置１０の情報とに基づき結像シミュレーションによって訓練画像を計算上で生成する。

撮像装置１０の情報は、訓練データを生成する前に撮像装置１０に対して計測を行い取得してもよい。例えば、撮像装置１０に対して一般的な波面収差計測法を用いて結像シミュレーションに用いるための収差データを取得する。

一方で、撮像装置１０の情報が未知のまま再構成を行うことも可能である。この場合には、訓練標本は実在する複数の標本２０３とし、これらの標本２０３に対して撮像装置１０を用いて図３（ａ）と同じ条件および手順で訓練画像を取得する。

この場合には、全ての訓練標本に関する標本データ、すなわち透過光の振幅および位相分布は既知でなければならない。この標本データは、一般的な波面収差計測法や表面形状計測手段によって得られたデータから生成してもよいし、人工的に作成した標本であればその設計値から生成してもよい。また、訓練標本は見かけ上必ずしも複数である必要はない。例えば、訓練標本として有効な要素が１つの標本の中に複数集積してもよいし、このような形態に限られることはない。

このように、撮像装置１０の情報が未知のままで逆問題を解くことができる、言い換えるとブラインド推定が可能である点に、本発明の特徴の１つがある。ブラインド推定のメリットは、種々の画像劣化要因に再構成精度が影響を受けない頑健性である。

画像劣化要因の一例としては、撮像装置１０や光学素子３０１の製造誤差に起因する性能のばらつき、撮像装置自体または環境に起因する振動や温度変化などが挙げられる。ブラインド推定が可能である理由は、標本２０３と観測像の関係を含む訓練データを用いることで、式（１）において撮像装置１０の情報を全て含む行列Ｇを直接再構成に用いることを回避しているためである。

もし、既に訓練データを用いて固有像Ｅ、変換データΦ、カーネル行列αが求められており、それらがＰＣ５０１または記憶装置５０３のいずれかに蓄積されているならば、以降のＳ２１３〜Ｓ２１４は行う必要がない。言い換えると、一度訓練データから定数行列（固有像および式（９）の行列Ｖ）を生成しておけば、撮像装置１０に変化がない限り未知の標本２０３および評価画像がいくら変化しても、Ｓ２１５以降の計算を反復なしに１回行うだけで再構成の処理が終わる。

Ｓ２１３では、標本データを用いて式（４）または（３）に基づきカーネル行列αを生成する。Ｓ２１４では、式（５）に基づき固有像Ｅを生成する。

Ｓ２１５（第１のステップ、第１の機能）では、結合係数演算部（第１の算出部）が式（７）または（８）などに基づき線形結合係数γを求める。

Ｓ２１７（第２のステップ、第２の機能）では、中間データ生成部（第２の算出部）が式（９）に基づき中間データφ（ｚ）を求める。Ｓ２１８（第３のステップ、第３の機能）では、変換部（第３の算出部）が中間データφ（ｚ）の逆写像であるｚを求める。

しかし、未知の標本２０３と訓練データとの組み合わせによっては再構成精度が大きく低下することもある。そのような場合には、γのノルムが異常に大きい値をとる。これを解決する１つの方法として、Ｓ２１５において求められた線形結合係数γのノルムが閾値を超える場合には、再構成に用いる訓練データを変えるという方法が挙げられる。この場合、判断部は結合係数のノルムが既定の閾値以下であるかどうかを判断する（Ｓ２１６）。この処理方法を、図３（ｃ）のフローチャートに示す。

Ｓ２１５において求められた線形結合係数γのノルムが閾値を超える場合には（Ｓ２１６のＮＯ）、Ｓ２１９に進み訓練データを交換する。具体的には、再構成に用いる訓練データの組数Ｔよりも多く訓練データをあらかじめ用意しておき、そのうちのＴ組だけを選択して再構成に用いてもよい。結像シミュレーションによって訓練画像を計算上で生成する場合には、この代わりに標本データを既定の規則に従い新たに生成し直して訓練画像を生成してもよい。また、訓練データを交換する方法はこれらに限られることはない。

本発明の特徴の１つは、再構成の途中の時点で再構成誤差が予測可能であり、なおかつ誤差が大きいと予測される場合には誤差を低減する手段（訓練データの交換）を自動的に実行できる点にある。これについては以下の実施例３において、具体例に基づき詳細を説明する。

次に、本発明の実施例１を、数値シミュレーションを用いて説明する。

以下では、シミュレーション条件を説明する。照明光学系１００から標本に照射される照明光の波長を０．５５μｍ、結像光学系３００の標本側の開口数を０．７０、照明光学系１００の開口数を０．４９（すなわちコヒーレンスファクタ０．７０）とする。

図４（ａ）に示すように、照明光学系の瞳面の透過光強度分布（または標本がない場合に結像光学系の瞳面に形成される光強度分布）は開口数０．４９に対応する円形境界の内部で一様である。

図４（ｂ）および（ｃ）に、結像光学系の瞳面に配置された光学素子による透過光の振幅と位相の変化の分布を示す。振幅分布はサンプリング点ごとに独立に０〜１の値を一様乱数で生成し、位相分布はサンプリング点ごとに独立に標準偏差２πラジアンのガウス乱数で生成した。また、結像倍率は説明の都合上１倍とするため、結像光学系３００の像側の開口数は０．７０であり、標本と像面の縮尺は同一とする。なお、現実の顕微鏡では数十倍〜数百倍の結像倍率が用いられるが、以降の議論の本質は何ら変わることがない。明視野顕微鏡は部分コヒーレント結像に従うことが知られているため、本実施例では一般的な２次元の部分コヒーレント結像の理論に基づくシミュレーションを行う。また、乾式顕微鏡を想定し、標本と結像光学系３００の間は屈折率１．００の空気で満たされているとする。

以上の条件は、以降の訓練標本および未知の標本の撮像において一切変化しないものとする。また、以降の全ての標本と像のサンプリング間隔は０．３０μｍとし、振幅は０から１の間の実数とし、位相はラジアン単位で示す。

図５（ａ）および（ｂ）に、１６０個の１１×１１画素の振幅および位相分布からなる標本データを縦に８個、横に２０個に並べて示す。各々の標本データは、透過率、位相、位置をランダムに決定した２個の開口の振幅および位相分布をベクトル化し、２値のランダム行列を乗じることで生成した見かけ上密な振幅および位相分布である。

図５（ｃ）に、この１６０個の標本データから前述の条件で計算上得られる１６０個の訓練画像（像強度分布）を、同様に並べて示す。また、この標本データから式（４）に従い１６０行１６０列のカーネル行列を生成し、対応する固有値が大きい順に１２０個の固有ベクトルを抽出し１６０行１２０列の行列αを算出する。

この行列αと図５（ｃ）の訓練画像とから式（５）に従い１２０個の固有像が生成される。図５（ｄ）に固有像を、見やすさのために絶対値の常用対数をとった上で縦に６個、横に２０個に並べて示す。図５（ｄ）の左方の５３個以外は１．００Ｅ−１０以下の輝度値であり、これらを用いなくても再構成精度に有意な変化は生じない。このことは、カーネル行列の固有値または固有ベクトルから必要な固有像の数Ｌが決定できることを示唆している。本実施例の場合には、Ｌは５３以上であれば十分であるといえるが、以降では固有像の数Ｌを１２０としたまま計算を行う。

一方、未知の標本２０３として、図６（ａ）および（ｂ）に示す１１×１１画素の振幅および位相分布を設定する。図６（ｃ）に、この未知の標本２０３から前記条件の明視野顕微鏡を介して得られる評価画像（像強度分布）のシミュレーション結果を示す。ここで、像強度分布が標本２０３とは全く異なる様相の分布をしているのは、図４に示した光学素子を用いているためである。

次に、図５に示した訓練データを用いて、図６（ｃ）の評価画像から未知の標本の振幅および位相分布を再構成する。まず、式（７）において正則化パラメータλを０として線形結合係数γを算出する。こうして得られたγを式（９）に代入し、得られた中間データφ（ｚ）を１２１行１２１列の行列に変形した後に特異値分解を行う。こうして得られた第１の特異値の平方根と第１の左特異ベクトルとの積を、１１行１１列の行列に変形することで、図７（ａ）および（ｂ）に示す再構成された未知の標本の振幅および位相分布を得る。これらは、図６（ａ）および（ｂ）の真の分布をほぼ完全に再構成している。この再構成の精度を定量化するために、式（１０）で定義されるＲＭＳＥ（Ｒｏｏｔ Ｍｅａｎ Ｓｑｕａｒｅ Ｅｒｒｏｒ）を用いる。

ここで、Ｎは画素数（本実施例においては１２１）、ｉは画素番号、ｘ_ｉは再構成された画素ｉの振幅または位相、ｘ_ｉ’は画素ｉの真の振幅または位相をそれぞれ表す。

図７（ａ）の図６（ａ）に対するＲＭＳＥは４．２９Ｅ−１２、図７（ｂ）の図６（ｂ）に対するＲＭＳＥは３．９８Ｅ−１１ラジアンであり、無視できる水準の誤差である。

このように、明視野顕微鏡を用いて取得した１枚の評価画像だけから標本２０３の振幅および位相分布を高精度で再構成することが、本実施例の方法により可能になる。

このようにして再構成された振幅および位相分布は、未知の標本２０３の３次元的な構造を把握する目的に利用することも可能である。簡単な例としては、位相分布に既定の定数を乗じることにより、屈折率が略一様な標本に対しては厚さ分布を推定することができる。

また、この振幅および位相分布を利用すれば、標本中の特定の構造を強調して描画するなどの従来にないレンダリングが可能になり、標本２０３の情報の見せ方の自由度が大きく拡張される。

さらには、再構成された分布は原理的に顕微鏡の画質劣化要因の影響が除去されているため、通常の明視野顕微鏡を用いて観察される画像に比べ解像力が向上した鮮明な画像が得られ、微細な構造を観察することが容易になる。なお、画質劣化要因とは具体的には、結像光学系の回折限界に起因するぼけ、撮像素子に起因するノイズや解像力の低下などである。

次に、本発明の実施例２を、数値シミュレーションを用いて説明する。観測ノイズ以外の全ての条件とデータは実施例１と同じとする。

図６（ｃ）に示す評価画像に観測ノイズとして加法的白色ガウスノイズを与える。各画素のノイズは互いに独立でありながら共通の統計分布に従い、その分布は平均が０、標準偏差が輝度の最大値の１．００％の正規分布である。本実施例では実施例１と異なり、図８に示す観測ノイズを加えた評価画像を用いて式（７）に基づき再構成を行う。

図９（ａ）および（ｂ）に、式（７）の正則化パラメータλを実施例１と同様に０とした場合の振幅および位相分布の再構成結果を示す。図９（ａ）の図６（ａ）に対するＲＭＳＥは１．０７Ｅ−１、図９（ｂ）の図６（ｂ）に対するＲＭＳＥは１．９２ラジアンである。一方、図９（ｃ）および（ｄ）に、式（７）の正則化パラメータλを１．００Ｅ−６とした場合の振幅および位相分布の再構成結果を示す。

図９（ｃ）の図６（ａ）に対するＲＭＳＥは７．８７Ｅ−３、図９（ｄ）の図６（ｂ）に対するＲＭＳＥは８．１８Ｅ−１ラジアンである。図とＲＭＳＥ値から分かるように、正則化を加えた方が明らかに真の分布に近い。

以上の結果が示すように、観測ノイズの影響を抑制するためには、式（６）の最小二乗問題における正則化が有効である。なお、本実施例では評価画像に観測ノイズを加えた場合だけを説明したが、訓練画像に同様に観測ノイズが加わっている場合にも同様に再構成ができる。

次に、本発明の実施例３を、数値シミュレーションを用いて説明する。図１０に示す未知の標本以外は全ての条件とデータは実施例１と同じとし、観測ノイズはないものとする。

図１０（ａ）および（ｂ）に、未知の標本の振幅および位相分布を、図１０（ｃ）に対応する評価画像を示す。この未知の標本は図５の訓練データとの相性が良くないために、その再構成結果は図１１（ａ）および（ｂ）に示す振幅および位相分布となり、誤差が比較的大きい。図１１（ａ）の図１０（ａ）に対するＲＭＳＥは５．７６Ｅ−２、図１１（ｂ）の図１０（ｂ）に対するＲＭＳＥは１．４７ラジアンである。

そこで、図３（ｃ）に示すフローに従い再構成を行う。すなわち、γのＬ２ノルムが閾値を超える場合には再構成を中断して訓練データを交換して再構成をやり直す。γのＬ２ノルムに関する閾値は１．００Ｅ＋６とする。即ち、判断部は、判定条件を満たさない場合に、複数の第１の標本を他と交換し、第１の算出部および第２の算出部に中間データの算出をやり直させる。この場合、判定条件は、前記結合係数のノルムが既定の閾値以下であるという条件である。

図１２に、γのＬ２ノルムが閾値以下になった時の訓練データおよび再構成された振幅および位相分布を示す。図１１の結果に伴うγのＬ２ノルムは６．３３Ｅ＋１４であったのに対し、図１２の結果に対しては閾値未満の１．０８Ｅ＋４であった。図１２において、（ａ）は訓練標本の振幅分布、（ｂ）は訓練標本の位相分布、（ｃ）は訓練画像、（ｄ）は固有像を、図５と同様の形式で示している。訓練データが実施例１と異なることに起因し、図１２（ｄ）の有意な固有像は右端の列を除く１１４個になっている。

また、図１３（ａ）および（ｂ）に、振幅および位相分布の再構成結果を示す。図１３（ａ）の図１０（ａ）に対するＲＭＳＥは２．６７Ｅ−１２、図１３（ｂ）の図１０（ｂ）に対するＲＭＳＥは４．４１Ｅ−１１ラジアンであり、実施例１と同等の精度である。以上の結果が示すように、再構成を確実に成功するために図３（ｃ）に示すフローが有効であり、再構成精度は再構成の途中で得られるγのＬ２ノルムの値から予測可能である。

本発明は、デジタル顕微鏡、画像処理装置、画像処理方法の分野に適用することができる。

１０…撮像装置、３００…結像光学系、５０１…コンピュータ（ＰＣ）、５０２…画像処理部

Claims (18)

1. 部分コヒーレントまたはコヒーレントな結像光学系が形成した既知の複数の第１の標本の光学像を光電変換して得られる複数の第１の画像から生成された基底の線形結合により、第２の画像を近似するための結合係数を算出する第１の算出部と、
   前記複数の第１の標本のデータの非線形写像によって得られる複数の複素量のデータと、前記第１の算出部によって算出された前記結合係数と、に基づいて中間データを算出する第２の算出部と、
   前記第２の算出部によって算出された前記中間データに基づいて未知の第２の標本に関する複素量のデータを算出する第３の算出部と、
   を有することを特徴とする画像処理装置。
2. 前記第２の画像は、前記結像光学系または前記結像光学系と同一の光学特性を有する結像光学系が形成した前記第２の標本の光学像を光電変換して得られ、前記基底は、前記第１の算出部により前記複数の第１の標本の画像の線形結合として生成されることを特徴とする請求項１に記載の画像処理装置。
3. 前記第１の算出部は、前記複数の第１の標本の各々に対応する複素量のデータどうしの内積を用いて定義されるカーネル行列の固有ベクトルを線形結合係数として、前記複数の第１の標本の画像を線形結合することを特徴とする請求項２に記載の画像処理装置。
4. 前記第１の算出部は、前記カーネル行列の固有値に基づき前記基底の数を決定することを特徴とする請求項３に記載の画像処理装置。
5. 前記第１の算出部は、前記複数の第１の画像をそのまま前記基底として出力することを特徴とする請求項２に記載の画像処理装置。
6. 前記中間データは、前記複数の第１の標本の各々に対応する複素量の行列とその複素共役とのクロネッカー積から得られる複数の行列の線形結合であることを特徴とする請求項１乃至５のうちいずれか１項に記載の画像処理装置。
7. 前記複数の第１の画像は、前記結像光学系または前記結像光学系と同一の光学特性を有する結像光学系が形成した複数の既知の第１の標本の光学像を光電変換して得られることを特徴とする請求項１乃至６のいずれか１項に記載の画像処理装置。
8. 前記複数の第１の画像は、前記複数の第１の標本の各々に対応する複素量のデータから計算によって生成されることを特徴とする請求項１乃至６のうちいずれか１項に記載の画像処理装置。
9. 判定条件を満たさない場合に、前記複数の第１の標本を他と交換し、前記第１の算出部および前記第２の算出部に中間データの算出をやり直させる判断部を更に有することを特徴とする請求項１乃至８のうちいずれか１項に記載の画像処理装置。
10. 前記判定条件は、前記結合係数のノルムが既定の閾値以下であるという条件であることを特徴とする請求項９に記載の画像処理装置。
11. 前記第１の算出部は、前記結合係数を算出するために前記基底からなる行列に関するＭｏｏｒｅ―Ｐｅｎｒｏｓｅの疑似逆行列を用いて最小二乗解を求めることを特徴とする請求項１乃至１０のいずれか１項に記載の画像処理装置。
12. 前記第１の算出部は、前記結合係数を算出するためにＴｉｋｈｏｎｏｖ正則化を用いて最小二乗解を求めることを特徴とする請求項１乃至１０項のいずれか１項に記載の画像処理装置。
13. 前記第１の算出部は、前記第２の画像に含まれるノイズの大きさの推定値に基づき正則化係数を決定することを特徴とする請求項１２に記載の画像処理装置。
14. 前記第３の算出部は、前記中間データからなる行列に対して特異値分解を行うことを特徴とする請求項１乃至１３のいずれか１項に記載の画像処理装置。
15. 標本の光学像を形成する部分コヒーレントまたはコヒーレントな結像光学系と、
    前記結像光学系が形成した前記標本の前記光学像を光電変換する撮像素子と、
    前記結像光学系が形成した複数の既知の第１の標本の光学像を光電変換することによって得られた複数の第１の画像から生成された基底の線形結合で、前記結像光学系が形成した未知の第２の標本の光学像を光電変換して得られた第２の画像を近似するための結合係数を算出し、前記複数の第１の標本のデータの非線形写像によって得られる複数の複素量のデータと、算出された前記結合係数と、に基づいて中間データを算出し、算出された前記中間データに基づいて前記第２の標本に関する複素量のデータを算出する画像処理部と、
    を有することを特徴とする撮像装置。
16. 前記結像光学系は、光の強度または位相の少なくとも一方の分布を変調する光学素子を瞳面近傍に備えることを特徴とする請求項１５に記載の撮像装置。
17. 部分コヒーレントまたはコヒーレントな結像光学系が形成した複数の既知の第１の標本の光学像を光電変換して得られる複数の第１の画像から生成された基底の線形結合で、第２の画像を近似するための結合係数を算出する第１のステップと、
    前記複数の第１の標本のデータの非線形写像によって得られる複数の複素量のデータと、前記第１の算出部によって算出された前記結合係数と、に基づいて中間データを算出する第２のステップと、
    前記第２の算出部によって算出された前記中間データに基づいて未知の第２の標本に関する複素量のデータを算出する第３のステップと、
    を有することを特徴とする画像処理方法。
18. コンピュータに、
    部分コヒーレントまたはコヒーレントな結像光学系が形成した複数の既知の第１の標本の光学像を光電変換して得られる複数の第１の画像から生成された基底の線形結合で、第２の画像を近似するための結合係数を算出する第１の機能、
    前記複数の第１の標本のデータの非線形写像によって得られる複数の複素量のデータと、前記第１の算出部によって算出された前記結合係数と、に基づいて中間データを算出する第２の機能、
    前記第２の算出部によって算出された前記中間データに基づいて未知の第２の標本に関する複素量のデータを算出する第３の機能、
    を実現させるためのプログラム。