JP2005197792A - 画像処理方法、画像処理装置、プログラム、記憶媒体及び画像処理システム - Google Patents

Abstract

【課題】 時間的に連続したデジタル画像の間の経時変化を検出した画像を強調することにより診断精度を高める。
【解決手段】 第１、第２の画像Ｆ１、Ｆ２を入力する(Ｓ１１０)。続いて、第１、第２の画像Ｆ１、Ｆ２を位置合わせする(Ｓ１２０)。この位置合わせは２つの画像Ｆ１、Ｆ２の間に存在する位置ずれを検出し、対応画素位置を解剖学的特徴を用いて、一方の画像を他方の画像に合わせるようにワーピング処理を施す。次に、位置合わせ完了画像Ｘ１、Ｘ２を用いて、独立成分分析処理を行う(Ｓ１３０)。主要成分判断では、独立成分分析処理における２つの出力画像のうち、何れが主要成分画像であり、何れが経時変化成分画像であるかを判断する(Ｓ１４０)。一方が主要成分画像として判断されれば、他方が経時変化成分画像となり、このようにして主要成分画像と経時変化検出画像を出力する(Ｓ１５０)。
【選択図】図５

Description

本発明は、ＣＲやＭＲＩなどの医用診断機器により得られた時間的に異なる同じ部位の複数の画像を位置合わせして差分処理を行う画像処理方法及び処理装置に関するものである。

近年、医用画像診断の分野におけるデジタル画像の利用が進んでいる。例えば、半導体センサを使用してＸ線画像を撮影する装置は、従来の銀塩写真を用いる放射線写真システムと比較して極めて広い放射線露出域に渡って画像を記録し、更に画像の保存や伝送において効率的なシステムを構築し易いという実用的な利点を有している。

また、医用画像をデジタル化することで、従来の銀塩写真では困難であった診断形態の可能性が生まれている。即ち、従来は患者の経過観察などで異なる時点で撮影されたＸ線画像を比較する際には、フィルムをシャウカステンに架けて比較読影することが一般的に行われている。この読影では、相違点を見付けることが重要であるが、２枚のフィルムの相違点を探し出すには、細部を注意深く観察できる観察力と豊富な経験とを有する必要がある。

ここでデジタル画像を用いれば、異なる時点で撮影された２枚のデジタル画像を正常な解剖学的構造が一致するように位置合わせして、差分処理を行うことにより差分画像を生成して出力し、この差分画像を基となった１対の画像と比較読影することにより、画像間の変化をより正確に把握することが可能となる。生成された差分画像には相違点だけが表示されているため、フィルムによる読影により容易に比較読影できるためである。これは、経時差分ＣＡＤ(Computed Aided Diagnosis)と呼ばれている読影を支援する技術である。

このような差分画像を生成する処理方法については、例えば特許文献１に開示されており、異なる時点で撮影された２枚の胸部Ｘ線画像を位置合わせし、差分画像を生成することができる。

米国特許５３５９５１３号

図１６はこのような差分画像の生成表示装置の構成図である。第１及び第２画像は異なる時点で撮影された特定部位の医用画像データであるが、異なる時点で撮影されたため両者の撮影条件は完全に一致することはない。例えばＸ線画像の場合に、Ｘ線の線量、管電圧、管電流などの撮影条件パラメータが一致せず、両者の画像の濃度は異なっており、また撮影被写体の位置や姿勢も異なっている。

被写体を撮影した２枚の画像データは、濃度補正部１において両者の画像信号の濃度値分布がほぼ等しくなるように補正され、位置合わせ部２において解剖学的な構造の位置関係を用いて補正され、差分演算部３において対応する画素間で差分処理が行われて差分画像が生成される。そして、差分画像は表示部４において第１又は第２画像と共に表示される。

しかしながら経時差分ＣＡＤは、経時変化を検出する医用画像処理方法においては、完全に満足できるものではない。差分処理の精度は、濃度補正部１において、両者の画像信号の濃度値分布を等しく補正する補正処理の精度に左右される。

例として、被写体を同じ位置から撮影した２枚の画像について考えてみると、第１画像は時刻ｔ１で撮影し、第２画像は時刻ｔ２で撮影したとする。簡単のため被写体の姿勢を同じとすると、第１画像と第２画像は共通部分や似ている部分が非常に多いが、細部と濃度とが必ず異なる。例えば、被写体が患者の胸部で、撮影方法がＸ線撮影であれば、第１画像と第２画像は横隔膜の位置などが異なり、またＸ線照射量が異なるため撮影画像の濃度が異なる。Ｘ線画像でない場合でも、被写体が人間の場合には第１画像と第２画像を全く同じに撮影することはできないと考えてよい。

第１画像の濃度をｘ１(ｘ，ｙ)、第２画像の濃度をｘ２(ｘ，ｙ)、経時変化分の濃度をΔ(ｘ，ｙ)とすると、これらは次の関係にある。
ｘ１(ｘ，ｙ)＝ｓ(ｘ，ｙ)
ｘ２(ｘ，ｙ)＝ａ・ｓ(ｘ，ｙ)＋Δ(ｘ，ｙ)

ｓ(ｘ，ｙ)は第１画像と第２画像と共通する構造物の濃度であり、ａは定数である。第１、第２画像のそれぞれの撮影条件は異なるため、定数ａは１とは限らない。濃度補正部１において第２画像の濃度が第１画像の濃度に等しくなるように定数ａを推定し、その推定結果をａ’とすると、差分処理結果ｈ(ｘ，ｙ)は次式となる。
ｈ(ｘ，ｙ)＝ｘ２(ｘ，ｙ)−ａ’・ｘ１(ｘ，ｙ)
＝ａ・ｓ(ｘ，ｙ)＋Δ(ｘ，ｙ)−ａ’・ｓ(ｘ，ｙ)
＝(ａ−ａ’)・ｓ(ｘ，ｙ)＋Δ(ｘ，ｙ)

定数ａを正しく推定している場合は、ａ’＝ａであるから、
ｈ(ｘ，ｙ)＝Δ(ｘ，ｙ)
であるが、定数ａの推定結果が正しくないと、
｜(ａ−ａ’)・ｓ(ｘ，ｙ)｜＞｜Δ(ｘ，ｙ)｜
となり、経時変化Δ(ｘ，ｙ)を検出できず、経時変化Δ(ｘ，ｙ)が左辺の(ａ−ａ’)・ｓ(ｘ，ｙ)に埋もれてしまう。また、細部構造に違いがない場合は、経時変化Δ(ｘ，ｙ)＝０であり、差分処理結果ｈ(ｘ，ｙ)は、０を出力するはずである。しかし、定数ａの推定結果が正しくないと次式となる。
ｈ(ｘ，ｙ)＝(ａ−ａ’)・ｓ(ｘ，ｙ)

つまり、被写体に病変などの変化がない場合であっても、医師は経時変化ありと判断してしまう。この判断は医師の診断ミスを誘発することになる。

本発明の目的は、上述の問題点を解消し、時系列的に異なる２つの画像間の経時変化を精度良く検出し得る画像処理方法及び処理装置を提供することにある。

上記目的を達成するための請求項１に係る本発明は、２つの画像を入力する画像入力工程と、前記２つの画像の位置合わせをする位置合わせ工程と、該位置合わせ工程で位置合わせをした画像を用いて独立成分分析処理を行う独立成分分析処理工程とを有することを特徴とする画像処理方法である。

また、請求項２に係る本発明は、前記独立成分分析処理工程において出力された処理結果を、主要成分画像と経時変化成分画像とに区別するための判断を行う主要成分判断工程を有することを特徴とする請求項１に記載の画像処理方法である。

請求項３に係る本発明は、前記独立成分分析処理工程は、前記位置合わせ工程で位置合わせを行った画像を二次元データ形式から一次元データ形式に変換する一次元化処理工程と、該一次元化処理工程で得られた一次元データから直流成分を減算する分析前処理工程と、該分析前処理工程で得られた一次元データから分離行列、独立成分係数行列及び独立成分を算出する分離行列算出工程と、該分離行列算出工程で算出した独立成分に前記直流成分を加算する分析後処理工程と、分析後処理工程で直流成分を加算した独立成分を二次元データ形式に変換する二次元化処理工程とを備えることを特徴とする請求項１又は２に記載の画像処理方法である。

請求項４に係る本発明は、前記主要成分判断工程は、前記独立成分分析処理から出力された処理結果と、前記独立成分分析処理工程で用いた前記位置合わせした画像と他方の画像との２つの画像のうちの一方との類似性を調べて前記主要成分画像を判断することを特徴とする請求項２又は３に記載の画像処理方法である。

請求項５に係る本発明は、前記主要成分判断工程は、前記独立成分分析処理から出力された処理結果をヒストグラム処理して、該ヒストグラムの分布幅が小さい方を前記経時変化成分画像と判断することを特徴とする請求項２又は３に記載の画像処理方法である。

請求項６に係る本発明は、前記主要成分判断工程で判断した前記経時変化成分画像をフィルタ処理により強調した画像を、前記画像入力工程に入力された画像又は前記独立成分分析処理工程で用いた前記位置合わせした画像又は前記主要成分判断工程で判断した前記主要成分画像と重ね合わせることにより、経時変化部分を強調した画像を生成する経時変化強調画像生成工程を有することを特徴とする請求項２〜６の何れか１つの請求項に記載の画像処理方法である。

請求項７に係る本発明は、前記２つの画像は異なる時点で撮影された同一部位の画像であることを特徴とする請求項１〜６の何れか１つの請求項に記載の画像処理方法である。

請求項８に係る本発明は、前記位置合わせ工程は前記２つの画像の対応画素位置を解剖学的特徴を用いて検出し、一方の画像を他方の画像に合わせるようにワーピング処理を施す工程であることを特徴とする請求項１〜７の何れか１つの請求項に記載の画像処理方法である。

請求項９に係る本発明は、２つの画像から独立成分分析処理を行う画像処理方法をコンピュータに実行させるためのプログラムであって、２つの画像を入力する画像入力工程と、前記２つの画像の位置合わせをする位置合わせ工程と、該位置合わせ工程で位置合わせをした画像を用いて独立成分分析処理を行う独立成分分析処理工程をコンピュータに実行させるためのプログラムである。

請求項１０に係る本発明は、２つの画像から独立成分分析処理を行う画像処理方法をコンピュータに実行させるためのプログラムを記録したコンピュータ読み取り可能なコンピュータ可読記憶媒体であって、前記２つの画像を入力する画像入力工程と、前記２つの画像の位置合わせをする位置合わせ工程と、該位置合わせ工程で位置合わせをした画像を用いて前記独立成分分析処理を行う独立成分分析処理工程をコンピュータに実行させるためのプログラムを記録したコンピュータ読み取り可能なコンピュータ可読記憶媒体である。

請求項１１に係る本発明は、２つの画像を入力する画像入力手段と、前記２つの画像の位置合わせをする位置合わせ手段と、該位置合わせ手段で位置合わせをした画像を用いて独立成分分析処理を行う独立成分分析処理手段とを有することを特徴とする画像処理装置である。

請求項１２に係る本発明は、請求項１１に記載の画像処理装置と、画像生成装置、表示装置、ファイルサーバのうちの少なくとも１つとがネットワークを介して接続されて成ることを特徴とする画像処理システムである。

本発明に係る画像処理方法及び処理装置によれば、独立成分分析処理により、経時変化を検出する際に重要となる、画像の濃度補正を高精度で行うことができ、経時変化を検出した画像を精度良く生成できる。

経時変化が検出された画像は、腫瘍、炎症等の病状の変化を示してので、医師や技師は精度良く生成した検出画像を見ることにより、経時変化を容易に認識できる。診断行為が長時間に及んで診断に費やす集中力が低下すると、病状変化を見逃す可能性が高くなるが、この可能性を低くする効果もあり、診断結果のばらつき程度を抑えることになり、これらの相乗効果により診断精度が向上する。

本発明を図示の実施例に基づいて詳細に説明する。

図１は本発明による画像処理装置の構成図であり、画像生成部１０から出力される２つの画像Ｆ１、Ｆ２は経時変化検出部２０を介して経時変化表示部３０に送信されるようになっている。画像生成部１０においては、センサ１１１を有する撮像部１１で生成された画像出力は、記憶部１２、選択部１３を経て経時変化検出部２０に接続されている。

経時変化検出部２０においては、画像生成部１０の２つの画像Ｆ１、Ｆ２を受信する位置合わせ部２１が設けられており、同時に２つの画像Ｆ１、Ｆ２は経時変化表示部３０に直接に送信されるようになっている。位置合わせ部２１で生成された位置合わせ画像Ｘ１、Ｘ２は独立成分分析処理部２２に送信されると共に、経時変化表示部３０に送信されている。独立成分分析処理部２２で得られた独立成分画像Ｙ１、Ｙ２、独立成分係数行列Ａは、主要成分判断部２３に送信され、主要成分判断部２３で生成された主要成分画像、経時変化成分画像及び独立成分係数行列Ａは経時変化表示部３０に送信されるようになっている。

経時変化表示部３０においては、経時変化検出部２０の出力は、経時変化強調画像生成部３１、表示画像選択部３２を経て表示部３３に接続されている。

画像生成部１０において、センサ１１１により得られ撮像部１１により生成された被写体の撮影画像は記憶部１２に記憶され、選択部１３により２つの画像が選択されてデジタル撮影画像Ｆ１、Ｆ２を出力する。

センサ１１１である撮影装置はデジタル画像信号を出力するものであればよく、例えばＦＰＤ(Flat Panel Detector)装置やＣＲ(Computed Radiography)装置やＣＴ(Computed Tomography)装置等のＸ線撮像装置や、ＭＲＩ(Magnetic Resonance Imaging)装置やＵＳ(UltraSonic)エコー装置などが用いられる。

記憶部１２、選択部１３は、例えば前述した撮影装置を制御するコンピュータに付随したハードディスク等の外部記憶装置と、マウス、キーボード等の入力指示装置を使用して実現することができる。即ち、ユーザがこれらの入力指示装置を用いて、ハードディスクに保存された画像を対話的に指定するようにしてもよいし、所定のルールに従って自動的に２枚の画像が選択されるようにしてもよい。

経時変化検出部２０は画像生成部１０から出力された撮影画像Ｆ１、Ｆ２を入力し、後述する主成分画像及び経時変化成分画像を出力する。経時変化検出部２０において、入力した２枚の撮影画像は位置合わせ部２１において相互の位置ずれを補正するための位置合わせが行われ、後続の独立成分分析処理部２２において分析が行われて、画像を構成する各成分が出力され、更に主要成分判断部２３において経時変化検出画像及び経時変化強調画像が生成、出力される。

図２は経時変化検出部２０の独立成分分析処理部２２の構成図である。この独立成分分析処理部２２は一次元化処理部２２１、分析前処理部２２２、分離行列算出部２２３、分析後処理部２２４、二次元化処理部２２５、前処理情報記憶部２２６、二次元情報記憶部２２７により構成されている。

位置合わせ部２１、独立成分分析処理部２２、主要成分判断部２３は、ソフトウェアとして実現することが可能であり、各部はコンピュータ上で動作するソフトウェアを構成するモジュールとして実現することができる。或いは、他の実現方法として、位置合わせ部２１、独立成分分析処理部２２、主要成分判断部２３の一部又は全てをハードウェアにより実装することも可能である。

経時変化表示部３０においては、経時変化表示部３０は入力される各画像を医師の読影に適するように処理、出力するものである。経時変化強調画像生成部３１は、経時変化検出部２０が生成した画像及び画像生成部１０が生成した画像を用いて経時変化強調画像を生成する。表示画像選択部３２は経時変化強調画像生成部３１が生成した画像を選択して、後続の表示部３３に入力する。表示部３３は表示画像選択部３２が出力した画像を表示する。

これらの機能を実施するためには、例えば図３に示すような形態の構成を用いることができる。図３において、コンピュータ４１がネットワーク４２を介して１つ又は複数の画像生成装置４３に接続されている。この画像生成装置４３としては、前述した画像生成部１０の各種の撮影装置が該当する。或いは、図示しない撮影装置により撮影された画像が予めファイルサーバ４４に保存されており、それを画像生成部１０と同等に扱うようにしても、全く同じように適用することができる。

図３において、経時変化検出部２０及び経時変化表示部３０は、コンピュータ４１及びそれに接続された周辺機器により実現することができ、一例として経時変化検出部２０及び経時変化表示部３０は、表示装置４５及びハードディスク４６に格納されたプログラムとして実現可能である。この場合に、図示しないユーザ入力により、ハードディスク４６に保存されたプログラムがバス４７を介してＲＡＭ４８に読み出され、ＣＰＵ４９がＲＡＭ４８内のプログラムを順次に実行して前述した機能を実現し、処理の結果、得られた画像が表示装置４５に表示される。ここで表示装置４５として、例えばＣＲＴモニタ、液晶ディスプレイ等を用いることができる。

また、本発明を実施するプログラムはハードディスク４６に保存されることに限定される必要はなく、外部記憶装置としてコンピュータ４１に接続された光磁気ディスク５０、或いはネットワークを介して接続されるファイルサーバ４４に記憶するようにしてもよい。或いは、本発明による経時変化検出部２０及び経時変化表示部３０は、部分的に或いは全てをコンピュータ４１に装着された専用のハードウェア５１として実現するようにしてもよい。

更に、本発明は図３に示すようなネットワーク４２を介した接続に限定されるものではなく、図４に示すような構成によることもできる。図４の構成においては、コンピュータ４１は画像生成装置４３と直接接続されており、画像を直接に入力する。また、前述したように画像生成装置４３がコンピュータ４１を備えている場合には、その中に前述した形態により撮影装置中に実現することも可能である。

図５、図６は本実施例における全体処理フローチャート図、図７は全体処理フローチャート図の位置合わせに関するサブルーチンのフローチャート図、図８は全体処理フローチャート図の独立成分分析処理に関するサブルーチンのフローチャート図、図９は全体処理フローチャート図の主要成分判断処理に関し相関値を用いるサブルーチンのフローチャート図、図１０は全体処理フローチャート図の主要成分判断処理に関しヒストグラムを用いるサブルーチンのフローチャート図である。

次に、各部の機能を図５〜図１０の各フローチャート図を参照して動作を説明する。先ず、第１の画像Ｆ１と第２の画像Ｆ２を入力する(ステップＳ１１０)。この画像は本実施例の場合に画像生成部１０において撮像部１１により生成される。撮像部１１により生成された画像は、蓄積用の記憶部１２に蓄積され、患者単位、部位単位に時系列的に記憶される。記憶部１２に蓄積された特定患者の特定部位の画像の中から、選択部１３は２つの画像Ｆ１、Ｆ２を選択する。このようにして、ステップＳ１１０における第１の画像Ｆ１と第２の画像Ｆ２が生成される。

続いて、第１の画像Ｆ１と第２の画像Ｆ２を位置合わせする(ステップＳ１２０)。この位置合わせは位置合わせ部２１に入力された２つの画像Ｆ１、Ｆ２に対し、２つの画像Ｆ１、Ｆ２の間に存在する位置ずれを検出し、２つの画像間Ｆ１、Ｆ２の対応画素位置を解剖学的特徴を用いて検出し、一方の画像を他方の画像に合わせるようにワーピング処理を施す。この位置合わせの動作は、先ず第１の画像Ｆ１と第２の画像Ｆ２を概略的にマッチングする(ステップＳ１２１)。その後に、第１の画像Ｆ１と第２の画像Ｆ２を局所的にマッチングする(ステップＳ１２２)。このマッチングは２段階に分けることにより、マッチング精度が高くなる。

マッチングの際には、位置合わせをする際の目印となる関心領域ＲＯＩ(ロイ、Region of interest)を設定する。図１１は設定された１組のＲＯＩを図示したものであり、第１の画像Ｆ１に対してテンプレートロイＴＲが、第２の画像Ｆ２に対してサーチロイＳＲが設定されている。第２の画像Ｆ２に対してテンプレートロイＴＲを設定して、第１の画像Ｆ１に対してサーチロイＳＲを設定してもよい。位置合わせのためのＲＯＩの設定方法及び位置合わせに関しては、特許文献２に詳細に述べられている。

特開平７−３７０７４号公報

マッチングはテンプレートロイＴＲ及びサーチロイＳＲを１組として、複数のＲＯＩの組が入力画像に対して設定し、テンプレートマッチングが行われる。このテンプレートマッチングは図１１に示すように、設定されたテンプレートロイＴＲ内の画像データとサーチロイＳＲ内の画像データとのマッチングを取りながら位置を変化させ、最も一致度の高い位置を求め、テンプレートロイＴＲの設定位置との差を２つの画像Ｆ１、Ｆ２の局所的な位置ずれ量として出力する。

図１１の場合に、第２の画像Ｆ２は第１の画像Ｆ１に比較して右下に平行移動しているため、シフトベクトルは図１２に示すようになる。シフトベクトルは画像上の或る点の座標(ｘ，ｙ)を新しい座標(ｘ’，ｙ’)にシフトさせるという意味があるので、このように呼ばれている。なお、マッチングの具体的な方法としては、例えば非特許文献１等に記載されているＳＳＤＡ法又は相互相関係数法などの公知の方法によることができる。

「画像解析ハンドブック」高木幹雄監修 東京大学出版会

ステップＳ１２４では、入力された複数のシフトベクトルに所定の補間モデルを適用し、補間モデルのパラメータを出力する。補間モデルのパラメータとは、補間モデルを一意に決定するための定数である。例えば、ｙ＝ａｘ²＋ｂｘ＋ｃのｘに関する２次式を補間モデルとして採用すれば、(ａ，ｂ，ｃ)という３つの定数が、補間モデルのパラメータとなる。補間モデルはこのような二次元多項式を使用してもよいし、またスプライン補間など他の補間モデルを使用してもよい。補間モデルの次数は、予め経時変化検出画像の画質評価等に基づいて適切な値が決定してもよいし、補間誤差が所定の値より小さくなるように適応的に決定することもできる。

算出された補間モデルのパラメータを用いて、画像Ｆ１、Ｆ２のうちの一方を他方に合わせる処理を行う(ステップＳ１２５)。この処理は画像変形を行うため、ワーピング(Warping)処理と呼ばれている。ワーピング処理は線形変換による画像変形でも、非線形変換による画像変形でもよいが、本実施例ではアフィン変換を行っている。このワーピング処理により、位置合わせ完了画像Ｘ１、Ｘ２が出力される。このようにして、位置合わせ処理(ステップＳ１２０)が完了するが、ＲＯＩ設定とシフトベクトル補間処理による位置合わせの計算は、複数回行い精度を高めるようにしてもよい。

次に、位置合わせ完了画像Ｘ１、Ｘ２を用いて、独立成分分析処理を行う(ステップＳ１３０)。独立成分分析処理部２２により、位置合わせ部２１から出力された位置合わせ完了画像Ｘ１、Ｘ２を独立成分分析(ＩＣＡ，Independent Component Analysis)で処理して、主要成分画像と経時変化成分画像の２つの画像を生成する。ここで、主要成分とは２つの画像Ｘ１、Ｘ２に共通した構造物であり、経時変化成分とは２つの画像Ｘ１、Ｘ２間の細部の違いである。独立成分分析は主要成分と経時変化成分が統計的に独立な信号であるときに、画像Ｘ１、Ｘ２から主要成分と経時変化成分をそれぞれ抽出する。

２つの画像Ｆ１、Ｆ２が医用画像の場合には、例えば主要成分画像は正常な解剖学的構造が相当し、経時変化成分画像は非正常な解剖学的構造が相当する。同一患者の同一部位を異なる時刻に撮影し、この撮影間隔期間内に非正常な解剖学的変化が生じた場合に、経時変化成分画像にその変化が現れることになる。

例えば、画像Ｆ１、Ｆ２が胸部画像であり、第１の画像Ｆ１は癌が発生する前の撮影画像であり、第２の画像Ｆ２は癌が発生した後の撮影画像であるとする。位置合わせ部２１で位置合わせした完了画像Ｘ１、Ｘ２とすると、主要成分画像は画像Ｘ１と画像Ｘ２に共通する胸部構造の画像であり、経時変化成分画像は癌の画像である。胸部構造と癌は独立な事象であるから、独立成分分析によってこれら２つの事象を抽出できる。

主要成分画像ｓ(ｘ，ｙ)と経時変化成分画像Δ(ｘ，ｙ)は、位置合わせした第１の画像Ｘ１＝ｘ１(ｘ，ｙ)と位置合わせした第２の画像Ｘ２＝ｘ２(ｘ，ｙ)を用いると、次の関係にある。
ｘ１(ｘ，ｙ)＝ａ₁₁・ｓ(ｘ，ｙ)＋ａ₁₂・Δ(ｘ，ｙ)
ｘ２(ｘ，ｙ)＝ａ₂₁・ｓ(ｘ，ｙ)＋ａ₂₂・Δ(ｘ，ｙ) …(１)

独立成分分析により、ｓ(ｘ，ｙ)、Δ(ｘ，ｙ)、係数ａ₁₁、ａ₁₂、ａ₂₁、ａ₂₂が出力され、この出力結果は後述の主要成分判断部２３で使用される。従来の経時差分ＣＡＤにおける差分処理では、濃度補正部１において、係数ａ₁₁とａ₂₁を正しく推定する必要があったが、独立成分分析ではこれらの値は自動的に算出されるため、濃度補正処理は不要である。

ここで、後の説明のために、独立成分係数行列Ａを次の式（２）（［数１］）のように定義しておく。

対象画像が胸部画像で、第１の位置合わせ完了画像Ｘ１が癌発生前の画像であり、第２の完了画像Ｘ２が癌発生後の画像であれば、式(１)において、ａ₁₂＝０と考えてよく、次の式（３）の通りとなる。
ｘ１(ｘ，ｙ)＝ａ₁₁・ｓ(ｘ，ｙ)
ｘ２(ｘ，ｙ)＝ａ₂₁・ｓ(ｘ，ｙ)＋ａ₂₂・Δ(ｘ，ｙ) …(３)

独立成分分析処理により、胸部構造である主要成分画像ｓ(ｘ，ｙ)と、癌に相当する経時変化成分画像Δ(ｘ，ｙ)とが出力される。対象画像は胸部画像である必要はなく、別の部位の画像であってもよい。また、経時変化成分画像は経時的な変化を表す画像であり、癌であるとは限らず、被検体の異常領域が経時変化成分画像として出力される。図１３は独立成分分析処理部２２における入力画像と出力画像の関係を表した説明図を示している。

独立成分分析処理部２２は主要成分画像ｓ(ｘ，ｙ)と経時変化成分画像Δ(ｘ，ｙ)の２つ画像を出力するが、その画像のうちの何れが主要成分画像であり、何れが経時変化成分画像であるかは分からない。そのため、後続の主要成分判断部２３で、何れが主要成分画像ｓ(ｘ，ｙ)であるかを判断する。

独立成分分析処理とは、互いに独立な信号が混在した観測信号から、信号の独立性以外の先見的知識なしに原信号を抽出するデータ分離手法であり、音声や画像処理の分野で注目されている新しい技術である。原信号と結合係数が共に未知であるにも拘らず、繰り返し学習によって結合係数を推定し、観測信号から原信号を分離・抽出することができる。

即ち、Ｎ個の独立な原信号を、
Ｓ(ｔ)＝［ｓ₁(ｔ)，ｓ₂(ｔ)，・・，ｓ_N(ｔ)］^T(ｔ＝０，１，２，・・)
としＮ個の観測信号を、
Ｘ(ｔ)＝［ｘ₁(ｔ)，ｘ₂(ｔ)，・・，ｘ_N(ｔ)］^T(ｔ＝０，１，２，・・)
と表記するとき、観測信号はＸ(ｔ)＝Ａ・Ｓ(ｔ)と記述できる場合において、原信号の独立性と観測信号Ｘ(ｔ)以外の情報を一切用いずに、原信号Ｓ(ｔ)を推定する。ここで、ＡはＮ・Ｎの実数行列であり、［ ］^Tは、行列の転置である。

行列Ａの逆行列Ｗが見付かれば、原信号はＳ(ｔ)＝ＷＸ(ｔ)という演算により復元可能である。独立成分分析では、各原信号ｓ₁(ｔ)が互いに統計的に独立な信号であると仮定して、行列Ｗを算出する。つまり、Ｙ(ｔ)＝ＷＸ(ｔ)という演算をして、得られる各信号Ｙ(ｔ)＝［ｙ₁(ｔ)，ｙ₂(ｔ)，・・，ｙ_N(ｔ)］^T(ｔ＝１，２，・・)が、互いに統計的に独立な信号になるようにＷを算出する。行列Ｗは原信号Ｓ(ｔ)を混合した信号であるＸ(ｔ)から原信号を分離するため、分離行列と呼ばれている。

分離行列Ｗを求める方法は、逐次的に学習しながら求める方法と一括して求める方法がある。ここでは、逐次的に学習する方法について説明する。また、統計的に独立かどうかを測る尺度は、相互情報量や非ガウス性などがあるが、ここでは相互情報量を採用する。相互情報量は非負の値を持ち、互いに統計的に独立な場合だけ、０となる性質を持つ尺度である。従って、信号Ｙ(ｔ)の各成分間の相互情報量を、なるべく０に近付けるような行列Ｗを探すことになる。

学習方法は、以下のステップで行われる。

ステップ（１）：分離行列Ｗの初期化を行う。

ステップ（２）：信号Ｙ(ｔ)＝ＷＸ(ｔ)を計算する。

ステップ（３）：分離行列Ｗを更新式Ｗ←Ｗ＋η・ΔＷに従って更新する。
ηは学習係数、ΔＷは更新の変化分、相互情報量をＩ(Ｙ)とすると、
(ΔＷ)_ij∝−{(∂Ｉ(Ｙ))／∂Ｗ_ij}Ｗ^TＷ
となる。これは、相互情報量Ｉ(Ｙ)に関する最急降下方向である。

ステップ（４）：εは収束条件とし、(ΔＷ)_ij＜εを満たすまで、ステップ（２）〜（３）を繰り返す。

ここで、相互情報量について補足する。相互情報量Ｉ(Ｙ)は信号Ｙ(ｔ)の同時分布確率密度関数をｐ(Ｙ)＝ｐ(ｙ₁，ｙ₂，・・・，ｙ_N)とおくと、次式で定義される。

先に述べたように、Ｉ(Ｙ)≧０であり、ｙ₁，ｙ₂，・・・，ｙ_Nが互いに独立な場合、つまり、

のとき、Ｉ(Ｙ)＝０となる。従って、相互情報量Ｉ(Ｙ)を最小化することにより、Ｙの各成分を独立にすることが可能になる。ステップ（３）では、相互情報量Ｉ(Ｙ)の最小値を最急降下法により求めている。

次に、相互情報量Ｉ(Ｙ)を分離行列Ｗで表現することを考えると、ステップ（３）で示したように、相互情報量Ｉ(Ｙ)が小さくなるような方向に分離行列Ｗを更新するときに必要になる。ＹのエントロピーＨ(Ｙ)は、
Ｈ(Ｙ)＝−∫ｐ(Ｙ) logｐ(Ｙ)ｄＹ
である。Ｙ＝ＷＸより、
ｐ(Ｙ)＝ｐ(Ｘ)／ det(Ｗ)
ｄＹ＝det(Ｗ)ｄＸ
となるので、ＹのエントロピーＨ(Ｙ)は次式となる。

Ｈ(Ｙ)＝Ｈ(Ｘ)＋ log｜det(Ｗ)｜
相互情報量Ｉ(Ｙ)はエントロピーＨ(Ｙ)を用いて、

となるから、結局は、

となり、相互情報量Ｉ(Ｙ)を分離行列Ｗの関数で表現できる。

次に、ステップ（３）における更新の変化分ΔＷを具体的に計算する。

ここで、最後の式の右辺の第１項は、非特許文献２を参照すると、次式で近似することができる。

ここでφ(ｙ_i)は活性化関数と呼ばれており、非線形関数であり、例えば次のような関数である。
φ(ｙ_i)＝ｙ_i ²

S. Haykin, Neural Networks, Prentice Hall, 1999.

また、次のような関数でもよい。
φ(ｙ_i)＝(１／２)ｙ_i ⁵＋(２／７)ｙ_i ⁷＋(１５／２)ｙ_i ⁹＋(２／１５)ｙ_i ¹¹

更に、det（Ｗ）はＷのｊ行ｋ列の要素をＷ_jkとするとき、余因子Ｗ_jkを用いて、次のように余因子分解できる。

従って、右辺の第２項の対数のＷ_jkによる微分は、

以上をまとめると、ΔＷは次式となる。
ΔＷ_jk＝{(Ｗ^-T)_jk−φ(ｙ_i)ｘ_k}Ｗ^TＷ

これを行列表現すると、次式となる。
ΔＷ＝{(Ｗ^-T)−φ(Ｙ)Ｘ^T}Ｗ^TＷ

ここで、Ｙ^T＝Ｘ^TＷ^Tという関係を使えば、
ΔＷ＝{(Ｗ^-T)−φ(Ｙ)Ｘ^T}Ｗ^TＷ
＝{Ｉ−φ(Ｙ)Ｘ^TＷ^-T}Ｗ
＝{Ｉ−φ(Ｙ)Ｙ^T}Ｗ
となり、ステップ（３）における具体的な更新の変化分ΔＷを得る。なお、上式のＩは単位行列であり、詳細は非特許文献３に記載されている。

甘利俊一，村田昇 共著，独立成分分析，臨時別冊・数理科学 「ＳＧＣライブラリ １８」，サイエンス社(２００２)

このようにして、位置合わせ完了画像Ｘ１、Ｘ２は、一次元化処理部２２１で二次元データ形式から一次元データ形式に変換される。変換後の完了画像Ｘ１、Ｘ２を次式のように表記する。なお、ｎは画像の画素数である。
Ｘ(ｎ)＝［ｘ₁（ｎ），ｘ₂（ｎ）］^T(ｎ＝１，２，・・，Ｎ)

ｘ₁（ｎ）とｘ₂（ｎ）の構成方法は、図１４(ａ)に示すように画像の左隅から順番に取り出して並べる方法でもよいし、図１４(ｂ)に示すようにジグザグに並べる方法でもよいが、ｘ₁（ｎ）とｘ₂（ｎ）の構成方法は同じでなければならない。変換で用いられた並べ方情報は、二次元情報記憶部２２７に保持される。この情報は分析処理後に、一次元データ形式から二次元データ形式に変換する際に参照される。

一次元データ形式に変換された後に、Ｘ(ｎ)は分析前処理部２２２で独立成分分析処理を精度良く行うために整形される。整形で使用された情報は、前処理情報記憶部２２６に保持されて、分析処理後の整形逆変換を行うときに使用される。本実施例では、ｘ₁（ｎ）の平均値とｘ₂（ｎ）の平均値が０に等しくなるように、直流成分値をそれぞれ差し引く。差し引いた直流成分値情報は、前処理情報記憶部２２６に保持しておく。

整形されたＸ(ｎ)は、独立成分分析処理の中心となる分離行列算出部２２３に入力される。分離行列算出部２２３は分離行列Ｗを算出して、独立成分係数行列Ａと独立成分Ｙ(ｎ)を出力する。

この分離行列Ｗは先に述べた学習ステップ（１）〜（４）により算出される。

ステップ(１)：分離行列Ｗの初期化を行う(ステップＳ１３１)。本実施例では、次のように初期化する。

ステップ(２)：独立成分Ｙ(ｔ)＝ＷＸ(ｔ)を計算する(ステップＳ１３２)。

ステップ(３)：分離行列Ｗを更新式に従って更新する(ステップＳ１３３)。更新式は前述したように、次の通りである。
Ｗ←Ｗ＋ηΔＷ
ΔＷ＝{Ｉ−φ(Ｙ)Ｙ^T}Ｗ
φ(ｙ_i)＝(１／２)ｙ_i ⁵＋(２／７)ｙ_i ⁷＋(１５／２)ｙ_i ⁹＋(２／１５)ｙ_i ¹¹
η＝１

ステップ(４)：(ΔＷ)_ij＜εを満たすまで、ステップ(２)〜(３)を繰り返す(ステップＳ１３４)。収束条件εについては、本実施例では０．０００１とする。なお、本実施例では学習係数ηを定数としたが、更新回数と共にηを小さくするような更新方法でもよい。

ステップ(５)：独立成分Ｙ(ｎ)＝［ｙ₁(ｎ)，ｙ₂(ｎ)］^T(ｎ＝１，２，・・・，Ｎ)を得る(ステップＳ１３５)。また、Ａ＝Ｗ^-1を計算することにより、次の独立成分係数行列Ａを得る。

分離行列算出部２２３により出力された独立成分Ｙ(ｎ)は、前処理情報記憶部２２６に保持されている情報を使用して、分析後処理部２２４で整形逆変換が行われる。本実施例の場合に、前処理情報記憶部２２６に保持されている直流成分値を、独立成分Ｙ(ｎ)に足し合わせることにより整形逆変換を行う。

分析後処理部２２４で整形逆変換された独立成分Ｙ(ｎ)は、二次元情報記憶部２２７に保持されている情報を使用して、二次元化処理部２２５で一次元データ形式から二次元データ形式に変換され、独立成分画像Ｙ１、Ｙ２として出力される。

ここで、従来の経時差分ＣＡＤとの構成の違いについて述べると、従来の経時差分ＣＡＤでは、第１の画像と第２の画像を入力とし、差分処理の結果として１つの画像を出力する。本実施例では、第１の画像と第２の画像を入力とし、独立成分分析処理の結果として２つの画像を出力する。また、被検体の異常領域である経時変化成分画像Δ(ｘ，ｙ)は、従来例による差分処理の結果と似たような画像ではあるが、同じではない。

主要成分判断部２３では、独立成分分析処理部２２における２つの出力画像のうち、何れが主要成分画像であり、何れが経時変化成分画像であるかを判断する(ステップＳ１４０)。

本実施例では、主要成分画像と経時変化変化画像の判断を相互相関値を用いて判断する。独立成分分析処理部４における処理結果により、次の関係式が得られる。
ｘ１(ｘ，ｙ)＝ａ₁₁・ｙ１(ｘ，ｙ)＋ａ₁₂・ｙ２(ｘ，ｙ)
ｘ２(ｘ，ｙ)＝ａ₂₁・ｙ１(ｘ，ｙ)＋ａ₂₂・ｙ２(ｘ，ｙ) …(４)

ｙ１(ｘ，ｙ)とｙ２(ｘ，ｙ)の何れが主要成分画像であるかを判断する。経時変化成分画像は主要成分画像に比較して信号値が非常に小さいので、式(４)の第１式において、左辺ｘ１(ｘ，ｙ)と第１項ａ₁₁・ｙ１(ｘ，ｙ)との相関値Ｒ₁₁の絶対値と、左辺ｘ１(ｘ，ｙ)と第２項ａ₁₂・ｙ２(ｘ，ｙ)との相関値Ｒ₁₂の絶対値とを比較し、大きい方を主要成分画像と判断する。つまり、
｜Ｒ₁₁｜＞｜Ｒ₁₂｜
であれば、ｙ１(ｘ，ｙ)が主要成分画像であると判断し、
｜Ｒ₁₁｜＜｜Ｒ₁₂｜
であれば、ｙ２(ｘ，ｙ)が主要成分画像であると判断する(ステップＳ１４１ａ、Ｓ１４２ａ)。相関の絶対値｜Ｒ₁₁｜と｜Ｒ₁₂｜がほぼ等しい場合には、独立成分分析処理部２２が主要成分画像と経時変化成分画像に正常に分離できなかったと判断する。

相関値の算出方法としては、本実施例では相互相関値を使用したが、これに限定されることはなく、他の相関値算出方法を使用してよい。ｙ１(ｘ，ｙ)とｙ２(ｘ，ｙ)のうち、一方が主要成分画像として判断されれば、他方が経時変化成分画像となる。このようにして、主要成分判断部２３は主要成分画像と経時変化検出画像を出力する(ステップＳ１５０)。

経時変化強調画像生成部３１は経時変化検出部２０が生成した画像を用いて、経時変化検出画像と経時変化強調画像を生成する(ステップＳ１６０)。経時変化検出画像ｈ１(ｘ，ｙ)、ｈ２(ｘ，ｙ)は、経時変化成分画像とその対応する係数を掛け合わせて生成する。ｙ２(ｘ，ｙ)が経時変化成分画像とすると、次の式を生成する。
ｈ１(ｘ，ｙ)＝ａ₁₂・ｙ２(ｘ，ｙ)
ｈ２(ｘ，ｙ)＝ａ₂₂・ｙ２(ｘ，ｙ)

更に、経時変化強調画像ｗ１(ｘ，ｙ)、ｗ２(ｘ，ｙ)、ｋ１(ｘ，ｙ)、ｋ２(ｘ，ｙ)、ｖ１(ｘ，ｙ)、ｖ２(ｘ，ｙ)を、次のように生成する。なお、ｃは十分に大きな定数である。
ｗ１(ｘ，ｙ)＝ｃ・ｈ１(ｘ，ｙ)
ｗ２(ｘ，ｙ)＝ｃ・ｈ２(ｘ，ｙ)
ｋ１(ｘ，ｙ)＝ｙ１(ｘ，ｙ)＋ｃ・ｈ１(ｘ，ｙ)
ｋ２(ｘ，ｙ)＝ｙ２(ｘ，ｙ)＋ｃ・ｈ２(ｘ，ｙ)
ｖ１(ｘ，ｙ)＝ｘ１(ｘ，ｙ)＋ｃ・ｈ１(ｘ，ｙ)
ｖ２(ｘ，ｙ)＝ｘ２(ｘ，ｙ)＋ｃ・ｈ２(ｘ，ｙ)

経時変化強調画像ｗ１(ｘ，ｙ)、ｗ２(ｘ，ｙ)は、経時変化検出画像を強調した画像である。経時変化部分が微小な変化のため、経時変化があるのかないのか分からない場合でも、経時変化強調画像を見れば、経時変化があるかどうかが一目で判断できる。

経時変化強調画像ｋ１(ｘ，ｙ)、ｋ２(ｘ，ｙ)は、同様に経時変化を強調して、主要成分画像の上に重ね書きした画像である。主要成分画像の上に経時変化部分が重ねて表示されているため、経時変化部分が主要成分部分のどの位置にあるかを判断できる。

経時変化強調画像ｖ１(ｘ，ｙ)、ｖ２(ｘ，ｙ)は、同様に被写体の撮影画像の上に経時変化部分が重ねて表示されているため、経時変化部分が被写体のどの位置にあるかを判断できる。

図１５に示す説明図は、経時変化強調画像が位置合わせ完了画像や経時変化検出画像や主要成分画像と、どのように関係しているかを表している。

画像表示選択部３２では、経時変化強調画像生成部３１が出力した画像を診断者が選択する。即ち、主要成分画像、経時変化成分画像、経時変化検出画像、経時変化強調画像、位置合わせ完了画像を複数又は１つだけ選択する。表示部３３では、画像表示選択部３２が選択した画像を同時に表示する。

実施例２として、主要成分判断部をヒストグラムによる判別する方式を説明する。ｙ１(ｘ，ｙ)とｙ２(ｘ，ｙ)の両者に対してヒストグラムを算出して、ヒストグラムの分布幅が小さい方を経時変化成分画像と判断する(ステップＳ１４１ｂ、Ｓ１４２ｂ)。主要成分画像の画素値は、最小値から最大値まで全ゆる値が分布しているが、経時変化成分画像は或る狭い範囲にだけ分布しているためである。ｙ１(ｘ，ｙ)とｙ２(ｘ，ｙ)のうち、一方が経時変化成分画像として判断されれば、他方が主要成分画像となる。このようにして、主要成分判断部２３は主要成分画像と経時変化検出画像を出力する(ステップＳ１５０)。

実施例３として、強調画像生成をハイパスフィルタ等による行う方式を説明する。

この実施例３では、経時変化強調画像ｗ１(ｘ，ｙ)、ｗ２(ｘ，ｙ)、ｋ１(ｘ，ｙ)、ｋ２(ｘ，ｙ)、ｖ１(ｘ，ｙ)、ｖ２(ｘ，ｙ)を次式のように生成する。
ｗ１(ｘ，ｙ)＝ｕ(ｘ，ｙ)＊ｈ１(ｘ，ｙ)
ｗ２(ｘ，ｙ)＝ｕ(ｘ，ｙ)＊ｈ２(ｘ，ｙ)
ｋ１(ｘ，ｙ)＝ｙ１(ｘ，ｙ)＋ｕ(ｘ，ｙ)＊ｈ１(ｘ，ｙ)
ｋ２(ｘ，ｙ)＝ｙ２(ｘ，ｙ)＋ｕ(ｘ，ｙ)＊ｈ２(ｘ，ｙ)
ｖ１(ｘ，ｙ)＝ｘ１(ｘ，ｙ)＋ｕ(ｘ，ｙ)＊ｈ１(ｘ，ｙ)
ｖ２(ｘ，ｙ)＝ｘ２(ｘ，ｙ)＋ｕ(ｘ，ｙ)＊ｈ２(ｘ，ｙ)

ここで、ｕ(ｘ，ｙ)は高周波数成分を強調するハイパスフィルタ関数であり、演算子＊は畳み込み演算を意味している。

経時変化強調画像ｗ１(ｘ，ｙ)、ｗ２(ｘ，ｙ)、ｋ１(ｘ，ｙ)、ｋ２(ｘ，ｙ)、ｖ１(ｘ，ｙ)、ｖ２(ｘ，ｙ)は、経時変化を強調して、経時変化部分がどの位置にあるのかを判断できる。

実施例１、２では、経時変化検出画像を定数倍することにより経時変化強調画像を生成したが、本実施例３では高周波数を強調するハイパスフィルタを使用して生成している。生成方法はこれらに限定されることなく、非線形フィルタなど各種フィルタを使用して生成してもよい。

画像表示選択部３２では、経時変化強調画像生成部３１が出力した画像を選択する。即ち、主要成分画像、経時変化成分画像、経時変化検出画像、経時変化強調画像、位置合わせ完了画像を複数又は１つだけ選択する。実施例１、２では診断者が所望の画像を選択したが、本実施例では、画像表示選択部３２が所定のルールにより画像を選択する。この所定のルールとは、例えば「初めて診断する患者は、詳細な情報が必要なので上記の画像を全て選択する」、「経過観察中の患者は、関心領域を把握しているので経時変化検出画像だけ選択する」などである。

実施例４として、表示形態の画像自動選択、順次表示等による別方式を説明する。この本実施例４では、統計的に独立かどうかを測る尺度として相互情報量を用いるが、４次のキュムラントや更に高次の統計量を用いてもよい。また、独立成分分析アルゴリズムとして、Fast ICAや、JADEを使用してもよい。これらのアルゴリズムは、先の非特許文献３に解説されている。

本実施例では、分離行列Ｗを求める方法は、逐次的に学習しながら求める方法を使用しているが、一括して求める方法でもよい。逐次的に学習しながら求める方法は、解が収束するまで時間を要するが、演算用計算領域が少量で済むという利点がある。一方、一括して求める方法は、解を求める時間は長くないが、演算用計算領域が多量に必要という問題がある。何れの方法を採用するかの判断は、画像サイズの大きさと計算機の能力に依存する。

分離行列Ｗは実施例１に記載した学習ステップ（１）〜（４）で算出される。ステップ(５)で［数１１］の独立成分係数行列Ａを得る。

分析前処理部２２２で直流成分値を差し引いた後に、Ｘ(ｎ)を無相関化するために白色化処理を施しておいてもよい。白色化処理を施すと、分離行列Ｗを少ないステップ数で収束する傾向がある。ただし、これは画像Ｆ１、Ｆ２のＳ／Ｎ比が大きい場合に限られている。

分離行列算出部２２３により出力された独立成分Ｙ(ｎ)は、前処理情報記憶部２２６に保持されている情報を使用して、分析後処理部２２４で整形逆変換が行われる。本実施例４の場合に、前処理情報記憶部２２６に保持されている直流成分値を独立成分Ｙ(ｎ)に足し合わせることにより、整形逆変換を行う。

分析後処理部２２４で整形逆変換された独立成分Ｙ(ｎ)は、二次元情報記憶部２２７に保持されている情報を使用して、二次元化処理部２２５において一次元データ形式から二次元データ形式に変換され、独立成分画像Ｙ１、Ｙ２として出力される。

ステップＳ１４０、ステップＳ１４１ａ、Ｓ１４２ａを実施例１と同様に経由して、ステップＳ１５０で主要成分判断部２３は主要成分画像と経時変化検出画像を出力する。

経時変化強調画像生成部３１はステップＳ１６０で、実施例１と同様に折り畳み演算を行って次式を生成する。
ｈ１(ｘ，ｙ)＝ａ₁₂・ｙ２(ｘ，ｙ)
ｈ２(ｘ，ｙ)＝ａ₂₂・ｙ２(ｘ，ｙ)

更に、経時変化強調画像ｗ１(ｘ，ｙ)、ｗ２(ｘ，ｙ)、ｋ１(ｘ，ｙ)、ｋ２(ｘ，ｙ)、ｖ１(ｘ，ｙ)、ｖ２(ｘ，ｙ)を、実施例３と同様に次式のように生成する。
ｗ１(ｘ，ｙ)＝ｕ(ｘ，ｙ)＊ｈ１(ｘ，ｙ)
ｗ２(ｘ，ｙ)＝ｕ(ｘ，ｙ)＊ｈ２(ｘ，ｙ)
ｋ１(ｘ，ｙ)＝ｙ１(ｘ，ｙ)＋ｕ(ｘ，ｙ)＊ｈ１(ｘ，ｙ)
ｋ２(ｘ，ｙ)＝ｙ２(ｘ，ｙ)＋ｕ(ｘ，ｙ)＊ｈ２(ｘ，ｙ)
ｖ１(ｘ，ｙ)＝ｘ１(ｘ，ｙ)＋ｕ(ｘ，ｙ)＊ｈ１(ｘ，ｙ)
ｖ２(ｘ，ｙ)＝ｘ２(ｘ，ｙ)＋ｕ(ｘ，ｙ)＊ｈ２(ｘ，ｙ)

実施例１、２では、経時変化検出画像を定数倍することにより経時変化強調画像を生成したが、本実施例４では実施例３と同様に、高周波数を強調するハイパスフィルタを使用している。

表示部３３では、画像表示選択部３２が選択した画像を表示する。実施例１、２、３では、画像表示選択部３２が選択した画像を同時に表示したが、本実施例４では、１枚ずつ順次に表示する。表示方法はこれに限定されることなく他の方法、例えば、眼鏡方式の表示装置を用いて、右目と左目に異なる画像を表示してもよい。

なお、本発明の目的は、実施形態の装置又はシステムの機能を実現するソフトウェアのプログラムコードを記憶した記憶媒体を、装置又はシステムに供給し、その装置又はシステムのコンピュータ（ＣＰＵ又はＭＰＵ等）が記憶媒体に格納されたプログラムコードを読み出して実行することによっても、達成されることは云うまでもない。

この場合に、記憶媒体から読み出されたプログラムコード自体が実施形態の機能を実現することとなり、そのプログラムコードを記憶した記憶媒体及び当該プログラムコードは本発明を構成することになる。

プログラムコードを供給するための記憶媒体としては、ＲＯＭ、フロッピ（登録商標）ディスク、ハードディスク、光ディスク、光磁気ディスク、ＣＤ−ＲＯＭ、ＣＤ−Ｒ、磁気テープ、不揮発性のメモリカード等を用いることができる。

また、コンピュータが読み出したプログラムコードを実行することにより、実施形態の機能が実現されるだけでなく、そのプログラムコードの指示に基づいて、コンピュータ上で稼動しているＯＳ等が実際の処理の一部又は全部を行い、その処理によって実施例１〜４の機能が実現される場合も、本発明の実施の態様に含まれることは勿論である。

更に、記憶媒体から読み出されたプログラムコードが、コンピュータに挿入された機能拡張ボードやコンピュータに接続された機能拡張ユニットに備わるメモリに書き込まれた後に、そのプログラムコードの指示に基づき、その機能拡張ボードや機能拡張ユニットに備わるＣＰＵなどが実際の処理の一部又は全部を行い、その処理によって実施形態の機能が実現される場合も、本発明の実施の態様に含まれることは云うまでもない。

このようなプログラム又は当該プログラムを格納した記憶媒体に本発明が適用される場合に、当該プログラムは例えば上述の図５〜図１０に示されるフローチャート図に対応したプログラムコードから構成される。

画像処理システムの構成図である。 独立成分分析処理部の構成図である。 実施例における経時変化検出部の構成図である。 実施例における経時変化表示部の構成図である。 全体処理フローチャート図である。 全体処理フローチャート図である。 位置合わせに関するサブルーチンのフローチャート図である。 独立成分分析処理に関するサブルーチンのフローチャート図である。 主要成分判断処理に関し相関値を用いるサブルーチンのフローチャート図である。 主要成分判断処理に関しヒストグラムを用いるサブルーチンのフローチャート図である。 １組のＲＯＩ設定の説明図である。 シフトベクトルに関する説明図である。 独立成分分析処理の説明図である。 独立成分分析処理の一部の説明図である。 経時変化強調処理の説明図である。 従来例の説明図である。

符号の説明

１０ 画像生成部
１１ 撮像部
１２ 記憶部
１３ 選択部
２０ 経時変化検出部
２１ 位置合わせ部
２２ 独立成分分析処理部
２３ 主要成分判断部
３０ 経時変化表示部
３１ 経時変化強調画像生成部
３２ 表示画像選択部
３３ 表示部
１１１ センサ
２２１ 一次元化処理部
２２２ 分析前処理部
２２３ 分離行列算出部
２２４ 分析後処理部
２２５ 二次元化処理部
２２６ 前処理情報記憶部
２２７ 二次元情報記憶部

Claims (12)

1. ２つの画像を入力する画像入力工程と、前記２つの画像の位置合わせをする位置合わせ工程と、該位置合わせ工程で位置合わせをした画像を用いて独立成分分析処理を行う独立成分分析処理工程とを有することを特徴とする画像処理方法。
2. 前記独立成分分析処理工程において出力された処理結果を、主要成分画像と経時変化成分画像とに区別するための判断を行う主要成分判断工程を有することを特徴とする請求項１に記載の画像処理方法。
3. 前記独立成分分析処理工程は、前記位置合わせ工程で位置合わせを行った画像を二次元データ形式から一次元データ形式に変換する一次元化処理工程と、該一次元化処理工程で得られた一次元データから直流成分を減算する分析前処理工程と、該分析前処理工程で得られた一次元データから分離行列、独立成分係数行列及び独立成分を算出する分離行列算出工程と、該分離行列算出工程で算出した独立成分に前記直流成分を加算する分析後処理工程と、分析後処理工程で直流成分を加算した独立成分を二次元データ形式に変換する二次元化処理工程とを備えることを特徴とする請求項１又は２に記載の画像処理方法。
4. 前記主要成分判断工程は、前記独立成分分析処理から出力された処理結果と、前記独立成分分析処理工程で用いた前記位置合わせした画像と他方の画像との２つの画像のうちの一方との類似性を調べて前記主要成分画像を判断することを特徴とする請求項２又は３に記載の画像処理方法。
5. 前記主要成分判断工程は、前記独立成分分析処理から出力された処理結果をヒストグラム処理して、該ヒストグラムの分布幅が小さい方を前記経時変化成分画像と判断することを特徴とする請求項２又は３に記載の画像処理方法。
6. 前記主要成分判断工程で判断した前記経時変化成分画像をフィルタ処理により強調した画像を、前記画像入力工程に入力された画像又は前記独立成分分析処理工程で用いた前記位置合わせした画像又は前記主要成分判断工程で判断した前記主要成分画像と重ね合わせることにより、経時変化部分を強調した画像を生成する経時変化強調画像生成工程を有することを特徴とする請求項２〜６の何れか１つの請求項に記載の画像処理方法。
7. 前記２つの画像は異なる時点で撮影された同一部位の画像であることを特徴とする請求項１〜６の何れか１つの請求項に記載の画像処理方法。
8. 前記位置合わせ工程は前記２つの画像の対応画素位置を解剖学的特徴を用いて検出し、一方の画像を他方の画像に合わせるようにワーピング処理を施す工程であることを特徴とする請求項１〜７の何れか１つの請求項に記載の画像処理方法。
9. ２つの画像から独立成分分析処理を行う画像処理方法をコンピュータに実行させるためのプログラムであって、２つの画像を入力する画像入力工程と、前記２つの画像の位置合わせをする位置合わせ工程と、該位置合わせ工程で位置合わせをした画像を用いて独立成分分析処理を行う独立成分分析処理工程をコンピュータに実行させるためのプログラム。
10. ２つの画像から独立成分分析処理を行う画像処理方法をコンピュータに実行させるためのプログラムを記録したコンピュータ読み取り可能なコンピュータ可読記憶媒体であって、前記２つの画像を入力する画像入力工程と、前記２つの画像の位置合わせをする位置合わせ工程と、該位置合わせ工程で位置合わせをした画像を用いて前記独立成分分析処理を行う独立成分分析処理工程をコンピュータに実行させるためのプログラムを記録したコンピュータ読み取り可能なコンピュータ可読記憶媒体。
11. ２つの画像を入力する画像入力手段と、前記２つの画像の位置合わせをする位置合わせ手段と、該位置合わせ手段で位置合わせをした画像を用いて独立成分分析処理を行う独立成分分析処理手段とを有することを特徴とする画像処理装置。
12. 請求項１１に記載の画像処理装置と、画像生成装置、表示装置、ファイルサーバのうちの少なくとも１つとがネットワークを介して接続されて成ることを特徴とする画像処理シ
    ステム。

Images