JP5753791B2 - ノイズ除去された所定の解像度の医療画像を提供する方法、所定の解像度のノイズ除去された所定の解像度の医療画像を提供するシステム - Google Patents

Description

本発明は、医療画像（以下、単に「画像」とも称する）の処理を行う方法とシステムに関し、特に、ノイズを低減し更に他の好ましい特徴を有する医療画像を生成する方法とシステムに関する。更に本発明は、医療画像内にノイズが不均一に分布することを考慮して、鮮明な詳細画像を得る非線形フィルタを用いるＣＴ画像の処理方法に関する。

非特許文献１は、画像のラプアシアン・ピラミッド分解方法（Lplacian pyramid decomposition）を用いる様々なカメラの焦点の設定で得られたシーンの２つのイメージを合成（融合：fuse）する方法を開示する。非特許文献２は、ラプアシアン・ピラミッド分解法の代わりに、ウェーブレット変換法（wavelet transform）を用いて様々な種類のセンサを用いて得られた同一領域の様々な画像を合成することを開示する。

非特許文献３は、ＣＴ、ＭＲＩ、ＰＥＴ、ＭＲＩ等の様々な画像生成方法を用いてあるいはウェーブレットピラミッド法を用いて得られた２つの医療画像を合成する方法を開示する。

米国特許公開第２００８／０３１０６９５号明細書 米国特許公開第２００８／０１１８１２８号明細書 米国特許公開第２００７／９５３４７７号明細書 米国特許公開第２００３／０９９４０５号明細書 韓国特許公開第２００５／００３１２１０号公報 日本特許公開第２０００／０５０１０９号公報 日本特許公開第２００１／３９８７４号公報 米国特許第６４５９７５５号公報 ヨーロッパ特許第１７７４８３７号明細書 国際公開パンフレット第２００７／０４７５９９号公報

E.H. Adelson, C.H. Anderson, J.R. Bergen, P.J. Burt and J.M. Ogden, "Pyramid methods in image processing", RCA Engineer, 29-6, Nov. 1984 Li Hui, B.S. Manjunath, and S.K. Mitra, "Multi-sensor image fusion using the wavelet transform," in Proceedings of IEEE International Conference on Image Processing, 1994 Yu Lifeng, Zu Donglin, Wan Weidong, and Bao Shanglian, "Multi-Modality Medical Image Fusion Based on Wavelet Ptramid and Evaluation," a report published by the Institute of Heavy Ion Physics, Peking University,Beijing 100871, China (2001) Hassam EI-Din Moustafa and Sameh Rehan, "Applying Image Fusion Techniques for Detection of Hepatic Lesions," Proceedings of the 6th WSEAS International Conference on Wavelet Analysis & Multirate Systems Bucharest, Romania, October 16-18, 2006, pages 41-44, Richard Alan Peters II, "A New Algorithm for Image Noise Reduction using Mathematical Morphology", IEEE Trans. Image Processing 4, 554-568(1995), A.J. Britten,M. Crotty, H.Kiremidjiam, A. Grundy, and E.j.Adam,"The addition of computer simulated noise to investigate radiation dose and image quality in images with spatial correlation of statistical noise: an example application to X-ray CT of tha brain" The British Journal of Radiology, 77(2004), 323-328. C. Tomasi and R. Manduchi, "Bilateral filtering for gray and color images", in Proceedings of the 6th International Conference in Computer Vision(ICCV), 1998, pp.839-846. by J. Weickert, "Coherence-Enhancing Diffusion Filtering", International Journal of Computer Vision, 31(2-3), pp. 111-127,1999, L.Rudin, S. Osher, and E. Fatemi, "Nonlinear total variation based noise removal algorithms," Physica D60, 259-268(1992). M. Mahmoudi, and G. Sapiro, "Fast image andvideo denoising vianonlocal means of similar neighborhoods,"IEEE, Signal Proc., vol.12, no. 12, pp. 839-842, Dec.2005 A. Heidarzadeh and A.n. Avanaki, "An Enhanced Nonlocal Means Algorithm for Image Denoising," 9th ISSPA, Fed. 2007, by E.H.Adelson, C.H. Anderson, J.P. Bergen, P.J. Burt and J.M. Ogden, "Pyramid Methods in Image Processing," RCA Engineer 29, no. 6, Nov-Dec1084, pp. 33-41, or in "Nonlinear Multiresolution Signal Decomposition Schemes- part I: Morphological Pyramids," John Goutsias, and Henk J.A.M. Heijmans, IEEE Transactions On Image Processing, vol.9, No.11, November 2000. A. buades, B. Coll, and J.M.Morel, "On image denoising methods2, Centre de Mathematiques et de leurs Application (CMLA) publication No. 2004-15,2004

非特許文献４は、ＣＴの合成結果とＭＲＩ画像の合成結果を様々な方法を用いて比較することを開示している。様々な方法には、例えばラプアシアンピラミッド法、ウェーブレット変換法、計算効率を上げたピクセルレベルイメージフィージョン法（Pixel-level Image Fusion method）、空間周波数に基づくマルチフォーカス技術（Multi-focus technique based on Spatial Frequency）等がある。

非特許文献５は、これらの微小（薄い）特徴を保存しながら、かつノイズを取り除く、morphological image cleaning algorithmを開示する。ここに開示された方法は、複数のスケールに基づいて残差画像を計算するがこれは、morphological size分配を用いて行い、残差画像内で、ノイズを含むあるいはノイズを出すと判定された様々な領域を廃棄して、微小特徴よりも小さなダイナミックレンジを提供する。

特許文献１は、ＭＲＩ画像を局所適用型非線形ノイズフィルタを用いてノイズを除去する方法を開示するが、これはＭＲＩ画像内に理論的に予測されるノイズレベルの空間的変動を考慮に入れて行っている。
特許文献２は、人工的に生成された所定量のノイズを有するシミュレートされた画像を生成する方法を開示する。

本発明の一実施例は、医療画像（ＣＴ画像）の処理し、所望の特徴を有する出力画像を生成する方法とシステムに関する。この所望の特徴の例は、ノイズの低減、ノイズの予測レベルと空間分布、予測解像度、ノイズや解像度が低減された場合でも保持される詳細と構造がある。
本発明の一実施例は、更に他の好ましい特徴を有する医療画像を生成する方法とシステムに関する。更に本発明は、医療画像内にノイズが不均一に分布することを考慮して鮮明な詳細画像を得る非線形フィルタを用いるＣＴ画像の処理方法に関する。

本発明の一実施例によれば、所定の解像度でノイズ除去された医療画像を提供する方法において、前記医療画像は、表示対象である部分の元の厚さより薄い複数のスライスからなる３次元医療画像であり、
（ａ）前記所定の解像度よりも高い解像度で得られた源医療画像をノイズ除去するステップと、
（ｂ）前記ステップ（ａ）の後の医療画像の解像度を、前記所定の解像度に戻すステップと、
を有し、
前記（ｂ）ステップは、前記複数の薄いスライスを非線形解像度低減手順を用いて合成し、厚いスライスの医療画像を得るステップを含み、
前記非線形解像度低減手順は、
（ｘ）前記薄いスライスを帯域通過成分に分解するステップと、
（ｙ）前記帯域通過成分の一部に対し、より大きな重み係数をより大きな振幅を有する薄いスライス成分に与えるようにして、前記薄いスライス成分を厚いスライス成分になるよう組み合わせるステップと、
（ｚ）前記薄いスライスの帯域通過成分を再結合して、所定の解像度の厚いスライスの医療画像を生成するステップと、
を有する。
本発明の一実施例によれば、前記（ｙ）ステップは、最低周波数成分を含む薄いスライス成分以外の帯域通過成分に対し行われ、
前記帯域通過成分に対し、前記厚いスライスに組み合わされる薄いスライスの全ての中で、前記薄いスライス成分の最大絶対値はしきい値以下である。
本発明の一実施例によれば、前記源医療画像は、スライス方向に配置された薄いスライスを有する３次元医療画像である。
本発明の一実施例によれば、前記非線形解像度低減手順は、帯域通過分解手順を含む。
本発明の一実施例によれば、前記帯域通過分解手順は、ラプラシアン・ピラミッド手順を含む。
本発明の一実施例によれば、前記帯域通過分解手順は、ウェーブレット変換手順を含む。
本発明の一実施例によれば、前記（ａ）ステップは、
（ａ１）医療画像内で検査中の各ピクセルに対し、サーチピクセルの組を選択するステップと、
（ａ２）前記の各サーチピクセルの近傍の特徴の値と、前記検査中のピクセルの近傍の対応する特徴の値を計算するステップと、
（ａ３）前記の各サーチピクセルの同一又は異なる近傍の他の特徴の値と、前記検査中のピクセルの同一又は異なる近傍の対応する値とを計算するステップと、
（ａ４）検査中の各ピクセルのノイズを低減した状態のグレイ・バリューを、そのサーチピクセルのグレイ・バリューに基づいて、計算するステップと、
を有する。
本発明の一実施例によれば、医療画像装置とコントローラとを有する所定の解像度のノイズ除去された医療画像を提供するシステムにおいて、
前記医療画像装置は、所定の解像度よりも高い解像度で３Ｄ医療画像に対するデータを生成し、
前記コントローラは、高解像度のデータから医療画像のスライスを再構成し、再構成された医療画像をノイズ除去し、前記ノイズ除去した医療画像の解像度を所定の解像度まで低減し、非線形解像度低減手順を用いて厚いスライスを形成するよう薄いスライスを合成し、
前記非線形解像度低減手順は、
（ｘ）前記薄いスライスを帯域通過成分に分解するステップと、
（ｙ）前記帯域通過成分の一部に対し、より大きな重み係数をより大きな振幅を有する薄いスライス成分に与えるようにして、前記薄いスライス成分を厚いスライス成分になるよう組み合わせるステップと、
（ｚ）前記薄いスライスの帯域通過成分を再結合して、所定の解像度の厚いスライスの医療画像を生成するステップと、
を有する。
本発明の一実施例によれば、除去された構造物をノイズ除去した医療画像に復元する方法において、
（ａ）源医療画像をノイズ除去しノイズ除去済み画像を生成するステップと、
（ｂ）前記ノイズ除去済み画像と源医療画像の差に基づいて、残差画像を得るステップと、
（ｃ）非局所スムージング・フィルタを前記残差画像に適用し、フィルタ処理された残差画像を得するステップと、
前記フィルタ処理された残差画像は、前記構造物の程度を前記残差画像内の位置の関数として示し、
（ｄ）前記ノイズ除去済み画像の中の各ボクセルのグレイ・バリューを、所定量だけ変えて、構造物を復元するステップと、
前記所定量は、前記残差画像内の対応する位置にある構造体の程度と、前記残差画像内又はフィルタ処理された残差画像内の対応する位置のグレイ・バリューに基づく、
を有する。
本発明の一実施例によれば、前記（ｄ）ステップは、前記ノイズ除去済み画像に、重み付けされた又は重み付けされていないフィルタ処理された残差画像を付加することを含む。
本発明の一実施例によれば、（ｅ）前記残差画像の中の位置の関数として構造体の程度の測定値を、前記フィルタ処理された残差画像から見出すステップを更に有し、
前記（ｄ）ステップは、前記ノイズ除去済み画像に、前記残差画像又は適用パラメータにより重み付けされたフィルタ処理された残差画像を付加するステップを含み、
前記適用パラメータは、位置の関数として構造の程度の値に局部的に基づく。
本発明の一実施例によれば、前記（ｅ）ステップは、構造体テンソル又はヘシアン・マトリックス（Hessian matrix）の固有値を見出すステップを含む。
本発明の一実施例によれば、前記（ａ）ステップは、
（ａ１）医療画像内で検査中の各ピクセルに対し、サーチピクセルの組を選択するステップと、
（ａ２）前記の各サーチピクセルの近傍の特徴の値と、前記検査中のピクセルの近傍の対応する特徴の値を計算するステップと、
（ａ３）前記の各サーチピクセルの同一又は異なる近傍の他の特徴の値と、前記検査中のピクセルの同一又は異なる近傍の対応する値とを計算するステップと、
（ａ４）検査中の各ピクセルのノイズを低減した状態のグレイ・バリューを、そのサーチピクセルのグレイ・バリューに基づいて、計算するステップと、
を有する。
本発明の一実施例によれば、前記ステップ（ｃ）で使用される非局所スムージング・フィルタ（nonlocal smoothing filter）は、非線形エッジ保存フィルタ（nonlinear edge-preserving filter）を含む。
本発明の一実施例によれば、前記非線形エッジ保存フィルタは、非線形異方性拡散フィルタ（nonlinear anisotropic filter）を含む。
本発明の一実施例によれば、前記非線形異方性拡散フィルタは、ベルトラミ・フロー・フィルタ（Beltrami flow filter）を含む。
本発明の一実施例によれば、前記非線形異方性拡散フィルタは、コヒレンス強化拡散フィルタ（coherence enhancement diffusion filter:CED）を含む。
本発明の一実施例によれば、前記グレイ・バリューの変化は、前記残差画像のグレイ・バリュー又はフィルタ処理された残差画像のグレイ・バリューに対しより感受性があり、前記構造体の程度が大きい場所では構造体の程度が小さい場所よりも、感受性がある。
本発明の一実施例によれば、取り除いた構造体が付加されたノイズ除去した医療画像を生成するシステムにおいて、
医療画像装置と、コントローラとを有し、
前記医療画像装置は、医療画像用の生データを獲得し、
前記コントローラは、前記生データから医療画像を再構成し、前記医療画像をノイズ除去し、前記ノイズ除去済み画像と源画像の差に基づいて残差画像を得て、ノンローカル・スムージング・フィルタを、前記残差画像に適用し、フィルタ処理された残差画像を得て、前記フィルタ処理された残差画像は、構造物の程度を前記残差画像内の位置の関数として示し、前記ノイズ除去済み画像の中の各ボクセルのグレイ・バリューを、所定量だけ変えて、構造物を復元し、
前記所定量は、前記残差画像内の対応する位置にある構造体の程度と、前記残差画像内又はフィルタ処理された残差画像の対応する位置のグレイ・バリューに基づく。
本発明の一実施例によれば、所定の振幅と分布を有するノイズを有する医療画像を生成する方法において、
（ａ）源医療画像からノイズを除去し、ノイズ除去した初期画像を生成するステップと、
（ｂ）前記除去したノイズの一部を、前記ノイズ除去した初期画像に付加するステップと、
を有し、ノイズの所定の振幅と空間分布を具備するノイズ除去済み画像を生成する
を有する。
本発明の一実施例によれば、前記（ｂ）ステップは、
（ｂ１）ノイズの所定の振幅と分布を達成する、空間的に変動する重み付けパラメータを選択するステップと、
（ｂ２）前記重み付けパラメータに基づいて、源画像とノイズ除去した初期画像の重み付けの和を生成するステップと
を有する。
本発明の一実施例によれば、前記（ｂ１）ステップは、源画像とノイズ除去した初期画像の局所標準偏差を計算するステップを含む。
本発明の一実施例によれば、前記（ａ）ステップは、取り除いた構造体を前記ノイズ除去した初期画像に復元するステップを含む。
本発明の一実施例によれば、前記（ｂ）ステップは、空間的に変動するパラメータで乗算した源画像を付加するステップを含み、前記パラメータは、ノイズの所定の振幅と分布に依存する。
本発明の一実施例によれば、医療画像装置とコントローラとを有する所定の振幅と分散を有するノイズを有する医療画像を生成するシステムにおいて、
前記医療画像装置は、医療画像用の生データを獲得し
前記コントローラは、前記生データから医療画像を再構成し、前記医療画像をノイズ除去し、除去したノイズの一部をノイズ除去した初期画像に付加し、ノイズの所定の振幅と分布を有するノイズ除去した医療画像を生成する。
本発明の一実施例によれば、医療画像をノイズ除去する方法において、
（ａ）画像内の所定の複数のボクセルに対し、前記ボクセルを包囲する大きな近傍内のノイズレベルを見積もるステップと、
（ｂ）選択されたボクセルに対する見積もられた前記ノイズレベルから、前記画像内の場所の関数として、ノイズレベルを見出すステップと、
（ｃ）前記画像を、ノイズレベルを場所の関数として用いるノイズ除去アルゴリズムを用いて、ノイズ除去するステップと、
を有し、前記（ａ）ステップは、
（ｉ）グレイレベルの所定の範囲内に入る、前記の大きな近傍内の一部のボクセルに対し、前記ボクセルの周りの小さな近傍内のグレイレベルの変動値を見出し、
（ｉｉ）変動値が見いだされたボクセル内で、前記グレイレベルの変動値が所定のクオンタイル以下であるボクセルのサブセットを見出すステップと、
（ｉｉｉ）前記ノイズレベルを、前記サブセットのボクセルに対する変化の程度の平均値に設定する。
本発明の一実施例によれば、既知のノイズレベルを医療画像の位置の関数として用いて、医療画像をノイズ除去する方法において、
前記医療画像の各ボクセルに対し、
（ａ）前記ボクセルに対する比較ボクセルの組を選択するステップと、
（ｂ）前記ボクセルに関連する近傍と前記比較ボクセルに関連する近傍との間の類似度を計算するステップと、
（ｃ）前記比較ボクセルのグレイレベルの重み付け平均値を見出すステップと、
前記各比較ボクセルの重み係数は、その近傍の類似度に依存し、前記ボクセルの位置におけるノイズレベルに依存する値又は比較ボクセルの位置におけるノイズレベルのいずれか一方又は両方に正規化され、
（ｄ）前記ボクセルのグレイレベルを、前記比較グレイレベルの重み付け平均値に従って変更するステップとを有する方法。
本発明の一実施例によれば、前記（ｂ）ステップは、
前記ボクセルの近傍の特徴値を計算するステップと、前記比較ボクセルの近傍の対応する特徴値を計算するステップとを含む。
本発明の一実施例によれば、前記（ｃ）ステップは、前記医療画像をノイズ除去するステップを含む。

本発明の一実施例によれば、医療画像をシャープにする方法は、
（ａ）医療画像のグレイバリューの局所変動を、医療画像の位置の関数として、測定するステップと、
（ｂ）シャープニング・フィルターを前記医療画像に適用するステップと
を有する。
このフィルターにより生成されたシャープさの程度は、測定された変動が大きい場所では、小さくなる。

本発明の一実施例による医療画像のノイズを除去する方法のフローチャート図。 図１Ａの方法により、医療画像を獲得し、医療画像のノイズを除去するシステムを表す図。 図１Ａの方法を用いて、医療画像の前処理に使用される帯域圧縮関数（range compression function）のプロット図。 ノイズのない二次元画像を表す図。 ノイズが付加され選択されたピクセルと近傍を有する図３Ａと同一の画像を表す図。 図３Ａと３Ｂは、図１Ａの方法で使用された医療画像内のノイズを除去するアルゴリズムをの効果の差を表す。 図３Ａの画像と、選択されたピクセルと、近傍と、図１Ａの方法で使用される画像中のノイズを低減するアルゴリズムに従って選択されたピクセルに類似する他のピクセルも表す図。 図３Ｂの画像と、選択されたピクセルと、近傍を表す図。図１Ａの方法で使用される画像中のノイズを低減するアルゴリズムに従って選択されたピクセルに類似する他のピクセルも表す図。 図１Ａの方法を用いて医療画像内のノイズを低減する方法のフローチャート図。 低Ｘ線ドーズ量で撮影したノイズを含むＣＴ画像を表す図。 図５の方法を用いてノイズを低減した図６Ａの画像を表す図。 図６Ａの画像に類似するが、高Ｘ線ドーズ量で撮影した低ノイズのＣＴ画像を表す図。

本発明は、医療画像の処理を行う方法とシステムに関し、特に、ノイズを低減し更に他の好ましい特徴を有する医療画像を生成する方法とシステムに関する。更に本発明は、医療画像内にノイズが不均一に分布することを考慮して、鮮明な詳細画像を得る非線形フィルタを用いるＣＴ画像の処理方法に関する。

本発明の一実施例は、医療画像のノイズを除去する方法に関する。この方法においては、ノイズ除去は、「高（空間）解像度」で実行される。この「高解像度」の源画像はＣＴにより得られる。その後ノイズ除去済み画像が、「低解像度」に、例えば非線形のダウンサンプリング・アルゴリズムで、変換される。「低（空間）解像度」とは、放射線医師（医師）が通常ＣＴ画像を検査する時の解像度である。この非線形アルゴリズムは、線形ダウンサンプリングアルゴリズムよりも高精細の画像を保持する。従来の画像のノイズ除去は、低解像度に変換された後行われていた。

本発明の一態様は、ある構造（some structure）をノイズ除去済み画像に戻す方法に関する。この方法においては、ノイズ除去のアルゴリズムは、ノイズを除去することに加えて、ある構造を画像から除去する。残差画像（即ち源画像とノイズ除去済み画像との差）は、フィルタ処理してノイズを除去し、エッジ（境界）を強化して、その全部あるいは一部が、ノイズ除去済み画像に戻される。本発明の一実施例においては、画像の各場所において、残差画像用に用いられる相対的な重み付けは、残差画像内の「構造の程度」（degree of structure）に、或いはその場所における源画像又はノイズ除去済み画像内の「構造の程度」に依存する。例えば、構造の程度が低い場所においては、残差画像には、殆ど重み付けは与えられないが、構造の程度が高い場所においては、残差画像には相当量の重みが与えられ、残差画像には、最大の重み付けが与えられる。「最大の重み付け」は、構造の程度が更に増加しても、全く変わらない。選択的事項として、残差画像に用いられる相対的重み付けは、源画像又は残差画像内の各場所におけるノイズのローカルレベルに依存する。ノイズのレベルが高い場所では、残差画像に対しては小さな重み付けが用いられる。

本発明の一実施例は、所定の振幅と分布のノイズを含むノイズ除去済み画像を生成する方法に関する。例えばノイズの所定の振幅と分布は、ＣＴ画像を見る際に、完全にノイズ除去された画像よりも画像がより自然に見えるように、放射線医師が期待するノイズの振幅と分布である。ノイズの所定の振幅と分布は、ノイズを画像に追加することにより行われ、空間的エンベロープが、ノイズの所望の振幅と分布を達成する。ノイズは、空間的に変動する重み付けパラメータを用いて、源画像とノイズ除去済み画像を平均化することで、付加される。

本発明の一実施例は、画像のノイズを除去する方法に関する。この方法においては、ノイズ除去アルゴリズムは、ノイズレベルの測定値を画像内の位置の関数として用いる。このノイズレベルは、所定のボクセル内に見出されるが、これはボクセルの周りの大きなウィンドウを見ることにより、行われる。この大きなウィンドウ内では、ある範囲のグレイ・バリュー内のボクセルのみが考慮される。この範囲は、画像化される組織に依存する。例えば、その組織内の中間範囲内のボクセルのみが考慮される。各ボクセルに対し、グレイ・バリュー内の変化の局所測定値、例えば標準偏差が、小さなウィンドウ内に見出される。変動のこれらの測定値のサブセットを取り出し、ある分位（quantile）以下の変動の測定値のサブセットと、サブセット内の変動の平均値を取り出して、大きなウィンドウに対する局所ノイズレベルを見出す。画像内の位置の関数としてのノイズレベルは、画像内に分散する複数のボクセルに対し、この手順を用いて見出すことができ、選択的事項として、補間手法を用いて画像全体をカバーする。画像内でノイズレベルが、位置の関数として分かると、それを他の手順と一緒に用いて、局所ノイズレベルに応じて、手順の性能を局部的に最適化したり改善したりする。

本発明の一実施例は、画像をノイズ除去する（画像のノイズを除去する）方法に関する。この方法においては、画像の位置の関数としてのノイズレベルを、ノイズ除去のアルゴリズムで用いる。ノイズ除去のアルゴリズムの一例は、アルゴリズムの非局所平均（nonlocal means:：ＮＬＭ）である。このアルゴリズムにおいては、所定のボクセルの周りのパッチと画像内の何処かの比較用ボクセルの周りのパッチとを比較し、その類似性測定値（similarity measure）を、２つのパッチに対し計算する。ノイズに起因する類似性測定値における予測変動値が計算されるが、これは、所定のボクセルの局所ノイズレベル、比較ボクセルの局所ノイズレベル、のいずれかあるいはその両方に依存して計算され、これを用いて類似性を正規化する。
従来のＮＬＭアルゴリズムにおいては、ボクセルの位置から独立した（無関係な）一定のノイズレベルを用いて、類似性を正規化している。比較用のボクセルの重み付けが、この正規化された類似性に基づいて見出され、比較ボクセルのグレイ・バリューの重み付けされた平均値が見出される。その後、所定のボクセルのグレイ・バリューを、比較ボクセルの平均グレイ・バリューに基づいて、変更する。位置の関数としてのノイズレベルは、上記の手順により見出される。

本発明の一実施例は、画像をシャープにする方法に関する。この方法においては、グレイ・バリュー内の局所変化の程度、例えばグレイ・バリュー内の局所標準偏差を、画像の位置の関数として見出す。sharpening filterをこの画像に適用するが、シャープにするレベルは、グレイ・バリュー内の局部変動のレベルに依存する。グレイ・バリュー内の局部変動のレベルが大きい時は、シャープにするレベルを低くする。その結果得られた画像は、シャープなエッジにおいては視覚的なアーチファクトが少ないか全くない。これは、線形のシャープニングフィルタが、所定の平均レベルのシャープさに対し、視覚的なアーチファクトを多く生成するのとは、対照的である。

一般的に、画像処理においては、用語「ピクセル」は二次元画像の要素に、用語「ボクセル」は三次元画像の要素に用いる。本発明は、二次元画像又は三次元画像のいずれかに用いられるために、用語「ピクセル」と用語「ボクセル」は、二次元又は三次元の画像に限定されるものではない。特に断りのない限り、ピクセルとボクセルは、一般的な用語として互換性を持って使用する。

医療画像の測定装置、例えばＣＴ又はＭＲＩにおいては、画像を表示する時は、画像密度（例、ＣＴ画像内のHounsfield Units(ＨＵ)の密度）を、特定の密度の範囲に対し、黒から白の範囲に渡る輝度レベル或いはグレイ・バリューにマッピングすることが一般的である。用語「グレイ・バリュー」は、画像の密度を意味し、ディスプレイ内の画像の密度を表す輝度レベル（bright level）を意味しない。この用語の使用は、画像の密度が、黒と白の間の輝度レベルにマッピングされる範囲外にある場合でも用いられる。用語「グレイ・バリュー」は、画像装置（方法）の画像密度の範囲外にある画像の密度も意味する。この画像の密度は、例えば負数を画像処理する中間ステップで生成される。更に用語「グレイ・バリュー」は、黒と白の画像の輝度レベルのみならず、カラー画像の色の変動（例、赤、緑、青）のレベル、更には強度又はカラー画像の輝度レベルと飽和レベルを意味する。ＣＴ又はＭＲＩの画像の様な医療画像においては、１個の密度変数のみが表示される。例えば、ＨＵ、又はＴ_１又はＴ_２の重み付け密度が表示され、これらはグレイスケールのディスプレイの輝度レベルにマッピングされる。この場合、「グレイ・バリュー」が特に使用される。しかし、本発明の方法は、医療画像あるいは白黒の画像に限定されない。

本発明のノイズ低減方法は、医療画像にとって特に有効である。その理由は、医療画像は比較的高いノイズレベルを有する傾向があるからである。ノイズレベルと画像獲得パラメータ（例、Ｘ線のドーズ量、ＭＲＩ獲得時間）の間にはトレードオフが存在し、ノイズレベルを低減すると、経済的あるいは安全性の問題が発生する。医療画像は、「照射（lighting）」の差異がないため、ピクセルの近傍の特徴（feature of neighborhood）は、真のグレイ・バリューを正しく表すことになる。このことは、医療画像は画像の様々な部位に対し繰り返される類似の構造を有し、これらの構造は、スケールが変わったりあるいは方向が変わったりするだけである。
本発明の方法は、イメージ獲得装置あるいはワークステーション（例ＣＴ、ＭＲＩ、超音波画像装置）に適用でき、更にイメージ作成ステーションあるいは離れた場所へのネットワーク接続にも適用できる。

ノイズ除去の概要
図１Ａに、本発明の一実施例により、ノイズを低減した画像を生成する手順を示す。この手順は、ＣＴ画像を例に説明するが、他の医療画像あるいは他のタイプの画像にも適用可能である。様々なタイプの医療画像更には医療画像でない通常画像は、様々な特徴を有し、この手順が適用可能である。医療画像（例、ＣＴ画像とＭＲＩ画像）は、通常の医療用でない画像が使用する照射には依存せず、ノイズのレベルと空間分布は、一定であり、あるイメージから別のイメージへの予測が可能である。

ステップ１００２において、高解像度でノイズを含む３Ｄ源画像Ｉが、医療用画像装置（例、ＣＴスキャナー）で得られる。画像を得てノイズ除去するシステム５００を、図１Ｂに示す。この画像は、画像装置５０２で患者５０４から得られた生データ（raw data）からコントローラ５０６を用いて、再構成される。このコントローラ５０６は、画像装置のコンピュータあるいは専用回路である。コントローラは、以下に説明するノイズ除去手順の一部あるいは全てを実行する。本発明の一実施例においては、個別のコンピュータあるいは回路が、ノイズ除去又は画像の再構成プロセスの様々な部分を実行するが、これらは全てコントローラと総称する。出力装置５０８は、例えばノイズ除去済み画像を表示するディスプレイモニターであり、入力装置５１０は、キーボード又はコンソールであり、これを用いて画像の処理又は画像化装置を制御する。

図１Ａにおいて、ステップ１００４において源画像Ｉは、非線形フィルタを用いてフィルタ処理し、その結果前処理済みの画像Ｃとなる。この前処理済みの画像Ｃは、ノイズ除去手順でパッチを比較するのに用いられる。通常のＣＴ画像において、最も興味ある点の密度範囲は、通常小さく例えば−５０ＨＵから＋１００ＨＵである。ただしフルレンジは−１０００ＨＵから＋４０００ＨＵである。選択的事項としての範囲圧縮手順は、源画像Ｉを非線形のフィルタ処理する前に、実行される。予測されるノイズモデルに従った画像の「プレ−ホワイトニング（pre-whitening）」も、選択的に、実行される。これに関しては、段落００３５で説明する。

ＣＴ画像においては、ノイズは空間的に均一ではない。その理由は、獲得機構とビーム・ハードニング（beam hardening ）が原因であり、このプロセスでは、Ｘ線の非単一エネルギービーム（non-monoenergetic beam）の平均エネルギーは、身体を通るに連れて増加するが、それは軟Ｘ線が吸収されやすいからである。ノイズ除去手順は、ＣＴ画像の各ポイントにおける局所ノイズ特性による空間依存性のあるパラメータを使用する。局所ノイズは、ステップ１００６で予測される。局所ノイズ特性は、画像の前処理あるいはそれに並行して見積もられる。一方で局所ノイズ特性は、前処理された画像の白色化されたノイズ(whitened noise)に基づき、源画像の着色されたノイズ(colored noise)に基づかない場合には、最初に前処理をするのが好ましい。これらの用語については、段落００３５で説明する。

ノイズ除去手順の局所平均（ＮＬＭ）は、ステップ１００８で、画像に対し実行される。一般的に、ＮＬＭ手順は、画像内の所定のボクセルをノイズ除去するが、これは所定のボクセルの周りのパッチに類似するパッチを見出すことにより、行われる。ボクセルのノイズ除去された値は、パッチの中心ボクセルの重み付けされた平均値として計算される。この重み付けは、パッチ間のある類似メトリック（similarity metric）に比例する。他のノイズ除去手順は、従来公知である。ＮＬＭ手順以外の手順も、ＮＬＭノイズ除去手順の代わりにあるいはそれと共に使用することができる。

フィーチャーベースのＮＬＭ手順においては、類似性基準は、パッチから抽出された様々な特徴に基づく。この特徴は、前処理済みの画像Ｃと源画像Ｉの両方に基づく。類似性基準は、ステップ１００６で見積もられたパラメータにも依存する。Ｄは、特徴ベースのＮＬＭ手順を実行した後のノイズ除去済み画像である。ＮＬＭ手順は、ステップ１００８で用いられるが、フィーチャーベースのＮＬＭ手順、あるいは上記文献に開示されたＮＬＭ手順のいずれかである。フィーチャーベースのＮＬＭ手順は、段落００４５で説明する手順であり、例えばＷＯ２００９／０８１４１０に記載したフィーチャーベースのＮＬＭ手順である。

ステップ１０１０において、ノイズ除去手順で画像から除去された構造は、画像に戻される。ノイズ除去ステップの後でも、残差画像（Ｒ＝Ｉ−Ｄ）からの更なる特徴を戻すこともできる。この残差画像は、除去したノイズの大部分を含むが、ある空間的に識別可能な構造も含む。Coherence Enhancement Diffusion（ＣＥＤ）フィルタは、残差画像の非線形平滑化（non-linear smoothing）を実行するが、これを用いてこれらの構造を見出す。選択的事項として、フィルタ処理された残差画像Ｓは、直接ノイズ除去済み画像Ｄには加えられず、パラメータαで乗算されて加えられる。このパラメータαは、Ｄの局部特性に依存する。Ｄ’＝Ｄ＋αＳは、回復後のノイズ除去済み画像である。

選択的事項として、ノイズが画像に加えられる（ステップ１０１２）。ノイズ除去済み画像は、ＣＴ画像で所定のレベルのノイズに慣れた放射線医師にとっては、不自然に見える。回復後のノイズ除去済み画像Ｄ’と源画像Ｉとの平均が、放射線医師にとって自然に見える。この様な画像の一例は、Ｅ＝βＤ’＋（１−β）Ｉである。ここでβは、得られた画像Ｅが所定の標準偏差のノイズパターンを有するよう、選択される。βは、画像が変わると変化する。選択的事項として、画像Ｄ’もセーブされ、更に後で使用可能となる。この場合、この追加されたノイズは、放射線医師が高解像度で画像の一部を拡大（ズーム）しようとする場合には、不都合である。

画像の解像度は、ステップ１０１４で減らされる。放射線医師は、低解像度の通常のＣＴ画像を検査し、ＣＴスキャナーにより得られた生の画像の解像度を検査しない。解像度を減らすことは、非線形解像度低減手順を用いて行われる。この解像度低減手順は、画像の線形フィルタリングとサンプリングよりも、より詳細な事項を保存する。

ノイズ除去方法の詳細
フローチャート１００の詳細を以下説明する。
画像の獲得
日常の診察においては、放射線医師が検査するＣＴ画像は、約２．５ｍｍの軸方向解像度（スライス厚さ）を有する。最新型のＣＴスキャナーの生の軸方向解像度は、約０．６ｍｍである。その為ＣＴスキャナーにより得られた解像度は、臨床医が見る解像度よりも遙かに高い。高解像度のデータは日常的には利用されない。その理由は、その画像は高いノイズレベルを有し、読み込むのに時間がかかるからである。臨床で見る各低解像度のスライスは、源画像内の複数のスライスから集められたデータを含む。

本発明のノイズ除去手順とそれに関連する画像手順は、高解像度のデータを用いるが、高解像度のデータは、ＣＴスキャナーで得られ、より多くの情報を含む。処理の完了後、低解像度の画像が形成され、放射線医師がそれを見る。

本発明の一実施例においては、臨床の軸方向解像度とオリジナルの軸方向解像度の差は、２つの方法で得られる。第１の方法は、データが豊富なセットを、画像処理アルゴリズムに入力する。高い情報内容のために、このアルゴリズムは、より正確な結果を生成できる。第２の方法は、処理された画像のダウンサンプリングが、非線形に実行され、低解像度の画像を生成し、源画像で得られた詳細情報を保存したままにしておく。

前処理された画像
前処理段階において、前処理された画像Ｃは、源画像Ｉから計算される。この前処理された画像を用いて、ノイズ除去手順でのパッチの比較を行う。基本的なノイズ除去手順は次式で記述される。

ここで左辺（ｘハット）は、ｉ番目のボクセルのノイズ除去された値であり、重みｗ（ｘ_ｉ，ｘ_ｊ）は、源画像Ｉのボクセルｘ_ｉの周りのパッチと、ボクセルｘ_ｊの周りのパッチの間のフィーチャーベースの距離の関数である。より一般的なノイズ除去アルゴリズムは、源画像Ｉのパッチと前処理された画像Ｃのパッチとを比較する。

前処理された画像Ｃは、源画像Ｉから次のステップを用いて構築される。

Ａ． ノイズのモデルによる色のついたノイズ（colored noise）のプレホワイトニング。
ＣＴノイズは、付属性着色ガウスノイズ（aditive colored Gaussian noise）即ちホワイトガウスノイズとして、モデル化され、フィルタによりフィルタ処理される。このフィルタは、ノイズカラーリングフィルタとも称する。画像をプレホワイトニングすることは、Wienerフィルタを適用することにより行われる。このWienerフィルタは、ノイズカラーリングの操作を逆にし、ホワイトノイズを含む画像を生成する。ホワイトノイズは、空間的に相関関係がない。この種のノイズは、様々な文献で広範囲に研究されており、色のついたノイズよりも容易に除去できる。
ノイズカラーリングフィルタは、均一のファントム（uniform phantom：体内に照射される放射線の線量を予測するコンピュータ制御によるモデル）のＣＴイメージからの共分散マトリックスを計算することにより、見積もられる。これは、非特許文献６に記載されている方法に類似する。

Ｂ． レンジ圧縮
ＣＴ画像のダイナミックレンジは、通常−１０００ＨＵ（空気）と、＋４０００ＨＵ（金属）戸の間にある。０ＨＵと１００ＨＵのレンジは、１０００ＨＵと１５００ＨＵのレンジよりもより重要である。その理由は、第１レンジは軟組織を表し、第２レンジは遙かに密度の高い組織（骨）を表すからである。更に骨のような硬い硬組織は、通常密度の遙かに幅の広いウィンドウを用いて検査され、これによりノイズは見えなくなる。
放射線医師は、０ＨＵから＋８０ＨＵの密度ウィンドウでもって、脳のＣＴ画像を検査し、肝臓のＣＴ画像を−１５ＨＵから＋１５５ＨＵの密度ウィンドウで検査する。画像Ｃに対しグレイレベルの変換（grey-level transformation）を用いるのが好ましい。これは軟組織の範囲に伸び、−１０００ＨＵから−２００ＨＵの範囲、そして更には＋３００ＨＵの範囲を圧縮する。例えばグレイレベルの変換は、源密度のｘの次の関数ｙである。

ａ，ｂはレンジ定数であり、ｙは変換された密度である。図２はレンジ圧縮関数ｙ（ｘ）のプロット２００を示す。
他のグレイレベル変換関数の形成しそれを用いることもできる。例えば適宜のグレイレベルの変換は、ＣＴ画像のヒストグラムを伸ばし、このＣＴ画像のヒストグラムを、選択的事項として、画像の興味ある部位のみを等価することにより得られる。そしてこれは周囲の空気なしで行うことができる。

Ｃ． 頑強な非線形ノイズ除去
ノイズ除去フィルタを画像Ｃに適用する。本発明によればノイズ除去フィルタが画像Ｃを過度に視覚的にスムーズにした場合でも、かつ解像度が幾分失われた場合でも、フィルタ処理済画像Ｃは、パッチを比較するためには、フィルタ未処理画像Ｃよりも良好である。画像Ｃは、類似性基準を評価するために、パッチを比較する目的のためだけに使用され、平均化操作は、源画像Ｉに対し実行され、フィルタ処理された画像Ｃに対しては実行されない。
画像Ｃは、双方向フィルタを用いてフィルタ処理される。これは非特許文献７に記載されている。更に画像Ｃは、Coherence Enhancement Diffusion(CED)filterを用いてフィルタ処理される。これに関しては非特許文献８に開示されている。これは画像を平滑にし、かつエッジを保存する。他の公知のノイズ除去フィルタも使用可能である。

局所ノイズレベルの見積もり
非局所平均のノイズ除去スキーム（non-local means denoising scheme）の重要なパラメータは、δ_Ｒである。これは、複数の近傍の間の重み付けを制御し、次式で与えられる。

ここでＣ_Ｘは、画像Ｃ内のボクセルｘの周囲の画像パッチを表し、ｄｐ（ｉ，ｊ）は、ボクセルｉとボクセルｊの間の空間距離である。Ｃ_ｉ−Ｃ_ｊの絶対値は、画像Ｃ内のボクセルｉとボクセルｊの周りのパッチの間の差の測定値である。

δ_Ｒの値は、ノイズの起因する測定値の変動の標準偏差と同じ桁である。これに関しては、非特許文献１３を参照のこと。この重み係数の選択は、以下の機能的観測から得られる。画像内の統計的ノイズの大きさにまで類似する２個の画像パッチは、大きな重み係数を採り、類似性を示す。これに対し、画像内のノイズから予測される差を超えた非類似性を表す２つのパッチは、低い重み係数を採る。

従来、画像内のノイズは同一に分布すると考えられていた。即ち画像内のあらゆる空間場所で類似の統計的特徴を有すると考えられていた。その為δ_Ｒは、画像ごとに一旦決まると、同一の値が、画像内の全ての空間位置に対し重みを計算するのに用いられていた。

しかし、ＣＴ画像内においてノイズは空間的に均一に分布していない。その原因は、獲得機構とトモグラフィ再構成アルゴリズムにある。本発明者は以下のことを見出した。画像全体でδ_Ｒを一定値にすることは、ＣＴ画像に対し十分機能しないこと。δ_Ｒを一定値にすることにより、ノイズのない画像領域が過剰に平滑(over-smoothing)になり、他方で、高いノイズレベルを有する画像領域のノイズ除去は不十分であること。代わりに、空間的に変化する値をδ_Ｒに用いること。例えば、ノイズに起因する測定値の局所標準偏差に従って、画像内の各ボクセルをδ_Ｒ用に適宜の値を割り当てること。

ノイズに起因する空間依存性のある標準偏差を見積もる方法を次に説明する。
まず、局所標準偏差あるいはグレイ・バリューの局所変動の別の測定値を、ノイズを含んだ源画像Ｉの各ボクセルに対し計算する。これは、考慮中のボクセルに関連する小さな近傍（例、３×３×３のボクセル）の標準偏差（又は他の偏差値）を計算することにより行われる。この近傍は、その中心に考慮中のボクセルを含む。
次に、遙かに大きな近傍（例、３３×３３×５のボクセル）を、画像内の各ボクセルの周りで検査する。各大きな近傍に対しては、ある密度範囲（例、−５０ＨＵから＋１５０ＨＵ）内にあるボクセルの一部あるいは全てが、抽出される。
最後に、抽出されたボクセルのサブサンプルの局所標準偏差の値は、大きなウィンドウ内にある抽出されたボクセルの局所標準偏差のある分位（例、０．３）よりも小さい局所標準偏差に対応するが、これを平均化して、ノイズに起因する局所標準偏差の予測値を生み出す。

大きな近傍の使用により、ノイズに起因する局所標準偏差の予測値は、隣接するボクセルの間で大きく変わることはないと思われる。選択的事項として、計算は、全ての画像のボクセルのサブサンプルに対してのみ行われ、他の場所にある予測値は、必要な時に、標準の補間技術で、計算された値を補間ことにより、得られる。
更に他の方法を用いて、局所ノイズレベルを予測することもできる。例えばファントムを画像化し、ファントムの均一な部分にあるグレーレベルの標準偏差から直接ノイズレベルを決定することもできる。
δ_Ｒの値は、ノイズに起因する局所標準偏差の予測値の一部（例、０．５４）に設定される。他のパラメータも、ノイズに起因する局所標準偏差の予測値に基づいて、ノイズ除去アルゴリズムあるいは他の関連するイメージ処理手順について用いることができる。
本発明の一実施例においては、ノイズ除去のアルゴリズムあるいは関連する手順で用いられるδ_Ｒ又は他のパラメータは、画像内の何処でも同じ値を有するが、ノイズに起因する局所標準偏差の平均値に依存する。これらは上記したように計算される。例えば、ある密度範囲内にあるボクセルのみが局所標準偏差を探し出すために考慮に入れられる。選択的事項として、ある値以下の局所標準偏差のみが、平均を計算する際に含まれる。

代表的なＮＬＭノイズ除去アルゴリズム
図３Ａは、ピクセルのアレイを含む二次元画像１００を示す。各ピクセルの数字は白と黒の間のグレイ・バリューにマッピングされている。ＣＴ画像内のグレイ・バリューは、画像化される対象物の実際の密度のマッピングを表す、単位はＨＵである。脳のＣＴ画像においては、画像は０（ゼロ）ＨＵは水の密度を表し、これは黒にマッピングされ７０ＨＵは白にマッピングされるようビジュアル化される。
画像１００は、明るい領域１０２と、暗い領域１０４と、その間の明確なシャープな境界線が示されている。

図３Ｂにおいて、画像１０８はノイズを付加した画像１００である。従来は、ノイズは、ピクセルのグレイ・バリューと近隣のピクセルのグレイ・バリューとを平均することにより、減らしていた。この場合、最近ピクセルに最大の重み（係数）を与えていた。これは微細な詳細領域のない均一の領域（例、画像１０８内の明るい領域１０２と暗い領域１０４）では十分機能する。しかしそれらの間の境界はぼけることになる。従来のノイズ低減方法は、バイラテラル・フィルタ（bilateral filter）を使用し、これは非線形フィルタであるが、ピクセルｉのグレイ・バリューＩ_ｉを、グレイ・バリューが似ている他のピクセルｊのグレイ・バリューＩ_ｊで平均化することにより、この問題を回避していた。例えば（ｘ_ｉ，ｙ_ｉ）にあるピクセルのｉに対し操作を行うと、重みＷ_ｊが（ｘ_ｉ，ｙ_ｉ）にある別のボクセルｊのグレイ・バリューに対し次式で与えられる。

ここでｄ_ｐは、２つのピクセルの間のユークリッド距離であり、（Ｉ_ｉ−Ｉ_ｊ）は、２つのピクセルの概算的な距離であり、それらが互いに似ている程度を表す。ピクセルのｉの新たなグレイ・バリューは次式で表される。

ここでＮは、ピクセルｉの周りのサーチウィンドウであり、加算は、このサーチウィンドウの中の全てのピクセルｊに対し行われる。

ノイズ低減に用いられる別の種類の非線形フィルタは、非特許文献９に記載されている。非局所平均フィルタにおいては、２つのピクセルの類似度は、この２つの近傍のピクセルごとの比較に基づく。
例えば、図３Ｂで番号１１０のピクセルｉのノイズレベルを低減するために、図３Ｂの番号１１２の近傍Ｍｉが、ピクセル１１０の周りであると定義される。
その後、他のピクセルｊが検索され、各サーチピクセルｊの周囲の同一サイズと同一形状の近傍Ｍ_ｊが検索され、平均自乗誤差ＭＳＥ（Ｍ_ｉ，Ｍ_ｊ）が、近傍１１２のピクセルと、サーチピクセルｊの近傍の対応するピクセルとの間で、見出される。ピクセル１１０の近傍の間で、平均自乗誤差が小さなサーチピクセルには、最も重い重み係数が与えられるが、これは、ピクセル１１０の低減したノイズのグレーレベルを得るために、サーチピクセルのグレイ・バリューを平均化する時である。
重み付けＷ_ｊは次式で与えられる。

ピクセルｉの新たな値は次式で与えられる。

図４Ａは、画像２００を示す。これは、図３Ｂの画像１０８と類似する。ピクセルの組２０２は、近隣近傍１１２又はピクセル１１０に似ている近傍を有する。各ピクセル２０２は、類似の近傍を有する。その理由は、ピクセル２０２は、明るい領域１０２と暗い領域１０４の間の境界から全て同一距離にあり、同一方向を向いているからである。
非線形フィルタを用いる他のノイズ低減方法においては、２個の近傍の類似性は、近傍内にある全てのピクセルの平均グレイ・バリュー、又は近傍内のピクセルのグレイ・バリューの傾斜方向に依存する。これに関しては非特許文献１０に記載されている。別の方法は非特許文献１１に開示されている。この別の方法においては、２つの近傍の類似性は、２つの近傍の二進のエッジマップの平均自乗誤差に依存する。これはCannyエッジディテクタを用いて決定され、かつ２つの近傍の元の源画像の平均自乗誤差に基づく。

図４Ｂは、画像２０４を示す。これは、図３Ｂの画像１０８と類似する。本発明の一実施例により、２つの近傍の間の類似性に対する別の基準を用いて、重みＷ_ｊを計算すると、より良好なサーチピクセル２０６が見出される。このサーチピクセル２０６は、ピクセル１１０の近傍１１２に十分類似する近傍を有する。
図２Ｂに示す実施例においては、以下に説明する基準は、近傍の相対方向に依存せず、暗い領域１０４からピクセル１１０と同一距離にある全てのピクセルは、この基準によれば、近傍１１２に極めて近い近傍を有する。高い重み係数を有するサーチピクセル２０６の拡大した組は、ＮＬＭ方法を用いた高い重み係数を有するサーチピクセル２０２に比較すると、更なるノイズ低減が可能となる。その理由は、グレイ・バリューを平均化するより多くのピクセルがあるからである。
本発明の一実施例よれば、２個の近傍の間の類似性基準は、２個の近傍の相対方向に依存する。
本発明の他の実施例においては、高い重み係数を有するサーチピクセルの数は、従来方法よりは多くはないが、ピクセル１１０の真のグレイ・バリューのよりよい予測値を提供する点から見て、サーチピクセルの質は良くなっている。

図５は、本発明の一実施例による画像のノイズ低減方法のフローチャート３００である。フローチャート３００の方法は、図４Ｂの方法を一般化したものであり、近傍間の類似性に対し様々な基準が示されている。ステップ３０２において、ノイズのある画像を得る。ノイズ低減アルゴリズムは、一時に１個のピクセルを検査する。ステップ３０４でｉ＝１のピクセルの組を検査する。ステップ３０６で、ｉ番目のピクセルが検討され、ステップ３０８で、ｉ番目のピクセルのフィーチャーベクトルＦ_１が見出される。フィーチャーベクトルは、複数のフィーチャーの値の番号順の組であり、それぞれ、考慮中のピクセルのグレイ・バリューあるいは関連する近傍の他のピクセルのグレイ・バリューに依存する。近傍は連続してる必要もなく、又ピクセルを囲む必要もなく、一方の側にある。考慮中のピクセルの座標軸は、（ｘ_ｉ，ｙ_ｉ）あるいは三次元では（ｘ_ｉ，ｙ_ｉ，ｚ_ｉ）であり、これらがフィーチャーとして取り扱われる。従来既知のフィーチャーの一例は、上記のバイラテラル・フィルタで使用されるピクセルｉのグレイ・バリューと、非局所平均フィルタで使用されるピクセルｉの周囲の所定サイズの近傍内の各ピクセルのグレイ・バリューを含む。従来既知の他のフィーチャーは、ピクセルｉの近傍の全てのピクセルの平均グレイ・バリューと、ピクセルｉの近傍のグレイ・バリューの傾斜方向と、ピクセルｉの近傍のバイナリ（二進）エッジマップの各ピクセルのグレイ・バリューを含む。これらはCanny edge detectorを用いて決定される。他にも様々なフィーチャーがあるが、これはＷＯ２００９／０８１４１０に記載されている。

ステップ３１０において、サーチピクセルのｊの組が検査され、ピクセルｉの特徴値に類似する特徴値を有するピクセルを見出す。ピクセルｉに最も類似するサーチピクセルｊのグレイ・バリューは、ピクセルｉのノイズのない真のグレイ・バリューの予測値に最も寄与する。最初にステップ３１０において、インデックスｊは１に設定される。ステップ３１２において、ｊ番目のサーチピクセルが考慮の対象となる。このサーチピクセルは、画像内の全てのピクセルあるいはピクセルｉを除く全てのピクセルを含む。別の構成として、サーチピクセルは、画像内のピクセルのサブセットのみを含む、例えば、ピクセルｉの周りのサーチウィンドウ内に入るピクセルのみを、サーチウィンドウ内でランダムに選択されたあるピクセルのみを、サーチウィンドウ内で通常の間隔にあるピクセルのみを、又はピクセルｉのグレイ・バリューに十分近いグレイ・バリューを有するピクセルのみを含む。選択的事項として、医療画像においては、画像は、公知のセグメンテーション技術を用いて、様々な種類の画像に分割される。サーチピクセルは、ピクセルｉと同じタイプの組織のピクセルからのみから選択される。

本発明の一実施例においては、ノイズ除去済み画像から選択されたサーチピクセルを用いることに加えてあるいはそれに換えて、サーチピクセルを、データベース内の他の画像から選択することもできる。例えば可能なサーチピクセルとその近傍をこの画像に類似すると期待される他の画像から得るディクショナリ（辞書）を予め編集することもできる。例えば画像が医療画像の場合には、ディクショナリ（辞書）は、所定のサーチピクセルとその近傍を有するが、これらは、同一の患者の身体の同一部分からの以前に取った医療画像からあるいは他の患者の身体の同一部分から得られたものである。

ステップ３１４において、特徴ベクトル（feature vector）Ｆ_２が各ピクセルｊに対し評価される。特徴ベクトルＦ_２は、特徴値の順番の付いた組であり、これは特徴ベクトルＦ_１の特徴値の１つに対応する。選択的事項として、Ｆ_１とＦ_２内の対応する特徴値は、ピクセルｉとピクセルｊの近傍内の対応するピクセルのグレイ・バリューを用いて、同じ方法で規定することもできる。本発明の一実施例においては、Ｆ_１とＦ_２内の対応する特徴値は別の方法でも定義できる。例えば、あるピクセルの周囲の近傍を他のピクセルの周囲の近傍に対し、異なる角度で方向付ける、あるいはサイズ換算することにより、特徴値を計算する際に必要によっては補間されたグレイ・バリューを用いることもできる。いずれの場合においても、Ｆ_１とＦ_２の対応する特徴値は、それらを比較するのに有意な方法で、あるいはピクセルｉとピクセルｊの間の概算距離測定値を計算する為に、その値の差を用いて行われる。これによりピクセルがどの程度似ているかを測定して、ノイズ低減に役立てる。

サーチピクセルｊは、検査中の画像からではなく、サーチピクセルの以前に記憶されたディクショナリ（辞書）から取り出されると、ピクセルｊに対する特徴ベクトルＦ_２又はその構成要素の一部は、ディクショナリに記憶され、それを使う毎に計算する必要はなくなる。同様にサーチピクセルｊが、別のピクセルｉのサーチピクセルとして以前に使用された場合には、その特徴ベクトルＦ_２は、メモリー内に記憶され、再度計算されることはない。選択的事項として、特徴ベクトルＦ_２を画像内の全てのピクセルに対し予め評価しておき、メモリー内に記憶しておくと、Ｆ_２は、サーチピクセルｊとピクセルｉに対しループを行う間、評価する必要はない。Ｆ_２の特徴値が、Ｆ_１の対応する特徴値と同様に規定される程度において、特徴ベクトルＦ_２あるいはその成分の一部は、メモリーから取り出され、再び計算されることはない。これはサーチピクセルｊが、ステップ３０６で検査されているピクセルｉとして以前に用いられた場合である。

ステップ３１６において、距離の測定値ｄ（Ｆ１，Ｆ２）が計算されるが、この値は、ピクセルｊがピクセルｉに似ているかを反映する概算距離であり、これは近傍のグレイ・バリューにより規定され、更にそれらの位置にも規定される。距離の測定値ｄは、特徴ベクトルＦ１，Ｆ２を形成する対応する特徴の値の差に依存する。特徴ベクトルＦ１とＦ２のｋ番目の成分（特徴値）が、Ｆ１（ｆ_１ ^１，ｆ_２ ^１−−−ｆ_ｋ ^１），Ｆ２（ｆ_１ ^２，ｆ_２ ^２−−−ｆ_ｋ ^２）で与えられると、距離の測定値は次式で定義される。

ここで、（α_１，α_２，…，α_ｋ）は、距離を計算する際に異なる特徴に対し使用される重みを与える重みベクトルである。パラメータβは、単位ユニットの正の値であるが、時に、２に等しく、この場合、ｄ（Ｆ_１，Ｆ_２）は直交成分のユークリッド距離になり、ピクセルｉとピクセルｊの特徴値の間の絶対距離に等しい。重みベクトル（α_１，α_２，…，α_ｋ）は、一般的なアルゴリズムを用いて見出される。このアルゴリズムは、ノイズ低減方法の効率を最大にする最適の重み付けベクトルを見出すアルゴリズムである。

ｄ（Ｆ_１，Ｆ_２）の別の表現は、異なる特徴値の間の相関を考慮に入れているが、例えば近傍の異なるピクセルのグレイ・バリューの間の相関を考慮に入れているが、ＷＯ２００９／０８１４１０号公報に開示されている。ｄ（Ｆ_１，Ｆ_２）は、クロスターム例えば、（ｆ_１ ^１−ｆ_１ ^２）（ｆ_２ ^１−ｆ_２ ^２）を含むが、異なる特徴値に相関関係がある場合には、異なる近傍の間の類似の程度に対する有効な測定値となる。

ステップ３１８において、ピクセルｊに対する重み係数Ｗ_ｊは、ｄ（Ｆ_１，Ｆ_２）から計算され、メモリーに記憶される。ピクセルｉの近傍とピクセルｊの近傍が、最もよく似ている場合には、即ちｄが小さい場合にはＷ_ｊは最大になり、ｄが大きい場合にはＷ_ｊは小さくなる。例えばＷ_ｊ＝ｅｘｐ（−ｄ^２／σ_Ｎ）である。特徴値がピクセルとその近傍のグレイ・バリューにのみ依存し、ピクセルの位置に依存しない場合には、重み係数Ｗ_ｊは、次式で規定される。Ｗ_ｊ＝ｅｘｐ（−ｄ^２／σ_Ｎ−ｄ_ｐ ^２／σ_ｐ）
ｄ_ｐは、ピクセルｉとピクセルｊの間の物理的距離例えばユークリッド距離である。ここσ_Ｎとσ_ｐは、概算距離ｄと、ピクセルｉとピクセルｊの間の空間距離ｄ_ｐが増加する時に、Ｗ_ｊが如何に速く減るかを決定するパラメータである。別の構成として、Ｗ_ｊは、ｄとｄ_ｐへの異なる依存性を有するが、ｄとｄ_ｐが大きな値になると減る。選択的事項として、計算時間を節約するためあるいは性能を上げるために、重み係数Ｗｊは、あるしきい値よりも小さい場合、或いはｄ又はｄ_ｐがあるしきい値より大きい場合には、０に設定される。

ステップ３２０において、サーチピクセルのインデックスｊを１増加させて、次のサーチピクセルを見る。ステップ３２２において、全てのサーチピクセルが検査されたか否かを決定する。検査されていない場合には、次のサーチピクセルがステップ３１２で候補に上げられる。全てのサーチピクセルが検査された時には、サーチピクセルｊのグレイ・バリューのＷ_ｊで重み付けされた平均値が、計算される。

ステップ３２６において、ノイズのないピクセルｉの真のグレイ・バリューを見積もる。この見積もりは、サーチピクセルのグレイ・バリューあるいは元のピクセルのグレイ・バリューに基づく。サーチピクセルは、それらが類似の特徴値に基づいてピクセルｉに類似する場合には、見積もられた真のグレイ・バリューに対し大きな影響を有する。例えば特徴値の類似性を用いて概算距離ｄ（Ｆ_１，Ｆ_２）を計算すると、各サーチピクセルｊは、ピクセルｉからの距離に基づいて、重み付け係数Ｗｊが割り当てられ、ピクセルｉの見積もられた真のグレイ・バリューは、サーチピクセルｊのグレイ・バリューの重み係数Ｗ_ｊで重み付けされた平均値から見出される。この平均値は、最も外れた値を除いた平均値あるいは他のタイプの平均値のいずれでもよい。

別の構成として、ピクセルｉの真のグレイ・バリューの見積もり方法は、以下のものとは異なる。サーチピクセルのグレイ・バリューの計算方法、サーチピクセルの特徴ベクトルＦ_２の計算方法、ピクセルｉの特徴ベクトルＦ_１との計算方法。例えばサーチピクセルは、クラスに分割され、各クラスが、その特徴ベクトルＦ_２のクラスタリングに基づいて様々な組織の種類を表し、ピクセルｉと同一クラスのサーチピクセルのみを用いて、ピクセルｉの真のグレイ・バリューを見積もる、あるいは、そのサーチピクセルが、ピクセルｉの見積もられた真のグレイ・バリューに対し大きな影響を持つ。別の方法として、ある基準で特徴ベクトルＦ_１に最も近い特徴ベクトルＦ_２を有する上位の少数のサーチピクセルｊを用いて、ピクセルｉの真のグレイ・バリュー値を見積もる。選択的事項として、寄与するサーチピクセルのグレイ・バリューの平均値を用いる代わりに、ピクセルｉの見積もられた真のグレイ・バリューが、サーチピクセルのグレイ・バリューに基づいて、ルックアップテーブルから見出される。

修正されたグレイ・バリューは、元のピクセルｉのグレイ・バリューとサーチピクセルの重み付けされた平均値の線形の組み合わせでもよい。選択的事項として、元のピクセルのグレイ・バリューは、明白には考慮に入れられていないが、選択的事項として、ピクセルｉそのものを別のサーチピクセルとして取り扱い、重み付けされた平均に入れてもよい。この場合、ピクセルｉの重み係数Ｗｊは、Ｆ２の特徴値がＦ１の対応する特徴値と同様に規定された場合には、１であるが、Ｆ２の特徴値が回転した或いは換算した近傍で別に定義される場合には、１未満である。
ここに記載したピクセルのグレイ・バリューは、必ずしも画像の元のグレイ・バリューではなく、変換された画像あるいはフィルタ処理された画像、例えばσが数ピクセルの幅以下のようなガウシアンフィルタ処理されたイメージのグレイ・バリューでよい。
ステップ３２８において、ピクセルのインデックスｉを１増加させ、ステップ３３０において、考慮すべきピクセルが残っているか否かを決定する。残っている場合には、次のピクセルｉがステップ３０６で考慮される。残ってない場合には、本発明の方法は、ステップ３３２で終わる。この場合ノイズを低減した画像は、ステップ３２６で見出された修正されたグレイ・バリュー値を用いて、出力される。

図６Ａは、ノイズを含む画像４００即ち頭をスライスしたＣＴ画像を示し、これにより図５の方法を示す。この画像は通常のものよりノイズを多く含む。その理由はｘ線の線量を減らしたからである。図６Ｂはノイズを低減した画像４０２を示す。これは図５の方法を用いて画像４００から得られたものであり、特徴と重み付けベクトルに付いては後述する。比較のために、図６Ｃは低ノイズの画像４０４である。これはこの種の画像に対し通常のｘ線量を用いて得たものである。ノイズを低減した画像４０２は、源画像４００よりもかなりノイズが少なく、より詳細に見ることができ、特に様々な組織の間でコントラストが低いような場所ではそうである。画像４０２の質は、低ノイズ画像４０４に近いが、画像４００以上である。

残差画像からの詳細の回復
Ｄは、源画像Ｉのノイズ除去したバージョンである。ノイズ除去した画像はより平滑にり、小さな構造物と器官は、平滑過ぎて見えなくなる傾向がある。ノイズ除去アルゴリズムを適用後、回復手順を用いて、ノイズ除去アルゴリズムにより一部除去された或いは平滑化された構造物を回復し強調する。
残差画像（残像）Ｒは、源画像Ｉとノイズ除去済み画像Ｄの差である、即ちＲ＝Ｉ−Ｄである。この残差画像は大部分がノイズであり、源画像から除去された付加的ノイズである。しかし空間の整合性を有する（小さな）構造体は、残差画像から非線形のエッジ補助フィルタを用いて、回復される。

本発明の一実施例において、ノイズ除去アルゴリズムは、類似の領域を平均化する。その結果、残差画像は、ファインスケールで検査した時は、大部分がノイズである。更に、残差画像は、隠れた情報を保有する。残差画像は、ファインスケールで検査した時には、ノイズに見えるが、ラージスケールで検査した時には、空間的整合性を示す。この様な隠れた構造物を暴露するために、非線形のエッジ保存フィルタ（non-linear edge-preserving fileter）を用いる。この非線形のエッジ保存フィルタは、３Ｄの残差画像内に２Ｄ表面の形態のエッジを探し出し、あるいは２Ｄの残差画像内に１Ｄ表面の形態のエッジを探し出し、これらのエッジに平行な残差画像を平滑化する。代表的な非線形のエッジ保存フィルタは、非線形違法性拡散フィルタ(nonlinear anisotropic diffusion filtes)、例えばBeltrami flow filters、coherence-enhancement diffusion(CED)filtersである。フィルタ処理された残差画像Ｓは、空間的に一貫性（整合性）のある構造体を表し、残留ノイズを平滑化してなくす。本発明の一実施例においては、エッジあるいは他の構造体はＳの中で強化（強調）されるが、ノイズは平滑化して取り除くことはできない。別の構成として、残差画像は、フィルタ処理せず、フィルタ処理していない残差画像Ｒをフィルタ処理された残差画像Ｓの代わりに用いることもできる。どちらの方法が使用されようとも、回復されたノイズ除去済み画像が得られる。回復されたノイズ除去済み画像は、残差画像Ｒ又はフィルタ処理された残差画像Ｓのグレイ・バリューに対し、第１場所では第２場所より感受性がある。この第１場所では、残差画像は大きな構造体を有し、第２場所では、残差画像は小さな構造体を有する。

フィルタ処理された残差画像Ｓをノイズ除去済み画像Ｄに戻し、平滑化されていた詳細構造物が戻された回復されたノイズ除去済み画像Ｄ’を得る。選択的事項として、画像は、Ｄにα・Ｓを加えることにより、過度に強調される。αは１を越える。選択的事項として、適合パラメータαを用いるが、パラメータαは、Ｓのローカルな構造体程度（structureness）に依存する。ある領域が一貫性のある構造体あるいは首尾一貫性（coherence）を有するかを見積もるために、あらゆる措置を講じることができる。例えば、Ｓの構造体テンソルの固有値又はＳのHessian matrixの固有値を用いて、各場所におけるパラメータαを決定することもできる。選択的事項として、残差画像Ｒ内の情報は、所定の場所でのみノイズ除去済み画像Ｄに回復される（戻される）。所定の場所では、構造体は残差画像内に見出され、その結果ノイズはノイズ除去済み画像には殆ど戻らない。
現場では、放射線医師は、少量のノイズが画像内に残るのを好み、そしてそれはより実際に近く見える。選択的事項として、これを達成するために、パラメータαは、残差画像内にノイズのみが表れるような領域でも、０以上に維持される。

Ｄ’＝Ｄ＋α・Ｓを用いると、誤った構造を、比較的平滑な領域で画像に付加することになる。選択的事項として、Ｓの密度の範囲が限定される。別の方法は、αをエンベロープと見なし、Ｓが各ロケーションにおける画像の密度を変化させる量を制限する。この方法においては、Ｓの切頭されたあるいは圧縮されたバージョンが、Ｄに付加される。
Ｄ’＝Ｄ＋α・truncate(S/α) 又は
Ｄ’＝Ｄ＋truncate(α).S
ここで、切頭関数（truncate function）は、例えばハイパーボリック・タンジェント（hyperbolic tangent（tanh）)関数、範囲圧縮関数ｙ（ｘ）、あるいは「しきい値迄は直線状に大きくなりでしきい値以上は一定であるような関数」である。αの局部値が、密度の許容範囲を決定する。この密度の許容範囲は、Ｓが、所定のボクセルでノイズ除去済み画像Ｄに付加される或いはそこから除かれる密度の範囲のことである。

選択的事項としてsharpening filterを、Ｄ又はＤ’に適用する。このフィルタは、適用されると、類似のエンベロープにより制限される。これは例えば画像の局所標準偏差（local standard deviation）に基づいて行われる。選択的事項としてsharpening filterは、標準偏差が小さい領域でのみ、或いはでより強く適用される。エンベロープがない場合にはsharpening filterは、視覚的なアーチファクトを有し、特にシャープなエッジ（例、頭蓋骨のエッジ）近傍で視覚的アーチファクトを有する。エンベロープを使用すると、更なるシャープ差は相対的なしきい値を越えることがなく、視覚的なアーチファクトが減るかあるいは完全に回避できる。本発明の一実施例においては、フィルタは、標準偏差が中間の領域で、強く適用され、標準偏差が高い領域では適用されない。標準偏差が高い領域では、フィルタは視覚的アーチファクトを誘導し、標準偏差が低い領域では、フィルタはノイズを増幅するだけである。

ＣＥＤスキームの新たな変形例は、ノイズ除去済み画像Ｄから構造情報を取りだし、同時にフィルタを残差画像Ｒに適用する。回復手順は３つのステップに分けられる。第１ステップは、画像の構造テンソル（structure tesnor）を計算し、第２ステップは、各ピクセルに対する拡散テンソル(diffusion tensor)を、構造テンソルの固有値分解（eigendecomposition）に基づいて計算し、第３ステップは、拡散ステップを行う。選択的事項として、繰り返しプロセスを、（１−λ）・Ｄ＋λ・Ｓから計算された所定の繰り返しに対し構造テンソルで、用いる。ここでＳは、前の繰り返しに続くフィルタ処理された残差画像であり、係数λは、繰り返しが継続する連れて１に徐々に近づく。

ノイズをノイズ除去済み画像に追加すること
ノイズ除去済み画像は、放射線医師には不自然に見えることがある。放射線医師はある量のノイズを持ったＣＴ画像を日頃検査しており、更に良好な画像状態の元でさえそのようにしている。より自然に見える画像を生成するために、Ｄ’とＩの重み付けされた平均画像を用いる。例えば、画像Ｅ＝βＤ’＋（１−β）Ｉがある。ここでβは空間的に変動する重み付けパラメータであり、βの選択は、得られた画像Ｅがノイズを含み、しかもこのノイズが、入力イメージＩの中のノイズと同様に分布し、その程度は小さいものとなるよう、行われる。回復手順が用いられない場合には、Ｄは、Ｄ’の代わりに用いてもよい。

一実施例においては、重み付けパラメータβは次のように計算される。まず局所標準偏差のマップが、Ｄ’とＩに対し計算される。標準偏差のマップは、ＳＴＤ_Ｄ’とＳＴＤ_Ｉとして示される。その後βに対する予備的な値を、次式に基づいて計算する。

ここでσ_Ｎは、、追加されたノイズの量を制御し、σ_Ｎ＝４であり、εは小さな正の定数で、例えばε＝１０^−３である。次にβは、範囲圧縮される例えばこれは上記で決定された範囲圧縮関数ｙ（ｘ）を用いて行われる。例えばａ＝０でｂ＝１である。最後に、βは、ガウスフィルタでフィルタ処理することにより平滑化される。

解像度低減
放射線医師は、日頃ＣＴ装置が得る源画像の解像度よりも低い解像度のＣＴ画像を検査している。一般的に検査されている画像のスライス厚さは、ＣＴ装置により再構成される生の画像のスライス厚さの４倍である。
選択的事項として、上記のノイズ除去アルゴリズムとそれに関連する手順が適用される源画像Ｉは、ＣＴ装置から得られる高解像度の残差画像である。放射線医師が検査するＣＴ画像用に使用される低解像度の画像を得るために、最終解像度低減手順が行われる。従来はこの様な解像度低減手順は、通常最初にローパスフィルタ（例、anti-aliasing filterを用いて）で画像をフィルタリングし、その後所望の解像度になるよう、サブサンプリングを行う。
本発明の一実施例においては、非線形解像度低減手順（non-linear resolution reduction procedure）が用いられるが、これは、従来の解像度低減手順よりも、より詳細な情報を保持する。例えばこの方法の一般論は非特許文献１１に記載されている。

上記の方法は、異なるカメラ焦点設定で撮られた２つの画像を合成して、１つの画像を得る方法である。この１つの画像は、２つの画像からの重要なデータ特徴を含む。同様な方法を用いて、複数のＣＴスライスを合成することもできる。この方法は、他の画像の寸法に沿って別々の方法で適用され、検査された画像の解像度を得る。この様な非線形のダウンサンプリング手順は、上記の文献に記載されているが、スライス方向に配列した複数のスライスを有する３つのノイズ除去済み画像に適用されると、この手順は、スライス方向には解像度を低減した画像が得られるが、しかしスライスに沿った方向の解像度は落ちていない。例えば、スライス方向に直行する面のスライスの面では、解像度は落ちていない。

最初に、帯域通過の分解（band-pass decomposition）を、高解像度の生のＣＴ画像のスライスに対し実行する。例えばこのスライスをラプラシアンピラミッド（Laplacian pyramid）に変換する。ラプラシアンピラミッドは、上記の文献に記載されている。その後、ｎ_Ｓの連続するスライスの各組（例、４枚のスライスの組）に対し、１枚の厚いスライスに対応する１個のピラミッド構造を、以下の手順で形成する。
ピラミッドのあるレベルに対しては、非線形の手順を用いて、ｎ_Ｓ枚のスライスを組み合わせて、厚いスライスのピラミッドのレベルを形成する。しかしこれはピラミッドの全てのレベルに対して行う必要はない。例えば、ピラミッドのトップ・レベル（画像の低周波数成分を含む）に対しては、厚いスライスが、ｎ_Ｓ個のラプラシアンピラミッドの上部レベルの画像の全てを平均化すること（線形手順）により、形成される。ピラミッドの他のレベルに対しては、非線形の手順を用いて、より厚いスライスを形成する。これは、ある場合には、大きい重み係数を薄いスライス（局部的に大きな振幅を有する）に加えることになる。例えば、厚いスライスのピラミッドの他のレベルは、各場所においてｎ_Ｓ個のピラミッドの内の１つのピラミッドの各値を取り込むことにより、形成される。この選択された値は、全てのｎ_Ｓ個のピラミッドの中で最も高い絶対値を有する値である。但しこれは、この値があるしきい値（例えば１０ＨＵ）より高くない場合である。この場合、全てのｎ_Ｓ個のピラミッドの平均値が取られる。
最後に、厚いＣＴスライスが、形成されたピラミッド構造を再構成することにより形成される。全ての厚いスライスが、非線形の手順で、低解像度の画像を形成する。この非線形の手順は、各スライスに沿って即ちスライス方向に直行する方向に沿って、高精細の詳細構造を保存する傾向がある。
本発明の一実施例においては、ラプラシアンピラミッドを帯域通過分解（band-pass decomposition）として用いることに加えてあるいはそれに代えて、あらゆるタイプの帯域通過分解を用いることができる。例えばウェーブレット変換を用いることができるが、これは非特許文献２に開示されている。更に非特許文献２に記載されている画像を合成する方法も用いることができる。
選択的事項として、この非線形スキームにより形成された低解像度の画像は、線形にサブサンプリングされた低解像度画像と、重み付け平均で、組み合わされる。

ＮＬＭノイズ除去アルゴリズムで使用される特徴のタイプ
幾つかのタイプの特徴が特徴ベクトルＦ_１，Ｆ_２で使用される。
本発明の一実施例においては、特徴値を計算することは、近傍のピクセルのグレイ・バリューの分布の特性を見出すことを含む。特徴値は、グレイ・バリューの分布のモーメント、あるいはそのモーメントの関数である。
分布の第１モーメントは平均値であり、第２モーメントは標準偏差であり、第３モーメントはスキュ（skew）である。分布のｋ番目のモーメント（ｋ＞１）は、次式で定義される。

ここで、Ｉ_ｎは、近傍のｎ番目のピクセルのグレイ・バリューである。加算は近傍内のＮ個のピクセルに対し行われ、Ｍ_１は第１モーメント即ちグレイ・バリューの平均値である。更に特徴値は、分布の順序統計量、分布の所定の割合に対応するグレイ・バリューである或いはそれに依存する、特徴値は、グレイ・バリューの中央値であり、これは５０％のパーセンタイルのグレイ・バリューである。別の構成として、異なるパーセンタイルのグレイ・バリューも使用することができる。例えば２５％，３７．５％，６２．５％，７５％のパーセンタイルも使用できる。選択的事項として、中間値のパーセンタイルを用いることもできる。例えば２５％と７５％の間の中間値も用いることもできる。これは、特徴値が、全体として近傍の特性であるが、近傍内の少ない外側のピクセルの特徴ではないという利点がある。選択的事項として、サーチピクセルが、他の画像のサーチピクセルを含むようなディクショナリから選択され、別の方法でグレイ・バリューが正規化された場合には、２個の画像のグレイ・バリューが、正規化され、それらが、順序統計量（order statistics）に基づいて特徴内で比較される。

近傍内のピクセルのグレイ・バリューの分布特性にのみ基づく特徴は、特に近傍が完全な四角形あるいは他の等方性形状の場合には、特徴値が画像内の構造の方向には無関係であるという利点がある。図４Ｂの画像２０４に対しこの様な特徴を用いることは、ピクセル１１０に近い特徴を有するピクセル２０６の様なピクセルの組を生成する。その理由は、特徴値は、暗い領域１０４からのピクセルの距離に大部分が依存し、暗い領域と明るい領域の間の境界の局部的な方向には依存しないからである。他方で、画像のある部分が、特定の方向例えば身体の組織のセグメンテーションマップから特定の方向に方向付けられたエッジあるいはテクスチャを有する場合には、構造体の方向性に敏感な特徴を用いるのが好ましい。

選択的事項として、特徴値は、近傍のピクセルの生のグレイ・バリューの分布から見出されるのではなく、画像が平滑化された後あるいは他の方法で処理された後のグレイ・バリューの分布から見出される。この様な特徴値を評価する前に画像を平滑化することは、次のような利点がある。即ち特徴値は近傍の画像の構造体の特性に多く依存し、近傍のノイズには敏感でないという利点である。選択的事項として、この様な平滑化処理あるいは画像処理は、特徴のあらゆる種類に対し特徴値を評価する前に、行い、グレイ・バリューの分布に依存する特徴に対しては行わない点である。平滑化は、例えばガウシアンフィルタ（Gaussian filter）、双務フィルタ（bilateral filter）、トータル変動フィルタ（total variation filter）を用いて、行われる。選択的事項として、平滑化は、特徴を得る為に使用される近傍の最大寸法のスケールで、あるいは近傍の最小寸法のスケールで、構造体の大部分を平滑化して除去しないような方法で行われる。例えばガウシアンフィルタが、幅パラメータσで用いられる場合には、パラメータσは、近傍の最大寸法更には近傍の最小寸法よりも小さいあるいは余り大きくはない値である。別の構成として、平滑化は、近傍の空間構造体の全てを効果的に平滑して除去するような方法で行われる。特徴値は、近傍内にない構造体の測定値であるが、近傍の周りのより大きなスケールの構造体の測定値、あるいは近傍の周りの平均傾斜の測定値である。

選択的事項として、画像はグレイ・バリューの分布を見出す前に、様々な方法で処理される。例えば、導関数演算子（derivative oprator）を画像に適用して、各ピクセルのグレイ・バリューを、特定の方向の画像の導関数（derivative）に比例する値で、あるいは画像の傾斜の大きさに比例する値で、置換する。これが行われると、近傍のピクセルの値の分布の平均は、近傍の平均傾斜の測定値となる。選択的事項として、画像が、傾斜を見出す前に、平滑化されるあるいは十分に平滑される。その結果、近傍のピクセルの大部分が同一の傾斜となり、特徴値をノイズに対し敏感とならない。

本発明の一実施例によれば、特徴値を計算することは、少なくとも近傍内の画像に変換又はフィルタ（transform or filter）を適用することを含む。この変換又はフィルタは、近傍の最大寸法とピクセルとの間のサイズ的に中間にあるスケールの範囲の構造物を選択する。別の構成として、変換又はフィルタは、他の方向よりもある方向に向いた構造物を選択する。この様に規定された特徴は、ある範囲内の寸法又は方向を有する期待される構造物（例、血管）をピックアップするのに役立つ。しかしノイズに起因する密度のより細かいスケールの変動は無視する。

本発明のこれらの特徴は、様々なフィルタと変換、あるいは線形と非線形の手法を採用する利点を有する。これらの利点は、コンピュータによる手書き文字の認識のようなアプリケーションで使用されるか、あるいは画像内の対象物を自動的に明確にし、テキストの記述によらずにサーチできる利点である。
この様な特徴は、近傍の画像のウェブレットフィルタ（waelet filter）への応答に依存する。例えばMeyer or Gabor filters, Laplacian and Gaussian Pyramidsあるいは様々なガウス分布フィルタあるいは従来公知の他の線形フィルタである。この様なフィルタは、特定の方向あるいは特定のスケールを有する近傍内の構造に対して、最も敏感である。
選択的事項として、フィルタは、近傍にのみ適用することもできる。別の構成として、フィルタは、近傍よりも広い領域にも適用できる。更には画像全体に対しても適用できる。特徴は、近傍のフィルタの応答に依存する。例えば、画像がフィルタ処理された後、近傍内のピクセルのグレイ・バリューに依存する。これらの選択的事項は、フィルタ又は変換を画像ピクセルに適用することに関連する他の種類にも応用可能である。

選択的事項として、特徴値は、２つの異なるサイズパラメータσ_１，σ_２において、画像のガウシアンフィルタあるいは他の平滑化フィルタに対する応答の差に依存する。こ２つののフィルタの差は、σ_１，σ_２の間のスケールの中間にある構造体を選択する傾向にあるが、これはフィルタが等方性の場合には、構造体の方向性には依存しない。かくして規定された特徴は、画像が様々な方向を向いた類似の構造を有する場合に、有益である。他の画像（特定の方向を向いたあるいは狭い範囲の方向を向いた構造体を有するとして既知である）に対しては、他の特徴の測定値は、構造体の方向に依存し、使用するのに利点がある。

本発明の一実施例においては、特徴は、非線形変換（例、形態のマルチスケール変換あるいは形態演算子）に対する近傍内の画像の応答に依存する。例えば特徴値は、非線形のマルチスケール変換を画像に特定のスケールパラメータを用いて適用した後、検査中のピクセルのグレイ・バリューあるいは近傍内の特定のピクセルに依存する。
選択的事項として、特徴値は、複数のスケールパラメータを用いたピクセルのグレイ・バリューに依存する。例えば２つのスケールパラメータ用のグレイ・バリューの差に依存する。形態的なマルチスケール変換の一例は、形態ウェーブレットあるいは形態ピラミッドを含み、これは非特許文献１１に記載されている。

別の構成として、特徴は、形態演算子（morphological operator）を適用した後、ピクセルのグレイ・バリューに依存する。形態演算子を画像に適用して、ある構造を強化あるいは抽出する。形態演算子の一例はtop hat transformである。即ち入力画像とその構造体要素による形態開口の間の差である。この様な演算は、暗い背景の中に明るい詳細画面を映し出し、構造を構成する要素のサイズが、検出された特徴のサイズを制御する。類似の演算子が、白色の背景に対し暗い構造体を抽出するようにすることもできる。

形状マッチあるいは画像形態の文献には様々な技術を含み、これらの技術を用いて、画像内の形状を特徴付け、近傍の特徴を規定するのに用いられる。これは形態変換あるいは演算子を適用した後に行われる。この一例はRubuner et al; the により導入されたEarth Mover's Distanceであり、Kantorovich-Wasserstein metricが画像の形態形成に用いられるがこれはHaker et al; the shape に記載されている。形状sinnatures defined by Osada et al;に記載されており、形状マッチングの方法はBelongie et al; all referenced aboveに記載されている。

本発明の一実施例においては、検査中のピクセルｉの特徴値とサーチピクセルｊの特徴値は、一方から他方に移行する位置変換（geometric transformation）で変更された近傍で計算される。例えば２個の近傍は、相対的に異なる方向、スケールあるいはその両方を有する。本発明の変形例として、近傍の１つは、互いにミラー反射される（勝手違いに配置される）。例えばピクセルｉの特徴値を見出すアルゴリズムは、ピクセルｉから＋ｘ方向にある距離離れた近傍のピクセルのグレイ・バリューを用い、代わりに、ピクセルｊの特徴値は、ピクセルｊから＋ｙ方向（９０度回転）の方向に同一距離離れたピクセル、あるいは−ｘ方向（反射）に同一方向離れたピクセル、＋ｘ方向に２倍離れた（スケール変換）ピクセル、＋ｙ方向に２倍離れた（回転＋スケール変換）ピクセルを用いる。
選択的事項として、回転したあるいはスケール変換された近傍のグレイ・バリューは、特徴値を計算する前に、補間され、特に回転角が９０度の整数倍でない場合（ピクセルのデカルト・グリッドの場合）、あるいはスケールファクタが整数でない場合は、そうである。ピクセルは三角形グリッドあるいは六角形グリッドで、あるいはより複雑なパターンに配置することもできる。

この様に規定された特徴用いることは、回転、スケール変換、反射がなく同一の特徴が規定された場合に、特に有効であり、様々な回転角や換算係数に対し有効であり、その結果得られた特徴は同一の重み付けを与える。これは画像内の構造の方向とスケールとは独立した距離の測定値であり、これはある範囲の方向とスケールに対し行われる。この方法は画像が異なる方向とスケールを有する類似の構造体を含む場合、あるいは互いにミラー状態のある場合特に有効である。

以上の説明は、本発明の一実施例に関するもので、この技術分野の当業者であれば、本発明の種々の変形例を考え得るが、それらはいずれも本発明の技術的範囲に包含される。特許請求の範囲の構成要素の後に記載した括弧内の番号は、図面の部品番号に対応し、発明の容易なる理解の為に付したものであり、発明を限定的に解釈するために用いてはならない。また、同一番号でも明細書と特許請求の範囲の部品名は必ずしも同一ではない。これは上記した理由による。用語「又は」に関して、例えば「Ａ又はＢ」は、「Ａのみ」、「Ｂのみ」ならず、「ＡとＢの両方」を選択することも含む。特に記載のない限り、装置又は手段の数は、単数か複数かを問わない。

図１
１００２：画像を獲得する。
１００４：画像を前処理する。
１００６：ローカルノイズを見積もる。
１００８：ノイズ除去アルゴリズムを適用する。
１０１０：構造を回復する。
１０１２：ノイズを付加する。
１０１４：解像度を低減する。


図５
３０２：画像を獲得する。
３０４：ｉ＝１
３０６：画像内のピクセルを見る。
３０８：ピクセルｉの特徴ベクトルＦ_１を見出す。
３１０：ｊ＝１
３１２：サーチウィンドウあるいは方向内サーチピクセルｊを見る。
３１４：サーチピクセルｊの特徴ベクトルＦ_２を見出す。
３１６：距離測定値ｄＦ_１Ｆ_２を見出す。
３１８：サーチピクセルｊの重み係数Ｗ_ｊを見出し記録する。
３２０：ｊ＝ｊ＋１
３２２：サーチピクセルで実行したか？
３２４：サーチピクセルの重み付けされた平均グレイ・バリューを見出す。
３２６：ピクセルｉの修正されたグレイ・バリューを見出す。
３２８：ｉ＝ｉ＋１
３３０：画像内のピクセルで実行したか？
３３２：終了

Claims (8)

1. 所定の解像度でノイズ除去された医療画像を提供する方法において、
   前記医療画像は、元のスライス厚さより薄い複数のスライスからなる３次元医療画像であり、
   （ａ）前記所定の解像度よりも高い解像度で得られた源医療画像をノイズ除去するステップと、
   （ｂ）前記ステップ（ａ）の後の医療画像の解像度を、前記所定の解像度に戻すステップと、
   を有し、
   前記（ｂ）ステップは、前記複数の薄いスライスを非線形解像度低減手順を用いて合成し、厚いスライスの医療画像を得るステップを含み、
   前記非線形解像度低減手順は、
   （ｘ）前記薄いスライスを帯域通過成分に分解するステップと、
   （ｙ）前記帯域通過成分に対し、より大きな重み係数をより大きな振幅を有する薄いスライス成分に与えるようにして、前記薄いスライス成分を厚いスライス成分になるよう組み合わせるステップと、
   （ｚ）前記薄いスライスの帯域通過成分を再結合して、所定の解像度の厚いスライスの医療画像を生成するステップと、
   を有する
   ことを特徴とするノイズ除去された所定の解像度の医療画像を提供する方法。
2. 前記ステップ（ｙ）は、最低周波数成分を含むスライス成分以外の帯域通過成分に対して行われ、
   前記厚いスライスに組み合わされる全ての薄いスライスの中で、前記薄いスライス成分の最大絶対値はしきい値以下である
   ことを特徴とする請求項１記載の方法。
3. 前記源医療画像は、スライス方向に配置された薄いスライスを有する３次元医療画像である
   ことを特徴とする請求項１−２のいずれかに記載の方法。
4. 前記非線形解像度低減手順は、帯域通過分解手順を含む
   ことを特徴とする請求項１記載の方法。
5. 前記帯域通過分解手順は、ラプラシアン・ピラミッド手順を含む
   ことを特徴とする請求項４記載の方法。
6. 前記帯域通過分解手順は、ウェーブレット変換手順を含む
   ことを特徴とする請求項４記載の方法。
7. 前記（ａ）ステップは、
   （ａ１）医療画像内で検査中の各ピクセルに対し、サーチピクセルの組を選択するステップと、
   （ａ２）前記の各サーチピクセルの近傍の特徴の値と、前記検査中のピクセルの近傍の対応する特徴の値を計算するステップと、
   （ａ３）前記の各サーチピクセルの同一又は異なる近傍の他の特徴の値と、前記検査中のピクセルの同一又は異なる近傍の対応する値とを計算するステップと、
   （ａ４）ノイズを低減した状態の検査中のピクセルのグレイ・バリューを、そのサーチピクセルのグレイ・バリューに基づいて、計算するステップと、
   を有する
   ことを特徴とする請求項１−６記載の方法。
8. 医療画像装置とコントローラとを有する所定の解像度のノイズ除去された医療画像を提供するシステムにおいて、
   前記医療画像装置は、所定の解像度よりも高い解像度で３Ｄ医療画像に対するデータを生成し、
   前記コントローラは、高解像度のデータから医療画像のスライスを再構成し、再構成された医療画像をノイズ除去し、非線形解像度低減手順を用いて、前記ノイズ除去した医療画像の解像度を所定の解像度まで低減し厚いスライスを形成するよう薄いスライスを合成し、
   前記非線形解像度低減手順は、
   （ｘ）前記薄いスライスを帯域通過成分に分解するステップと、
   （ｙ）前記帯域通過成分の一部に対し、より大きな重み係数をより大きな振幅を有する薄いスライス成分に与えるようにして、前記薄いスライス成分を厚いスライス成分になるよう組み合わせるステップと、
   （ｚ）前記薄いスライスの帯域通過成分を再結合して、所定の解像度の厚いスライスの医療画像を生成するステップと、
   を有する
   ことを特徴とする所定の解像度のノイズ除去された所定の解像度の医療画像を提供するシステム。

Images