JP5946800B2 - 磁気共鳴イメージング装置、画像処理装置、画像処理方法、及び画像処理プログラム - Google Patents

Description

本発明は、磁気共鳴イメージング（Ｍａｇｎｅｔｉｃ Ｒｅｓｏｎａｎｃｅ Ｉｍａｇｉｎｇ、以下、ＭＲＩ）装置で取得した画像を処理する技術に関する。特に、拡散尖度イメージング技術に関する。

ＭＲＩ装置は、主にプロトンの核磁気共鳴現象を利用した医用画像取得装置である。ＭＲＩ装置は、非侵襲に任意の断面の撮像が可能であり、形態情報の他、血流や代謝機能などの生体機能に関する情報を取得可能である。ＭＲＩで撮影される重要な画像の一つに、拡散強調画像（ＤＷＩ：Ｄｉｆｆｕｓｉｏｎ Ｗｅｉｇｈｔｅｄ Ｉｍａｇｅ）がある。拡散強調画像は、生体組織に含まれる水分子の自己拡散を強調した画像であり、撮像時は、拡散速度に応じた信号を取得するため、ランダム運動する核スピンに対して位相分散による信号強度低下を誘起するＭＰＧ（Ｍｏｔｉｏｎ Ｐｒｏｂｉｎｇ Ｇｒａｄｉｅｎｔ）パルスを、リフォーカス用高周波磁場パルスの前後に印加し、エコーを取得する。

位相分散による信号強度低下は、ＭＰＧパルスの印加方向に拡散する核スピンによって引き起こされるため、ＭＰＧパルスの印加方向を制御することで任意方向の拡散情報を取得できる。また、拡散強調度は、ＭＰＧパルスの印加強度と印加時間に関するパラメータである拡散因子（ｂ値）にて調整でき、ｂ値が高いほど拡散強調度の高い画像を取得できる。

水分子の空間的な拡散分布を計測する手法として、ＤＴＩ（Ｄｉｆｆｕｓｉｏｎ Ｔｅｎｓｏｒ Ｉｍａｇｉｎｇ）がある。ＤＴＩでは、水分子の空間的な拡散分布が正規分布の３次元楕円拡散モデルに従うものと仮定し、その異方性比率ＦＡ（Ｆｒａｃｔｉｏｎａｌ Ａｎｉｓｏｔｒｏｐｙ）を算出することで、白質の神経走行路の構造を解析する。ＤＴＩのパルスシーケンスは、ＭＰＧパルスの印加方向を変えながら拡散強調画像ＤＷＩを取得するパルスシーケンスを繰り返すように構成される。

また、近年、細胞膜や細胞内小器官などによる拡散運動の制限度合いを強調する手法として、水分子の空間的な拡散分布を非正規分布の拡散モデルと仮定した拡散尖度イメージングＤＫＩ（Ｄｉｆｆｕｓｉｏｎ Ｋｕｒｔｏｓｉｓ Ｉｍａｇｉｎｇ）が提案されている。この手法は、正規分布の拡散モデルと仮定したＤＴＩと比し、組織変性や細胞増殖に伴う微細な構造変化を捉える手法として期待されている。ＤＫＩのパルスシーケンスは、ｂ値を変更しながらＤＴＩのパルスシーケンスを繰り返すように構成される。

一般に、ＤＫＩで得た画像の解析（以後、ＤＫＩ解析）では、ＭＰＧパルス印加方向が等しく、かつ、ｂ値が異なるシーケンスにより得た拡散強調画像において、画素ごとに非線形の最小二乗フィティング処理を実行し、ＭＰＧパルスの印加方向毎の拡散に関するパラメータとして拡散係数と尖度係数を推定する（例えば、非特許文献１）。そして、各係数の空間分布を描写するため、例えば、拡散テンソルおよび尖度テンソルの成分を計算し、第一主成分方向や主成分に垂直な方向の拡散係数と尖度係数等を算出する。

Ｆａｌａｎｇｏｌａ ＭＦ、Ｊｅｎｓｅｎ ＪＨ、Ｂａｂｂ ＪＳ、Ｈｕ Ｃ、Ｃａｓｔｅｌｌａｎｏｓ ＦＸ、Ｍａｒｔｉｎｏ ＡＤ、Ｆｅｒｒｉｓ ＳＨ、Ｈｅｌｐｅｒｎ ＪＡ著、"Ａｇｅ−ｒｅｌａｔｅｄ ｎｏｎ−Ｇａｕｓｓｉａｎ ｄｉｆｆｕｓｉｏｎ ｐａｔｔｅｒｎｓ ｉｎ ｔｈｅ ｐｒｅｆｒｏｎｔａｌ ｂｒａｉｎ" Ｊｏｕｒｎａｌ ｏｆ Ｍａｇｎｅｔｉｃ Ｒｅｓｏｎａｎｃｅ Ｉｍａｇｉｎｇ ２８，２００８，ｐ１３４５−１３５０

ＤＫＩ解析において高画質な画像を短時間で取得するには、拡散に関するパラメータ推定における非線形の最小二乗フィッティングの安定化および高速化が課題となる。非特許文献１に記載の手法では、計算を安定化させるため、前処理として、拡散強調画像全体に平滑化フィルタを適用する。従って、拡散に関するパラメータ推定の元となる拡散強調画像にボケが発生し、結果として得られるパラメータ画像の品質にも影響を与える。また、一般的な計算安定化手法として、制限付きの非線形最小二乗フィッティングがあるが、制限付きの場合、初期値を複数回変更し、計算を繰り返す必要があり、処理時間が長くなる。

本発明は、上記事情に鑑みてなされたもので、ＤＫＩ解析において、高速に高品質の画像を得る技術を提供することを目的とする。

本発明は、ＤＫＩ解析において、ＭＰＧパルスの印加方向毎の拡散に関するパラメータの推定時に、最小二乗フィッティングと制約処理とを分離し、最小二乗フィッティングにおいて制約条件の範囲外となった画素の画素値についてのみ補正処理を行う。そして、当該画素について、補正後の画素値から拡散に関するパラメータを再推定し、得られた拡散に関するパラメータを用いて、パラメータ画像を生成する。

本発明によれば、ＤＫＩ解析において、高速に高品質の画像を得ることができる。

第一の実施形態のＭＲＩ装置のブロック図である。 第一の実施形態の計算機の機能ブロック図である。 拡散強調画像取得シーケンスのパルスシーケンスを説明するための説明図である。 第一の実施形態の画像処理のフローチャートである。 第二の実施形態のパラメータ推定処理のフローチャートである。 第三の実施形態の画素値補正処理のフローチャートである。 第四の実施形態の計算機の機能ブロック図である。 第四の実施形態の画像処理のフローチャートである。 第四の実施形態のパラメータ補正処理のフローチャートである。 第五の実施形態の計算機の機能ブロック図である。 第五の実施形態の画像処理のフローチャートである。

＜＜第一の実施形態＞＞
以下、本発明を適用する第一の実施形態について説明する。本発明の実施形態を説明するための全図において、特に断らない限り、同一機能を有するものは同一符号を付し、その繰り返しの説明は省略する。

まず、本実施形態の磁気共鳴イメージング（ＭＲＩ）装置について説明する。本実施形態のＭＲＩ装置１００は、静磁場に置かれた被検体にスライス傾斜磁場を印加すると同時に特定の周波数をもつ高周波磁場を印加し、撮影したい断面内の核磁化を励起させる。励起された核磁化に位相エンコード傾斜磁場およびリードアウト傾斜磁場を印加することにより、平面位置情報を与え、核磁化が発生する核磁気共鳴信号（エコー）を計測する。ｋ空間と呼ばれる計測空間が充填されるまで核磁気共鳴信号の計測を繰り返し行い、ｋ空間に充填された信号を、逆フーリエ変換により画像化（撮影）する。

図１は、これを実現する、本実施形態のＭＲＩ装置１００の典型的な構成を示すブロック図である。本実施形態のＭＲＩ装置１００は、静磁場を発生するマグネット１０１と、傾斜磁場を発生する傾斜磁場コイル１０２と、被検体（生体）１０３に高周波磁場パルス（以下、ＲＦパルス）を照射するＲＦコイル１０７と、被検体１０３から発生するエコー信号を検出するＲＦプローブ１０８と、マグネット１０１の発生する静磁場空間内で被検体（例えば、生体）１０３を載置する寝台（テーブル）１１５と、を備える。

さらに、本実施形態のＭＲＩ装置１００は、傾斜磁場コイル１０２を駆動する傾斜磁場電源１０５と、ＲＦコイル１０７を駆動する高周波磁場発生器１０６と、ＲＦプローブ１０８で検出したエコー信号を受信する受信器１０９と、傾斜磁場電源１０５と高周波磁場発生器１０６とに命令を送り、それぞれ傾斜磁場および高周波磁場を発生させるとともに、検波の基準とする核磁気共鳴周波数を受信器１０９にセットするシーケンサ１０４と、検波された信号に対して信号処理（演算処理）を施す計算機１１０と、計算機１１０での処理結果を表示する表示装置１１１と、同処理結果を保持する記憶装置１１２と、操作者からの指示を受け付ける入力装置１１６と、を備える。

以上の構成を有するＭＲＩ装置１００では、シーケンサ１０４の制御により、ＲＦパルスがＲＦコイル１０７を通じて被検体１０３に印加されるとともに、スライス選択や位相エンコードなどの位置情報をエコー信号に与えるための傾斜磁場パルスが傾斜磁場コイル１０２によって印加される。また、被検体１０３から発生した信号はＲＦプローブ１０８によって受波され、検波された信号は計算機１１０に送られ、ここで画像再構成などの信号処理が行われる。なお、記憶装置１１２には、信号処理の結果だけでなく、必要に応じて、検波された信号自体、撮影条件等を記憶させてもよい。

また、ＭＲＩ装置１００は、静磁場均一度を調節する必要があるときには、シムコイル１１３と、シムコイル１１３を駆動するシム電源１１４と、をさらに備えてもよい。シムコイル１１３は、複数のチャネルからなり、シム電源１１４から供給される電流によりにより静磁場不均一を補正する付加的な磁場を発生する。静磁場均一度調整時にシムコイル１１３を構成する各チャネルに流す電流は、シーケンサ１０４により制御される。

また、計算機１１０は、検波された信号に対する信号処理（演算処理）だけでなく、ＭＲＩ装置１００全体の動作の制御も行う。例えば、シーケンサ１０４に指示を与え、上述の各部の制御を実現し、計測、撮像を実行させる。

なお、シーケンサ１０４への指示は、予めプログラムされたタイミング、強度で各部が動作するよう、行われる。上記プログラムのうち、特に、高周波磁場、傾斜磁場、信号受信のタイミングや強度を記述したものはパルスシーケンスと呼ばれる。計測は、パルスシーケンスとこれを制御するために必要な撮影パラメータとに従って行われる。パルスシーケンスは、予め作成され、記憶装置１１２に保持され、撮影パラメータは、操作者からユーザインタフェースを介して入力される。

このため、本実施形態の計算機１１０は、図２に示すように、所定のパルスシーケンスに従って、計測を行い、画像を再構成する撮像部２１０と、再構成された画像に対し、演算処理を施し、パラメータ画像を得る画像処理部２２０と、を備える。なお、本明細書において、パラメータ画像は、画素毎に計算された拡散係数および尖度係数の平均値や拡散テンソルおよび尖度テンソル等を用いて拡散係数および尖度係数の空間分布を描写した画像とする。

このパルスシーケンスは、目的に応じて種々のものが知られている。本実施形態では、このパルスシーケンスとして、傾斜磁場として拡散に伴う信号変化を核磁気共鳴信号に付加するＭＰＧ（ｍｏｔｉｏｎ ｐｒｏｂｉｎｇ ｇｒａｄｉｅｎｔ）パルスを印加する拡散強調撮像シーケンスを用いる。具体的には、ＤＫＩのパルスシーケンス（ＤＫＩパルスシーケンス）を用いる。ＤＫＩパルスシーケンスは、拡散強調画像を得るＤＷＩのパルスシーケンスを、ＭＰＧパルスの印加方向とｂ値とを変化させて繰り返し実行するものである。なお、拡散強調画像は、生体組織に含まれる水分子の自己拡散を強調した画像である。急性期脳梗塞の発症直後の病変描出が可能であり、Ｔ１強調画像やＴ２強調画像とは異なるコントラスを示すことが知られている。

図３に、ＤＷＩのパルスシーケンスの一例として、ｓｓ−ＤＷＥＰＩ（ｓｉｎｇｌｅ ｓｈｏｔ Ｄｉｆｆｕｓｉｏｎ ｗｅｉｇｈｔｅｄ Ｅｃｈｏ Ｐｌａｎａｒ Ｉｍａｇｉｎｇ）のパルスシーケンス３００を示す。

ｓｓ−ＤＷＥＰＩのパルスシーケンス３００の動作は以下の通りである。スライス方向（ｚ）の傾斜磁場パルス３０１の印加とともに、対象となる核磁化の共鳴周波数ｆ₀にてＲＦパルス３０２を印加し、対象物体内のあるスライスの核磁化に核磁気共鳴現象を誘起する。次に、核磁化のランダムな運動を計測するためのＭＰＧパルス３０３を印加する。このとき、ＭＰＧパルス３０３を印加する軸は、核磁化のランダム運動において計測したい方向に印加する。なお、本図においては、一例として、ｘ軸、ｙ軸、ｚ軸の３軸方向に印加する場合を示す。しかし、印加方向は、これらに限定されない。

次に、スライス傾斜磁場パルス３０４の印加とともに共鳴周波数ｆ₀の磁化リフォーカス用ＲＦパルス３０５を印加し、静磁場の不均一によって分散した磁化の位相を再収束させる。次に、ＭＰＧパルス３０３と同方向にＭＰＧパルス３０６を印加して、ＭＰＧパルス３０３によって分散した動きのない核磁化の位相をリフォーカスさせる。

次に、磁化の位相に位相エンコード方向（ｙ）の位置情報を付加するための位相エンコード傾斜磁場パルス３０７、ディフェーズ用リードアウト傾斜磁場パルス３０８を印加した後、リードアウト方向（ｘ）の位置情報を付加するためのリードアウト傾斜磁場パルス３０９を印加しながら核磁気共鳴信号（エコー）３１０を計測する。

その後は、位相エンコード傾斜磁場パルス３１１を印加し位相エンコード量を変化させ、リードアウト傾斜磁場パルス３１２を反転させ、再び核磁気共鳴信号（エコー）３１０を取得することを繰り返し、一枚の画像を得るのに必要なエコーを取得する。各エコー３１０は、ｋ空間に配置され、２次元逆フーリエ変換によって画像再構成される。

本実地形態の撮像部２１０は、ＤＫＩパルスシーケンスとして、上記ｓｓ−ＤＷＥＰＩパルスシーケンス３００を、前記ＭＰＧパルス３０３の印加方向およびｂ値を変化させて繰り返し実行して複数の拡散強調画像を、ＭＰＧパルス３０３の印加方向およびｂ値毎に得る。なお、ｂ値は、上述のように、ＭＰＧパルス３０３および３０６の印加強度および時間に関するパラメータである。従って、ｂ値を変更することにより、ＭＰＧパルス３０３および３０６の印加強度または時間が変化する。

例えば、ＭＰＧパルスの印加方向を異なるＭ方向に変化させ、ｂ値を異なるＮ個の間で変化させる場合、撮像部２１０は、ｓｓ−ＤＷＥＰＩパルスシーケンス３００を、Ｍ×Ｎ回実行する。本実施形態では、例えば、ＭＰＧパルス３０３の印加方向は、１５軸以上の異なる方向とする。また、ｂ値は、異なる３つ以上（ｂ＝０を含む）の値を用いる。

本実施形態の画像処理部２２０は、得られた複数の拡散強調画像から、画素毎に、拡散に関するパラメータ（拡散パラメータ）として、見かけの拡散係数ＡＤＣおよび見かけの尖度係数ＡＫＣを推定し、それを用いてパラメータ画像を生成する。本実施形態では、拡散パラメータの推定の際、制約条件付きの非線形最小二乗フィッティング処理を行う代わりに、制約条件無しの非線形最小二乗フィッティング処理を行い、得られた拡散パラメータが制約条件を満たさないものについてのみ、画素値を補正し、再度、拡散パラメータを推定する。また、本明細書では、見かけの拡散係数ＡＤＣおよび見かけの尖度係数ＡＫＣを、それぞれ、単に、拡散係数ＡＤＣおよび尖度係数ＡＫＣと呼ぶ。

本実施形態の画像処理部２２０は、これを実現するため、図２に示すように、得られた画像（拡散強調画像）上で関心領域を設定する関心領域設定部２２１と、ＭＰＧパルス３０３の印加方向を等しくし、ｂ値を変化させてパルスシーケンス３００を実行することにより得た複数の画像（拡散強調画像）の関心領域の画素の画素値を用い、画素毎に、拡散に関するパラメータである拡散パラメータを推定するパラメータ推定部２２２と、推定された拡散パラメータが、予め定めた制約条件の範囲内であるか否かを、判別する判別部２２３と、判別結果が制約条件の範囲外である画素の画素値を補正する画素値補正部２２４と、画素毎の拡散パラメータを用いて、パラメータ画像を生成する画像生成部２２５と、を備える。パラメータ推定部２２２は、制約条件無しの最小二乗フィッティングにより拡散パラメータを推定するとともに、補正後の画素値を用いて拡散パラメータを再推定し、画像生成部２２５は、画素値が補正された画素については、再推定された拡散パラメータを用いる。

まず、本実施形態の画像処理部２２０の各部による、画像処理の流れについて、説明する。ここでは、撮像部２１０が、ＤＫＩパルスシーケンスを、ＭＰＧパルス３０３の印加方向およびｂ値を変化させて実行することにより得た、複数の拡散強調画像から、パラメータ画像を得る。

図４は、本実施形態の画像処理の処理フローである。ここでは、ＭＰＧパルス３０３の印加方向数をＭ、異なるｂ値の数をＮとする（Ｍ、Ｎは２以上の整数である。以下、本明細書において同様とする）。本処理は、撮像部２１０がＤＫＩパルスシーケンスに従って、所定枚数（ここでは、Ｍ×Ｎ枚）の拡散強調画像を得た後、ユーザからの指示を受け、開始する。なお、撮像部２１０が所定枚数の拡散強調画像を得たことを契機に開始してもよい。

はじめに、関心領域設定部２２１は、取得した拡散強調画像上で関心領域を設定する（ステップＳ１１０１）。ここで、得られた関心領域内の画素数をＰとする（Ｐは、１以上の整数とする。以下、本明細書内において同様とする。）。

次に、画像処理部２２０の各部は、Ｍ×Ｎ枚の拡散強調画像の関心領域の各画素について、取得時のＭＰＧパルス３０３の印加方向毎に、以下のステップＳ１１０３からステップＳ１１０７の処理を実行する（ステップＳ１１０２、ステップＳ１１０３）。

まず、パラメータ推定部２２２は、取得する際の、ＭＰＧパルス３０３の印加方向がｍ（ｍは１以上Ｍ以下の整数）で、ｂ値の異なる複数（Ｎ枚）の拡散強調画像の、画素ｐ（ｐは１以上Ｐ以下の整数）の画素値を用いて、当該画素ｐの、拡散パラメータとして拡散係数ＡＤＣおよび尖度係数ＡＫＣを推定する（ステップＳ１１０４）。次に、判別部２２３は、推定された拡散パラメータが、予め定めた制約条件の範囲内であるか否かを判別する（ステップＳ１１０５）。

推定された拡散パラメータが、制約条件の範囲外である場合、画素値補正部２２４は、Ｎ枚の拡散強調画像の当該画素ｐの画素値を補正し（ステップＳ１１０６）、パラメータ推定部２２２は、補正後の画素値を用いて、当該画素ｐの拡散パラメータ（拡散係数ＡＤＣおよび尖度係数ＡＫＣ）を再推定する（ステップＳ１１０７）。

各ＭＰＧパルス３０３印加方向について、全関心領域の画素Ｐの拡散パラメータの推定が完了すると（ステップＳ１１０８、Ｓ１１０９）、画像生成部２２５は、各ＭＰＧパルス３０３印加方向の、画素毎の推定された拡散パラメータから、平均拡散画像および平均尖度画像などのパラメータ画像を作成する（ステップＳ１１１０）。

なお、ステップＳ１１０２より前に、全拡散強調画像の関心領域の各画素の画素値に対し、画素値補正部２２４がステップＳ１１０６で行う補正処理を施し、記憶装置１１２に格納しておいてもよい。この場合、上記ステップＳ１１０６では、補正処理を行う代わりに、記憶装置１１２から、該当する画素の補正後の画素値をリロードする。

次に、上記各部の処理について、説明する。

まず、関心領域設定部２２１による関心領域設定の手法について説明する。関心領域は、以降の画像処理演算を行う対象となる画素により構成される領域である。すなわち、拡散強調画像における背景以外の領域である。関心領域設定部２２１は、例えば、判別分析法などにより、関心領域として抽出する画素を特定する。具体的には、分離度が最大となる閾値を求め、画素値が閾値以上となる画素を関心領域として抽出する。なお、操作者が、インターフェイスを介して、閾値、あるいは、関心領域そのものを指定するよう構成してもよい。

パラメータ推定部２２２における、拡散パラメータの算出方法について説明する。算出する拡散パラメータは、上述のように、拡散係数ＡＤＣおよび尖度係数ＡＫＣである。

ここでは、ＭＰＧパルス３０３の印加方向数をＭとし、印加方向毎に連番が付されているものとする。また、ｍ番目の印加方向にＭＰＧパルス３０３が印加される場合、ＭＰＧパルス３０３の印加方向をｍと呼ぶ。また、ｂ値の数はＮとする。

ＭＰＧパルス３０３の印加方向をｍ、ｂ値をｂとしたとき、ＤＫＩ解析において、拡散強調画像の所定の画素の信号強度Ｓ（ｍ，ｂ）は、式（１）で表わされる。

ここで，Ｓ₀は、ｂ値が０の信号強度（＝Ｓ（ｍ，０））、ＡＤＣ_mは、印加方向ｍの拡散係数、ＡＫＣ_mは、印加方向ｍの尖度係数である。

パラメータ推定部２２２は、当該画素の拡散係数ＡＤＣ_mおよび尖度係数ＡＫＣ_mを、Ｎ枚の拡散強調画像の当該画素の画素値から、準ニュートン法や、Ｌｅｖｅｎｂｅｒｇ−Ｍａｒｑｕａｒｄｔ法などの制約条件なしの非線形最小二乗フィッティングにて推定する。Ｎ枚の拡散強調画像は、印加方向が等しく、ｂ値の異なるＮ個のＭＰＧパルス３０３を印加することにより得たものである。

判別部２２３は、推定された拡散パラメータ（拡散係数ＡＤＣ_m、尖度係数ＡＫＣ_m）が、制約条件を満たすか否かを判別する。制約条件は、各拡散パラメータについて、その値の範囲で定められる。制約条件の範囲は、予め定められ、記憶装置１１２に保持される。例えば、拡散係数ＡＤＣ_mについては０以上（０≦ＡＤＣ _m ）、尖度係数ＡＫＣ_mについては０以上３以下（０≦ＡＫＣ_m≦３）などの範囲とする。

これは、拡散係数ＡＤＣ_mおよび尖度係数ＡＫＣ_mは、それぞれ、物理的に、０≦ＡＤＣ_m、−２≦ＡＫＣ_m≦３の範囲で定義されるためである。従って、拡散係数ＡＤＣ_mについては、０より小さいものは、物理的に有りえない。また、尖度係数ＡＫＣ_mについても、−２より小さいもの、および、３より大きいものも、物理的に有りえない。

さらに、尖度係数ＡＫＣ_mについては、所定の画素の信号強度（画素値）Ｓ（ｍ，ｂ）の拡散モデルを、式（２）で表わされる正規分布に従うモデルと仮定すると、ｂ値の増加に伴い、拡散係数ＡＤＣ_mが低下する現象が観察される。

これは、一定以上の拡散速度をもつ組織では、ｂ値が大きくなると、核磁化の位相分散が飽和し、それ以上信号強度が減少しなくなるためである。従って、ｂ値の増加に伴って拡散係数ＡＤＣ_mが増加する状態を表す−２≦ＡＫＣ_m＜０の範囲も、尖度係数ＡＫＣ_mの範囲外と考えてよい。従って、上述の範囲が、拡散係数ＡＤＣ_mおよび尖度係数ＡＫＣ_mの満たすべき値の範囲（制約条件）として定められる。

画素値補正部２２４における画素値の補正処理について説明する。補正には、例えば、平滑化フィルタを用いる。本実施形態では、補正対象とする画素について、当該画素近傍の予め定めた範囲の画素の画素値を用いて、平滑化フィルタを適用し、補正を行う。使用する平滑化フィルタは、例えば、ガウシアンフィルタ、バイラテラルフィルタ、メジアンフィルタなどとする。あるいは、ｋ空間上で使用するローパスフィルタやフェルミフィルタ（Ｆｅｒｍｉ Ｆｉｌｔｅｒ）などを用いてもよい。

次に、画像生成部２２５におけるパラメータ画像の生成処理の一例を説明する。本実施形態の画像生成部２２５は、推定された各画素の拡散係数ＡＤＣ_mおよび尖度係数ＡＫＣ_mから、以下の手法で、パラメータ画像として、平均拡散、平均尖度、軸方向の尖度、周方向の尖度を計算し、それぞれの値を画素値とする画像を生成する。

拡散係数ＡＤＣ_mおよび尖度係数ＡＫＣ_mは、ＭＰＧパルス３０３の印加方向をｍ、２階の拡散テンソルＤと４階の尖度テンソルＷの各要素をそれぞれＤ_ij，Ｗ_ijklとすると、以下の式（３）、（４）で表わされる。

ここで、ＭＤ（Ｍｅａｎ Ｄｉｆｆｕｓｉｖｉｔｙ）は、平均拡散である。このＭＤは、拡散テンソルＤにおける固有値をλ₁、λ₂、λ₃（λ₁＞λ₂＞λ₃）とすると，以下の式（５）で定義される。

また、ＭＰＧパルス３０３の印加軸数をＭ、拡散テンソルＤの固有ベクトルを（ｅ₁，ｅ₂，ｅ₃）とし、Ｗ^T _ijkl、Ｇ₁、Ｇ₂をそれぞれ式（６）、（７）、（８）で定義すると、平均尖度ＭＫ（Ｍｅａｎ Ｋｕｒｔｏｓｉｓ）、軸方向の尖度Ｋ_A（Ａｘｉａｌ Ｋｕｒｔｏｓｉｓ）、および、周方向の尖度Ｋ_R（Ｒａｄｉａｌ Ｋｕｒｔｏｓｉｓ）は、以下の式（９）、（１０）、（１１）で表わされる。

ただし、平均尖度ＭＫ（Ｍｅａｎ Ｋｕｒｔｏｓｉｓ）、軸方向の尖度Ｋ_A（Ａｘｉａｌ Ｋｕｒｔｏｓｉｓ）、周方向の尖度Ｋ_R（Ｒａｄｉａｌ Ｋｕｒｔｏｓｉｓ）の算出方法はこれに限定されない。

本実施形態の画像生成部２２５では、上記パラメータ推定部２２２が、Ｍ個のＭＰＧパルス３０３の印加方向それぞれについて、画素毎に推定した拡散パラメータ（拡散係数ＡＤＣおよび尖度係数ＡＫＣ）を用い、上記式（３）、式（４）に従って、各テンソルの成分を計算する。テンソル成分は、連立方程式を作成し、逆行列を求めて演算することにより算出する。このとき、ステップＳ１１０５で範囲外と判別された画素については、ステップＳ１１０７で再推定された拡散パラメータを用いる。

次に，上記式（５）、式（９）、式（１０）、式（１１）にて平均拡散、平均尖度、軸方向の尖度、周方向の尖度を計算し、それぞれの値を画素値とする画像を、パラメータ画像として生成する。

生成した各パラメータ画像は、表示装置１１１に送られ操作者に提示される。

なお、本実施形態の計算機１１０は、ＣＰＵ、メモリ、記憶装置等を備え、ＣＰＵが記憶装置に格納されたプログラムをメモリにロードして実行することにより、上記各機能を実現する。各機能の全てまたはその一部は、ＭＲＩ装置１００とは独立して設けられる汎用の情報処理装置であって、ＭＲＩ装置１００とデータの送受信が可能な情報処理装置により実現されていてもよい。特に、画像処理部２２０は、画像処理装置として独立していてもよい。

以上説明したように、本実施形態のＭＲＩ装置１００は、所定のパルスシーケンスに従って、計測を行い、画像を再構成する撮像部２１０と、前記再構成された画像に対し、演算処理を施し、パラメータ画像を得る画像処理部２２０と、を備え、前記パルスシーケンスは、傾斜磁場として拡散に伴う信号変化を核磁気共鳴信号に付加するＭＰＧ（ｍｏｔｉｏｎ ｐｒｏｂｉｎｇ ｇｒａｄｉｅｎｔ）パルス３０３を印加する拡散強調撮像シーケンスであり、前記撮像部２１０は、前記ＭＰＧパルス３０３の印加方向およびｂ値を変化させて前記パルスシーケンスを実行して複数の前記画像を取得し、前記画像処理部２２０は、前記画像上で関心領域を設定する関心領域設定部２２１と、前記ＭＰＧパルス３０３の印加方向を等しくし、ｂ値を変化させて前記パルスシーケンスを実行することにより得た複数の前記画像の前記関心領域の画素の画素値を用い、画素毎に、拡散に関するパラメータである拡散パラメータを推定するパラメータ推定部２２２と、推定された前記拡散パラメータが、予め定めた範囲内であるか否かを、判別する判別部２２３と、前記判別結果が前記範囲外である画素の画素値を補正する画素値補正部２２４と、画素毎の前記拡散パラメータを用いて、前記パラメータ画像を生成する画像生成部２２５と、を備え、前記パラメータ推定部２２２は、制約条件無しの最小二乗フィッティングにより前記拡散パラメータを推定するとともに、前記補正後の画素値を用いて前記拡散パラメータを再推定し、前記画像生成部２２５は、前記画素値が補正された画素については、再推定された前記拡散パラメータを用いる。また、前記最小二乗フィッティングのモデル関数は、非線形関数とする。

すなわち、本実施形態によれば、ＤＫＩパルスシーケンスにより得た複数の拡散強調画像から、制約条件なしの非線形最小二乗フィッティングにより、拡散係数および尖度係数を推定する。そして、予め、制約条件として定めた範囲外の値が算出された画素についてのみ、平滑化フィルタなどにより、補正を行い、拡散係数および尖度係数を推定し直す。制約条件範囲から外れた画素についてのみ、得られた計測データが、ＳＮ比低下などによる計測誤差による影響を受けていると判断する。

このように、本実施形態の手法は、関心領域の全画素の画素値に対して、平滑化フィルタを適用する従来手法に比べ、平滑化処理などの補正処理を行う画素が限られるため、精度を維持したままパラメータ画像におけるボケの発生を抑制できるため、それにより得られるパラメータ画像の画質も向上する。

また、拡散パラメータの推定を、制約条件なしの非線形最小二乗フィッティングにより行うため、制約条件付きの非線形最小二乗フィッティングを用いる手法と比べても、初期値を変更する回数が最小限に抑えられるため、計算速度が向上する。

したがって、本実施形態によれば、ＤＫＩの解析において、得られる画像のボケを抑制し、かつ、計算速度を向上させることができる。

なお、ステップＳ１１０７における再推定時は、制約条件有り、無し、いずれの非線形最小二乗フィッティングを用いてもよい。制約条件有りの非線形最小二乗フィッティングを用いる場合、制約条件は、上述の、各拡散パラメータについて予め定められた制約条件を用いる。

また、制約条件有りの非線形最小二乗フィッティングを行う場合は、ステップＳ１１０４で推定された拡散パラメータを初期値に用いてもよい。このように構成することにより、非線形の最小二乗フィッティングにおける初期値に適切なものが設定されることになるため、計算を素早く収束させることができる。

＜＜第二の実施形態＞＞
次に、本発明の第二の実施形態について説明する。第一の実施形態では、拡散係数および尖度係数を、制約条件無しの非線形最小二乗フィッティングにより算出する。しかし、本実施形態では、この非線形最小二乗フィッティングの処理を改良し、拡散係数および尖度係数をより高速に算出する。

本実施形態のＭＲＩ装置は、基本的に第一の実施形態のＭＲＩ装置１００と同様の構成を有する。ただし、上述のように、パラメータ推定時の算出手法が異なる。従って、パラメータ推定部２２２の処理が異なる。以下、本実施形態について、第一の実施形態と異なる構成に主眼をおいて説明する。

本実施形態のパラメータ推定部２２２は、ＭＰＧパルス３０３の印加方向は同じで、ｂ値の異なるパルスシーケンスで得たＮ枚の拡散強調画像の当該画素の画素値から、当該画素の拡散係数ＡＤＣ_mおよび尖度係数ＡＫＣ_mを、最小二乗フィッティングにより推定する。このとき、近似する関数（モデル関数）を一次関数（線形）とする。すなわち、線形の最小二乗フィッティングにより、各係数を推定する。

近似する一次関数は、ＭＰＧパルス３０３の印加方向が等しく、ｂ値のみ異なるパルスシーケンスで得られた複数の拡散強調画像の所定の画素の画素値のペアを複数作成し、作成した複数のペアを用いて、当該画素の画素値のｂ値による変化を表す一次関数を作成する。

図５は、本実施形態のパラメータ推定部２２２のパラメータ推定処理の処理フローである。ここでは、ＭＰＧパルス３０３の印加方向数をＭとし、異なるｂ値の個数をＮ個とする。印加方向については、第一の実施形態同様、印加方向毎に連番が付されているものとし、Ｎ個のｂ値をそれぞれｂ_n（ｎ＝１，２，・・・・，Ｎ，ｂ＝０があれば含める。）とする。

ＭＰＧパルス３０３の印加方向をｍ、ｂ値をｂ_nとしてパルスシーケンスを実行し、得られた拡散強調画像の所定の画素の信号強度をＳ（ｍ，ｂ_n）（ｎ＝１，２，・・・・，Ｎ）とする。本実施形態のパラメータ推定部２２２は、以下の処理を、各画素について行う。以下、ＭＰＧパルス３０３の印加方向をｍ、ｂ値をｂ_nとしてパルスシーケンスを実行して得た拡散強調画像を、ＭＰＧパルス３０３の印加方向がｍ、ｂ値がｂ_nの拡散強調画像と呼ぶ。

はじめに、パラメータ推定部２２２は、ＭＰＧパルス３０３の印加方向がｍの、Ｎ個の拡散強調画像の、所定の画素ｐ間で、_NＣ₂組の画素値のペアを作成する（ステップＳ１２０１）。

式（１）は、１組の画素値ペア（ｂ値がｂ_sの時の信号強度（画素値）Ｓ（ｍ，ｂ _s ），および，ｂ値がｂ_tの時の信号強度（画素値）Ｓ（ｍ，ｂ _t ））を用い、以下の式（１２）のように変形できる。

本実施形態では、式（１２）で表される一次関数を、最小二乗フィッティングに用いるモデル関数とする。式（１２）より、一次関数ｙ＝Ａｘ＋Ｂにおけるデータ点（ｘ，ｙ）、傾きＡ、および切片Ｂは、以下の式（１３）から式（１６）で表される。

パラメータ推定部２２２は、作成した_NＣ₂組の画素値のペアを用いて、式（１３）、式（１４）よりフィッティング計算用のデータ点（ｘ_k，ｙ_k）（ｋ＝１，２，・・・，_NＣ₂）を計算する（ステップＳ１２０２）。

次に、パラメータ推定部２２２は、傾きＡおよび切片Ｂを式（１７）、式（１８）より算出する（ステップＳ１２０３）。

そして、パラメータ推定部２２２は、式（１７）および式（１８）により得られた傾きＡおよび切片Ｂの値を式（１５）、式（１６）に代入し、拡張係数ＡＤＣ_mおよび尖度係数ＡＫＣ_mを算出する（ステップＳ１２０４）。

なお、本実施形態の画像処理部２２０による画像処理の流れは、ステップＳ１１０４のパラメータ推定部２２２による推定処理以外は、第一の実施形態と同様である。また、本実施形態では、ステップＳ１１０７における、補正後の画素値を用いての再推定においては、上記線形最小二乗フィッティングを用いても、非線形最小二乗フィッティングを用いてもよい。また、非線形最小二乗フィッティングにおいては、制約条件の有る無し、を問わない。

以上説明したように、本実施形態のＭＲＩ装置は、第一の実施形態と同様に、拡散強調撮像シーケンスを、ＭＰＧパルス３０３の印加方向およびｂ値を変化させて実行し、複数の拡散強調画像を得る撮像部２１０と、画像処理部２２０と、を備え、画像処理部２２０は、関心領域設定部２２１と、パラメータ推定部２２２と、判別部２２３と、画素値補正部２２４と、画像生成部２２５とを備える。そして、本実施形態のパラメータ推定部２２２が最小二乗フィッティングに用いるモデル関数は、一次関数とする。

従って、本実施形態によれば、第一の実施形態同様、制約条件の範囲外の拡散パラメータ（拡張係数および尖度係数）が推定された画素についてのみ補正処理を行うため、高い品質のパラメータ画像を得ることができる。

非線形の最小二乗フィッティングは、繰り返し演算が必要となる。一方、線形の最小二乗フィッティングは、繰り返し演算が不要である。本実施形態では、この繰り返し演算の不要な線形の最小フィッティングにより、パラメータ（拡張係数および尖度係数）の推定を行う。従って、これらの演算を、さらに高速に行うことができる。従って、操作者へのパラメータ画像の提示速度を速めることができ、操作性が向上する。

なお、例えば、撮像部２１０により得られた拡散強調画像のＳＮが良い場合、画像処理部２２０による画像処理において、ステップＳ１１０５における判別処理を省いてもよい。すなわち、はじめに、関心領域設定部２２１にて計測された拡散強調画像において、画像処理を実行する画素を決定する。そして、パラメータ推定部２２２にて、上記手法により、ＭＰＧパルス印加方向毎に、画素値から拡散係数および尖度係数を算出する。最後に、画像生成部２２５にて、推定した拡散係数および尖度係数から、平均拡散画像および平均尖度画像などのパラメータ画像を作成する。この場合、補正処理がなく、平滑化フィルタを使用しない。このため、より高速に計算を完了させる効果が期待できる。

また、本実施形態においても、ステップＳ１１０７において制約条件有りの非線形最小二乗フィッティングを実行する際、ステップＳ１１０４で推定された拡散パラメータを初期値に用いてもよい。本実地形態では、ステップＳ１１０７において、線形の最小二乗フィッティングを用いることで高速にかつ真値に近い初期値を得ているため、再推定時に非線形の最小二乗フィッティングを使用したとしても、高速に演算できる。

また、このときステップＳ１１０５における判別処理を省き、ステップＳ１１０４で推定された拡散パラメータを初期値として全画素でステップＳ１１０７の制約条件有りの非線形最小二乗フィッティングを実行してもよい。これにより従来の初期値の設定なしに制約条件有りの非線形最小二乗フィッティングを用いるより繰り返し演算が減少し、高速に計算できる。

＜＜第三の実施形態＞＞
次に、本発明の第三の実施形態について説明する。本実施形態では、第二の実施形態同様、線形最小二乗フィッティングにより拡散パラメータを推定する。このとき、制約条件外とされた画素値について、第二の実施形態では、平滑化フィルタなどにより補正する。一方、本実施形態では、フィッティング関数の傾きを大きくする画素値を取り除く（間引く）ことにより、補正する。

本実施形態のＭＲＩ装置は、基本的に第二の実施形態のＭＲＩ装置と同様の構成を有する。しかし、上述のように、制約条件の範囲外とされた画素値の補正処理が異なるため、画素値補正部２２４による処理が異なる。以下、本実施形態について、第二の実施形態と異なる構成に主眼をおいて説明する。

本実施形態の画素値補正部２２４は、上述のように、最小二乗フィッティングに用いる一次関数の傾きＡを大きくすることに大きく寄与するデータを取り除く（間引く）。ＤＫＩ解析では、判別部２２３による制約条件の範囲内であるか否かの判定において、尖度係数ＡＫＣ_mが負の値（ＡＫＣ_m＜０）となることにより、制約条件の範囲外と判定される場合が比較的多い。そこで、本実施形態では、画素値補正部２２４による補正処理として、該当する画素値を特定し、それを間引く処理を実行する。間引く対象は、パラメータ推定部２２２によるパラメータ推定時に、尖度係数ＡＫＣ_mを負の値にすることに寄与するデータ点である。

以下、本実施形態の画素値補正部２２４による補正処理を説明する。図６に、本実施の形態における画素値補正部２２４による補正処理の流れを示す。

はじめに、線形の最小二乗フィッティング用に計算したデータ点（ｘ，ｙ）を用い、以下の式（１９）により、評価値Ｅを計算する（ステップＳ１３０１）。

評価値Ｅは、データ点（ｘ，ｙ）と、各データ点のｘ値の平均値ｘ₀およびｙ値の平均値ｙ₀による中心点（ｘ₀，ｙ₀）との差の内積である。この評価値Ｅは、第二の実施形態のパラメータ推定処理の線形最小二乗フィッティングにおいて、データ点（ｘ，ｙ）が、フィッティングする一次関数の傾きに与える影響の指標となる。

式（１５）より、尖度係数ＡＫＣ_mが０以上（ＡＫＣ_m≧０）となるのは、傾きＡが０以下（Ａ≦０）の場合である。式（１９）で表される評価値Ｅが大きいと、尖度係数ＡＫＣ_mが負（ＡＫＣ_m＜０）となる可能性が高い。このため、ここでは、評価値Ｅが最大となるデータ点を抽出する（ステップＳ１３０２）。

次に、抽出した各データ点とデータ点の中心点（ｘ₀，ｙ₀）との差の総和を用い、条件処理を行う（ステップＳ１３０３）。ここでは、総和が０となるか否かを判別する。そして、総和が０となる場合は補正処理を終了する。

ステップＳ１３０３において、総和が０でない場合、抽出したデータ点を間引き（ステップＳ１３０４）、残りのデータ点により、評価値Ｅを再計算する（ステップＳ１３０５）。

評価値Ｅの総和が負（総和＜０）になるか、あるいは、計算に必要な最小のデータ点数（＝２）となるまで、間引き処理（ステップＳ１３０２からステップＳ１３０５の処理）を繰り返す（ステップＳ１３０６）。ステップＳ１３０６の条件を満たした場合、間引き処理を終了する。

なお、画素値補正部２２４による補正処理以外の画像処理部２２０による画像処理は、第二の実施形態と同様である。すなわち、画像処理部２２０は、撮像部２１０がＤＫＩパルスシーケンスを実行することにより得た複数の拡散強調画像から、線形の最小二乗フィッティングを用いて拡散係数ＡＤＣ_mおよび尖度係数ＡＫＣ_mを算出する。そして、得られた結果について、制約条件の範囲内であるか否かを判別し、範囲外の画素値に対し、上記補正処理を行い、尖度係数ＡＫＣ_mを負とするものを除き、除いたデータ点から拡散係数ＡＤＣ_mおよび尖度係数ＡＫＣ_mを再推定する。そして、最終的に得られた拡散係数ＡＤＣ_mおよび尖度係数ＡＫＣ_mから、パラメータ画像を生成する。

以上説明したように、本実施形態のＭＲＩ装置は、第二の実施形態と同様に、拡散強調撮像シーケンスを、ＭＰＧパルス３０３の印加方向およびｂ値を変化させて実行し、複数の拡散強調画像を得る撮像部２１０と、画像処理部２２０と、を備え、画像処理部２２０は、関心領域設定部２２１と、パラメータ推定部２２２と、判別部２２３と、画素値補正部２２４と、画像生成部２２５とを備える。そして、本実施形態のパラメータ推定部２２２が最小二乗フィッティングに用いるモデル関数は、一次関数とする。このとき、画素値補正部２２４は、拡散パラメータを範囲外とする画素を取り除くことにより、画素値の補正を行う。

このように、本実施形態によれば、パラメータの推定に線形の最小二乗フィッティングを用いる。このため、第二の実施形態同様、高品質なパラメータ画像を高速に得ることができる。さらに、本実施形態では、補正処理において平滑化フィルタを使用しないため、平滑フィルタ使用によるパラメータ画像のボケを抑制することができる。これにより、ＤＫＩ解析において、さらに高画質な画像を提供できる。

＜＜第四の実施形態＞＞
次に、本発明の第四の実施形態を説明する。本実施形態では、制約条件の範囲外となった画素値を補正後、再推定を行ったパラメータについても、制約条件の範囲内であるか否かを判別し、範囲外である場合は、パラメータ自体を補正する。

本実施形態のＭＲＩ装置は基本的に第一の実施形態のＭＲＩ装置１００と同様の構成を有する。ただし、上記パラメータの補正処理を行うため、本実施形態の画像処理部２２０は、図７に示すように、第一の実施形態の構成に加え、パラメータ補正部２２６を備える。また、判別部２２３は、パラメータが推定される毎に、制約範囲内であるか否かを判別する。以下、本実施形態について、第一の実施形態と異なる構成に主眼をおいて説明する。

まず、本実施形態の画像処理部２２０による画像処理の流れを説明する。図８は、本実施形態の画像処理の処理フローである。第一の実施形態同様、ＭＰＧパルス３０３の印加方向数をＭ、異なるｂ値の数をＮとし、関心領域内の画素数をＰとする。

本図に示すように、関心領域を設定し、設定した関心領域内で、ＭＰＧパルス３０３の印加方向毎に、画素毎に、拡散パラメータ（拡散係数ＡＤＣ_mおよび尖度係数ＡＫＣ_m）を算出し、制約条件の範囲内であるか否かを判別し、範囲外であれば、画素値を補正する処理までは、第一の実施形態と同様である。

本実施形態では、パラメータ推定部２２２は、補正後の画素値を用いて、拡散パラメータを再推定する（ステップＳ１４０１）。このとき、本実施形態のパラメータ推定部２２２は、制約条件無しの非線形最小二乗フィッティングにより、拡散パラメータの再推定を行う。

その後、判別部２２３は、再推定された拡散パラメータが、制約条件の範囲内であるか否かを判別する（ステップＳ１４０２）。そして、範囲内であれば、ステップＳ１１０８へ移行する。一方。範囲外であれば、パラメータ補正部２２６は、パラメータ自体を補正する（ステップＳ１４０３）。補正方法は、例えば、制約条件内で当該拡散パラメータに最も近い値に置き換えるものとする。

以後の処理は、第一の実施形態と同様である。

ステップＳ１４０３のパラメータ補正部２２６によるパラメータ補正処理の流れを図９に示す。まず、推定された拡散係数ＡＤＣ_mの値が０未満（ＡＤＣ_m＜０）である場合（ステップＳ１４１１）、拡散係数ＡＤＣ_mの値を０に置き換える（ステップＳ１４１２）。推定された尖度係数ＡＫＣ_mの値が０未満（ＡＫＣ_m＜０）である場合（ステップＳ１４１３）、尖度係数ＡＫＣ_mの値を０に置き換える（ステップＳ１４１４）。推定された尖度係数ＡＫＣ_mの値が３より大きい（ＡＫＣ_m＞３）場合（ステップＳ１４１５）、尖度係数ＡＫＣ_mの値を３に置き換える（ステップＳ１４１６）。

なお、本処理において、上記各判別順序は問わない。

なお、パラメータ補正部２２６による補正方法は、上記手法に限定されない。再推定した前記拡散パラメータを、予め定めた補正値だけ補正するものとしてもよい。このとき、その補正値は、各画素について算出する全ての拡散パラメータが制約条件の範囲内に収まるよう決定される。

一般に、式（１）における信号強度（画素値）Ｓ（ｍ，ｂ）は、ｂ値が増加するに伴って単調減少する。しかしながら、ｂ値に応じた信号強度Ｓ（ｍ，ｂ）の変化を、下に凸の二次関数モデルで近似すると、近似にあたり、計測に使用するｂ値の範囲内で極値をとることがある。このような場合、ｂ値が増加するに従って、信号強度Ｓ（ｍ，ｂ）が増加に転じる。本実施形態のパラメータ補正部２２６では、このような拡散係数および尖度係数が算出されないようにパラメータを補正する。

変化させるｂ値の中で最大のｂ値をｂ_maxとすると、ｂ_max以下で信号強度Ｓ（ｍ，ｂ）が極値を持たない、すなわち、信号強度Ｓ（ｍ，ｂ）が極値を取るときのｂ値が、最大ｂ値ｂ_max以上となればよい。従って、拡散係数ＡＤＣ_mおよび尖度係数ＡＫＣ_mは、以下の式（２０）で制限される。

制約条件として、式（２０）を使用する場合、上記のパラメータ補正処理のように単純に値を置き換えることができない。従って、この場合、パラメータ補正部２２６は、制約条件範囲内に値が収まるように、拡散係数ＡＤＣ_mおよび尖度係数ＡＫＣ_mを、同じ補正値Δで補正する。

補正値Δは、以下の式（２１）の解のうち、小さい方の値とする。具体的には、以下の式（２２）で与えられる。

以上説明したように、本実施形態のＭＲＩ装置は、第一の実施形態と同様に、拡散強調撮像シーケンスを、ＭＰＧパルス３０３の印加方向およびｂ値を変化させて実行し、複数の拡散強調画像を得る撮像部２１０と、画像処理部２２０と、を備え、画像処理部２２０は、関心領域設定部２２１と、パラメータ推定部２２２と、判別部２２３と、画素値補正部２２４と、画像生成部２２５とを備える。
そして、本実施形態のＭＲＩ装置は、再推定された拡散パラメータを補正するパラメータ補正部２２６をさらに備え、前記パラメータ推定部２２２は、制約条件無しの最小二乗フィッティングにより、補正後の画素値から前記拡散パラメータを再推定し、前記判別部２２３は、再推定された前記拡散パラメータが、制約条件の範囲内であるか否かを判別し、前記パラメータ補正部２２６は、前記判別結果が範囲外である場合、再推定した前記拡散パラメータを補正する。

本実施形態によれば、上記各実施形態同様、高速で高品質のパラメータ画像を提供できる。また、本実施形態によれば、画素値補正後に再推定されたパラメータが、制約条件の範囲外であっても、それを補正し、制約条件内に収める。従って、補正後の画素値によるパラメータの再推定に、高速な処理を用いることができる。従って、さらに高速に画像処理を行うことができる。

なお、本実施形態では、ステップＳ１１０４のパラメータ推定処理およびステップＳ１４０１のパラメータ再推定処理において、第二の実施形態のように線形の最小二乗フィッティングを用いてもよい。また、ステップＳ１１０６の補正処理において、第三の実施形態の補正処理を行ってもよい。

＜＜第五の実施形態＞＞
次に、本発明の第五の実施形態を説明する。本実施形態では、得られた拡散強調画像の各画素値の信頼性を算出し、ユーザに提示する。

本実施形態のＭＲＩ装置は基本的に第一の実施形態のＭＲＩ装置１００と同様の構成を有する。ただし、本実施形態の画像処理部２２０は、上記処理を行うため、図１０に示すように、第一の実施形態の構成に加え、パラメータ画像算出時の計算の信頼性を算出する信頼度演算部２２７を備える。以下、本実施形態について、第一の実施形態と異なる構成に主眼をおいて説明する。

画像処理部２２０において、ＤＫＩパルスシーケンスによって得られた複数の拡散強調画像から、拡散パラメータ（拡散係数および尖度係数）を推定し、それを用いてパラメータ画像を生成する処理の流れは、第一の実施形態と同様である。

本実施形態の信頼度演算部２２７は、推定された前記拡散パラメータの中の、前記範囲内と判別されたパラメータ数を算出し、算出結果から得た情報を信頼度表示としてユーザに提示する。信頼度として、例えば、各ＭＰＧパルス３０３の印加方向の拡散係数および尖度係数の中で、画素値の補正を行わずに計算された画素数の割合を算出する。算出は、ステップＳ１１０５において、制約条件の範囲内の画素をカウントすることにより行う。ステップＳ１１０５において、制約条件の範囲内であれば、信頼度演算部２２７は、当該ＭＰＧパルス３０３印加方向の拡散強調画像の、当該画素について、制約条件内であることを示すフラグを設定し、フラグが設定された画素数をカウントする。

本実施形態の画像処理部２２０による画像処理の流れを図１１に示す。第一の実施形態同様、ＭＰＧパルス３０３の印加方向数をＭ、異なるｂ値の数をＮとし、関心領域の画素数をＰとする。

本図に示すように、パラメータ推定部２２２は、関心領域設定部２２１により関心領域として設定された領域内で、ＭＰＧパルス３０３の印加方向毎に、画素毎に、パラメータ（拡散係数ＡＤＣ_mおよび尖度係数ＡＫＣ_m）を推定し、判別部２２３は、推定されたパラメータが制約条件の範囲内であるか否かを判別する。

そして、制約条件範囲内であれば、信頼度演算部２２７は、当該ＭＰＧパルス印加方向ｍの、当該画素ｐの信頼度フラグＦｒ（ｍ，ｐ）を設定する（ステップＳ１５０１）。そして、ステップＳ１１０８へ移行する。

なお、ステップＳ１１０５において、制約条件範囲外であった場合の処理は、第一の実施形態と同様である。

そして、全関心領域の画素、全ＭＰＧパルス３０３印加方向について、パラメータの推定を終えると、第一の実施形態同様、画像生成部２２５は、推定した拡散係数ＡＤＣおよび尖度係数ＡＫＣから、平均拡散画像および平均尖度画像などのパラメータ画像を作成する。本実施形態では、このとき、信頼度演算部２２７は、信頼度表示を作成する（ステップＳ１５０２）。作成されたパラメータ画像および信頼度表示は、画像処理部２２０により表示装置１１１に表示される。

本実施形態では、信頼度表示として、例えば、関心領域の各画素ｐについて、制約条件内の画素数の、全ＭＰＧパルス印加方向数Ｍに対する割合を算出する。なお、算出する信頼度は、画素数そのものであってもよい。このとき、算出された割合（あるいは、画素数）毎に予め対応づけた色などで表示してもよい。

以上説明したように、本実施形態のＭＲＩ装置は、第一の実施形態と同様に、拡散強調撮像シーケンスを、ＭＰＧパルス３０３の印加方向およびｂ値を変化させて実行し、複数の拡散強調画像を得る撮像部２１０と、画像処理部２２０と、を備え、画像処理部２２０は、関心領域設定部２２１と、パラメータ推定部２２２と、判別部２２３と、画素値補正部２２４と、画像生成部２２５とを備える。
そして、推定された拡散パラメータの中の、制約条件範囲内と判別されたパラメータ数を算出し、算出結果から得た情報を、信頼度表示としてユーザに提示する信頼度演算部２２７をさらに備える。

このように、本実施形態によれば、上記各実施形態同様、制約条件無しの最小二乗フィッティングにより拡散パラメータを推定し、制約条件外のものについてのみ補正を行う。このため、上記各実施形態と同様の効果を奏する。さらに、本実施形態では、補正無しで得られた画素数の割合を算出してユーザに示す。

上述のように、例えば、平滑化フィルタなどを用いて画素値の補正を行うことにより、画質は劣化する。関心領域内の各画素について、補正無しで、制約条件の範囲内の拡散パラメータを推定できた画素数の割合を示すことにより、ユーザは、パラメータ画像における、各画素の劣化の度合いを、大まかに把握できる。すなわち、計測データに対する計測誤差の影響を把握できる。

従って、本実施形態によって、得られたパラメータ画像を用いて診断するにあたり、画素毎の信頼度を把握できるため、所望の部位が、どの程度の信頼度で計算されているかを考慮しながら診断を行うことができ、診断能が向上する。

なお、ＭＰＧパルス３０３印加方向毎に、推定された拡散パラメータが制約条件の範囲内となる画素数を算出し、関心領域内の全画素数Ｐに対する割合を信頼度として算出してもよい。ここでも、画素数そのものを信頼度としてもよい。このとき、当該ＭＰＧパルス印加方向ｍの、補正無しの画素数をカウントするカウンタＣＴ_mを設け、上記ステップＳ１５０１において、フラグを設定するとともに、あるいは、フラグの設定に替えて、当該印加方向ｍの、補正無し画素数をカウントするよう構成してもよい。

さらに、ＭＰＧパルス３０３印加方向毎に、かつ、画素毎に、推定された拡散パラメータが制約条件範囲内であるか否かの情報を、信頼度として算出し、表示してもよい。

また、逆に、制約条件の範囲外と判別された画素数をカウントするよう構成してもよい。この場合、ステップＳ１１０５において、制約条件の範囲外と判別された場合、次の画素の処理に移るまでの間に、当該ＭＰＧパルス印加方向の拡散強調画像の、当該画素について、制約条件外であることを示すフラグを設定する。そして、このフラグを用い、各画素の、補正を行った印加方向数、あるいは、印加方向毎の、補正を行った画素数を算出し、信頼度表示としてもよい。制約条件内、制約条件外、それぞれにフラグを設定してもよい。

なお、本実施形態の手法は、第一の実施形態に限定されず、第二から第三のいずれの実施形態に組み合わせてもよい。

１００：ＭＲＩ装置、１０１：マグネット、１０２：傾斜磁場コイル、１０３：被検体、１０４：シーケンサ、１０５：傾斜磁場電源、１０６：高周波磁場発生器、１０７：ＲＦコイル、１０８：ＲＦプローブ、１０９：受信器、１１０：計算機、１１１：表示装置、１１２：記憶装置、１１３：シムコイル、１１４：シム電源、１１５：寝台、１１６：入力装置、２１０：撮像部、２２０：画像処理部、２２１：関心領域設定部、２２２：パラメータ推定部、２２３：判別部、２２４：画素値補正部、２２５：画像生成部、２２６：パラメータ補正部、２２７：信頼度演算部、３００：ｓｓ−ＤＷＥＰＩパルスシーケンス、３０１：傾斜磁場パルス、３０２：ＲＦパルス、３０３：ＭＰＧパルス、３０４：スライス傾斜磁場パルス、３０５：磁化リフォーカス用ＲＦパルス、３０６：ＭＰＧパルス、３０７：位相エンコード傾斜磁場パルス、３０８：ディフェーズ用リードアウト傾斜磁場パルス、３０９：リードアウト傾斜磁場パルス、３１０：エコー、３１１：位相エンコード傾斜磁場パルス、３１２：リードアウト傾斜磁場パルス

Claims (15)

1. 所定のパルスシーケンスに従って、計測を行い、画像を再構成する撮像部と、
   前記再構成された画像に対し、演算処理を施し、パラメータ画像を得る画像処理部と、を備え、
   前記パルスシーケンスは、傾斜磁場として拡散に伴う信号変化を核磁気共鳴信号に付加するＭＰＧ（ｍｏｔｉｏｎ ｐｒｏｂｉｎｇ ｇｒａｄｉｅｎｔ）パルスを印加する拡散強調撮像シーケンスであり、
   前記撮像部は、前記ＭＰＧパルスの印加方向およびｂ値を変化させて前記パルスシーケンスを実行して複数の前記画像を取得し、
   前記画像処理部は、
   前記ＭＰＧパルスの印加方向を等しくし、ｂ値を変化させて前記パルスシーケンスを実行することにより得た複数の前記画像の画素の画素値を用い、画素毎に、拡散に関するパラメータである拡散パラメータを推定するパラメータ推定部と、
   推定された前記拡散パラメータが、予め定めた範囲内であるか否かを、判別する判別部と、
   前記判別結果が前記範囲外である前記画素の画素値を補正する画素値補正部と、
   画素毎の前記拡散パラメータを用いて、前記パラメータ画像を生成する画像生成部と、を備え、
   前記パラメータ推定部は、制約条件無しの最小二乗フィッティングにより前記拡散パラメータを推定するとともに、前記判別結果が前記範囲外の前記画素については、前記補正後の画素値を用いて前記拡散パラメータを再推定し、
   前記画像生成部は、画素値が補正された前記画素については、再推定された前記拡散パラメータを用いること
   を特徴とする磁気共鳴イメージング装置。
2. 請求項１記載の磁気共鳴イメージング装置であって、
   前記最小二乗フィッティングのモデル関数は、非線形関数であること
   を特徴とする磁気共鳴イメージング装置。
3. 請求項１記載の磁気共鳴イメージング装置であって、
   前記最小二乗フィッティングのモデル関数は、一次関数であること
   を特徴とする磁気共鳴イメージング装置。
4. 請求項３記載の磁気共鳴イメージング装置であって、
   前記拡散パラメータが、拡散係数および尖度係数であり、
   前記画素値補正部は、前記パラメータ推定部によるパラメータの推定を行う際に用いる画素のうち、前記尖度係数が前記範囲外であり、かつ、負の値を示す画素を取り除くことにより、前記補正を行い、
   前記パラメータ推定部が、前記画素値補正部により取り除かれた画素以外の画素を用いて前記拡散パラメータを再推定すること
   を特徴とする磁気共鳴イメージング装置。
5. 請求項１記載の磁気共鳴イメージング装置であって、
   再推定された前記拡散パラメータを補正するパラメータ補正部をさらに備え、
   前記パラメータ推定部は、制約条件無しの最小二乗フィッティングにより、前記補正後の画素値から前記拡散パラメータを再推定し、
   前記判別部は、再推定された前記拡散パラメータが前記範囲内であるか否かを判別し、
   前記パラメータ補正部は、前記判別結果が前記範囲外である場合、再推定した前記拡散パラメータを補正すること
   を特徴とする磁気共鳴イメージング装置。
6. 請求項５記載の磁気共鳴イメージング装置であって、
   前記パラメータ補正部による補正は、再推定した前記拡散パラメータを前記範囲内の最も近い値に置き換えるものであること
   を特徴とする磁気共鳴イメージング装置。
7. 請求項５記載の磁気共鳴イメージング装置であって、
   前記パラメータ補正部による補正は、再推定した前記拡散パラメータを予め定めた補正値だけ補正するものであり、
   前記補正値は、１つの画素について算出する全ての拡散パラメータが前記範囲内に収まるよう決定されること
   を特徴とする磁気共鳴イメージング装置。
8. 請求項１記載の磁気共鳴イメージング装置であって、
   推定された前記拡散パラメータの中の、前記範囲内と判別されたパラメータ数を算出し、算出結果から得た情報を信頼度表示としてユーザに提示する信頼度演算部をさらに備えること
   を特徴とする磁気共鳴イメージング装置。
9. 請求項１記載の磁気共鳴イメージング装置であって、
   前記拡散パラメータは、拡散係数および尖度係数であり、
   前記パラメータ画像は、平均拡散、平均尖度、軸方向の尖度、周方向の尖度のいずれかを画素値とする画像であること
   を特徴とする磁気共鳴イメージング装置。
10. 請求項１記載の磁気共鳴イメージング装置であって、
    前記画素値補正部は、平滑化フィルタにより前記補正を行うこと
    を特徴とする磁気共鳴イメージング装置。
11. 請求項１記載の磁気共鳴イメージング装置であって、
    前記パラメータ推定部は、再推定時は、推定した当該画素の前記拡散パラメータを初期値とし、制約条件有りの最小二乗フィッティングにより前記拡散パラメータを再推定すること
    を特徴とする磁気共鳴イメージング装置。
12. 請求項１記載の磁気共鳴イメージング装置であって
    前記画像処理部は、さらに前記画像上で関心領域を設定する関心領域設定部を有し、
    前記パラメータ推定部は、前記関心領域の画素の画素値を用いること
    を特徴とする磁気共鳴イメージング装置。
13. 傾斜磁場として拡散に伴う信号変化を核磁気共鳴信号に付加するＭＰＧ（ｍｏｔｉｏｎ ｐｒｏｂｉｎｇ ｇｒａｄｉｅｎｔ）パルスを印加する拡散強調撮像シーケンスを、前記ＭＰＧパルスの印加方向を等しくし、ｂ値を変化させて実行することにより磁気共鳴イメージング装置において得た複数の画像の画素の画素値を用い、画素毎に、拡散に関するパラメータである拡散パラメータを推定するパラメータ推定部と、
    推定された前記拡散パラメータが、予め定めた範囲内であるか否かを、判別する判別部と、
    判別結果が前記範囲外である前記画素の画素値を補正する画素値補正部と、
    画素毎の前記拡散パラメータを用いて、パラメータ画像を生成する画像生成部と、を備え、
    前記パラメータ推定部は、制約条件無しの最小二乗フィッティングにより前記拡散パラメータを推定するとともに、前記補正後の画素値を用いて前記拡散パラメータを再推定し、
    前記画像生成部は、画素値が補正された前記画素については、再推定された前記拡散パラメータを用いること
    を特徴とする画像処理装置。
14. 傾斜磁場として拡散に伴う信号変化を核磁気共鳴信号に付加するＭＰＧ（ｍｏｔｉｏｎ ｐｒｏｂｉｎｇ ｇｒａｄｉｅｎｔ）パルスを印加する拡散強調撮像シーケンスを、前記ＭＰＧパルスの印加方向を等しくし、ｂ値を変化させて実行することにより磁気共鳴イメージング装置において得た複数の画像の画素の画素値を用い、画素毎に、制約条件無しの最小二乗フィッティングにより拡散に関するパラメータである拡散パラメータを推定し、
    推定された前記拡散パラメータが、予め定めた範囲内であるか否かを、判別し、
    判別結果が前記範囲外である前記画素の画素値を補正し、
    前記補正後の画素値を用いて前記拡散パラメータを再推定し、
    画素値が補正された前記画素については、再推定された前記拡散パラメータを、その他の前記画素については、推定された前記拡散パラメータを用いて、前記パラメータ画像を生成すること
    を特徴とする画像処理方法。
15. 傾斜磁場として拡散に伴う信号変化を核磁気共鳴信号に付加するＭＰＧ（ｍｏｔｉｏｎ ｐｒｏｂｉｎｇ ｇｒａｄｉｅｎｔ）パルスを印加する拡散強調撮像シーケンスを、前記ＭＰＧパルスの印加方向を等しくし、ｂ値を変化させて実行することにより磁気共鳴イメージング装置において得た複数の画像の画素の画素値を用い、画素毎に、制約条件無しの最小二乗フィッティングにより拡散に関するパラメータである拡散パラメータを推定し、
    推定された前記拡散パラメータが、予め定めた範囲内であるか否かを、判別し、
    判別結果が前記範囲外である前記画素の画素値を補正し、
    前記補正後の画素値を用いて前記拡散パラメータを再推定し、
    画素値が補正された前記画素については、再推定された前記拡散パラメータを、その他の前記画素については、推定された前記拡散パラメータを用いて、前記パラメータ画像を生成すること
    をコンピュータに実行させる画像処理プログラム。

Images