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

Description

本発明は、磁気共鳴イメージング技術に関し、特に、計算によって被検体パラメータを推定する方法に関する。

磁気共鳴イメージング（ＭＲＩ）装置は、被検体を横切る任意の平面内の水素原子核に核磁気共鳴を起こさせ、発生する核磁気共鳴信号からその平面内における断層像を撮影する医用画像診断装置である。一般的には、撮影面を特定するスライス傾斜磁場を印加すると同時にその面内の磁化を励起させる励起パルス（高周波磁場パルス）を与え、これにより励起された磁化が収束する段階で発生する核磁気共鳴信号（エコー）を得る。磁化に位置情報を与えるため、励起からエコーを得るまでの間に、断層面内で互いに垂直な方向の位相エンコード傾斜磁場とリードアウト傾斜磁場を印加する。

エコーを発生させるためのパルスと各傾斜磁場は、あらかじめ設定されたパルスシーケンスに基づいて印加されるようになっている。このパルスシーケンスは、目的に応じて種々のものが知られている。例えば、グラディエントエコー（ＧＥ）タイプの高速撮影法は、そのパルスシーケンスを繰り返して作動させ、繰り返しごとに位相エンコード傾斜磁場を順次変化させることにより、１枚の断層像を得るために必要な数のエコーを順次計測していく方法である。

ＭＲＩ装置では、このように計測したエコー信号に、フーリエ変換等の演算を施すことにより被検体の画像を再構成する。再構成された画像の画素値（ピクセル値）は、組織のプロトン密度、縦緩和時間Ｔ１や横緩和時間Ｔ２などの被検体自体の特性の他に、上述したパルスシーケンスの種類やその繰り返し時間（ＴＲ）や高周波磁場パルスである励起パルスの強度及び位相などの撮影条件（撮影条件の個々の要素を撮影パラメータという）、及び、磁場強度などの装置側の条件によって決まる。ピクセル値を決めるこれらの各要素は、被検体、撮影条件及び装置のそれぞれに分けて、被検体パラメータ、撮影パラメータ、装置パラメータと呼ばれる。

このように再構成画像自体は、被検体パラメータ等の値（定量値）をそのまま表すものではないが、撮影パラメータ、被検体パラメータ、装置パラメータとピクセル値の関係を表す関数(信号関数)がわかっていれば、それを用いた計算によって定量値を推定し、定量値を画素値とする画像（計算画像）を得ることが可能である。この信号関数は撮影シーケンスに依存する。特許文献１には、信号関数を数値シミュレーションによって求め、複数の元画像の撮影パラメータとそのピクセル値に対する信号関数の最小二乗フィットを求めることで、被検体パラメータあるいは装置パラメータの値をピクセル値とする計算画像を算出する手法が提案されている。

この方法では、被検体パラメータと装置パラメータの定量値推定誤差を最小にするため、誤差伝播の法則を用いて複数の元画像の撮影条件（撮影パラメータ）を決定する。誤差伝播の法則を用いた撮影パラメータの最適化は、与えられた撮影パラメータの全ての組合せについて、対象とする被検体パラメータあるいは装置パラメータに対する定量値推定誤差を計算し、それが最小になる組合せを探索する。

特開２０１１－０２４９２６号公報

特許文献１の誤差伝播の法則を用いた撮影パラメータの決定の過程では、特定の撮影パラメータにおいて画像に発生するアーチファクトについては考慮されていない。そのため、血流や脳脊髄液の流れ、体の動きなどの影響によってアーチファクトが発生しやすい撮影パラメータが選ばれることがある。そうした場合、元画像にアーチファクトが発生してしまい、そこから計算された被検体パラメータと装置パラメータにもアーチファクトが発生したり、アーチファクトがノイズのように作用してＳＮ比が低下したりしてしまうことがある。また、流れや動きによるアーチファクトは一定ではないため、撮影ごとに定量値が変動してしまうこともある。

本発明は、上記事情に鑑みてなされたもので、被検体パラメータや装置パラメータの計算画像を生成する際に、流れや体動のアーチファクトの影響を抑えて安定的に定量値マップを取得する技術を提供することを課題とする。

流れや体動によるアーチファクトは、撮影条件によって発生しやすさに違いがある。本発明では、あらかじめどの撮影条件でアーチファクトが発生しやすいかを調べておき、アーチファクトが発生しやすい撮影条件を除いて定量値マップ算出用の撮影パラメータセットの最適化を実施する。

具体的には、本発明のＭＲＩ装置は、予め定めた撮影条件および予め定めた撮影シーケンスに従って、静磁場の中に置かれた被検体に高周波磁場および傾斜磁場を印加して、前記被検体から発生するエコー信号を計測する計測部と、計測した前記エコー信号から再構成画像を作成する画像再構成部と、撮影条件を変えて前記撮影シーケンスを実行してそれぞれ得られた再構成画像と撮影シーケンス毎に定まる信号関数とを用いて被検体に依存する被検体パラメータあるいは装置に依存する装置パラメータの分布を推定するパラメータ推定部と、パラメータ推定部が用いる複数の再構成画像を取得するための複数の撮影条件の組み合わせを決定する撮影条件決定部と、を備える。撮影条件決定部は、再構成画像の画質を評価し、複数の撮影条件の組み合わせから除外すべき撮影条件の組み合わせを探索する撮影条件探索部を備え、当該除外すべき撮影条件の組み合わせを除外して、複数の撮影条件の組み合わせを決定する。

上述した撮影条件決定部の機能は、ＭＲＩ装置本体に備えることも可能であるが、機能の一部をＭＲＩ装置とは別の画像処理装置で実行することも可能である。

また本発明の画像処理方法は、ＭＲＩ装置で異なる撮影条件で取得した複数の再構成画像を処理し、被検体に依存する被検体パラメータあるいは装置に依存する装置パラメータの分布を推定する画像処理方法であって、撮影条件毎に再構成画像の画質を評価し、ＭＲＩ装置によるパラメータ推定に用いる撮影条件を探索するステップを含む。

本発明によれば、撮影条件（撮影パラメータ）毎に得られる再構成画像の画質も考慮して撮影パラメータを決定することにより、アーチファクトが発生しにくい撮影パラメータだけを用いて最適化を実施することができる。これによって安定的に定量値マップを取得できる。

本発明の実施形態のＭＲＩ装置の概略構成を示すブロック図である。 本発明の実施形態の計算機の機能ブロック図である。 本発明の実施形態における処理フローである。 本発明の実施形態におけるシーケンス図（a）と、k空間データ（b）を説明する図である。 本発明の実施形態における信号関数の一部を示す図である。 本発明の実施形態における撮影パラメータセット最適化対象の組織とそのＴ１値、T2値である。 撮影パラメータ探索処理の詳細を示す処理フローである。 本発明の実施形態における撮影条件による画質の違いを説明する図である。 本発明の実施形態による撮影条件の探索を説明する図である。 （Ａ）、（Ｂ）は、本発明の実施形態におけるマップの画質改善効果を説明する図である。 （Ａ）、（Ｂ）は、本発明の実施形態におけるマップの画質改善効果を説明する図である。

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

まず、本実施形態のＭＲＩ装置について説明する。図１は、本実施形態のＭＲＩ装置１０の概略構成を示すブロック図である。ＭＲＩ装置１０は、大きく分けて、被検体から核磁気共鳴信号であるエコー信号を計測する計測部１００と、エコー信号を用いた演算や計測部１００の制御を行う計算機１１０とを備える。また計算機１１０に付随する装置として、ディスプレイ１２０、入力デバイス１３０、及び記憶媒体１４０などを備える。

計測部１００は、予め定めた撮影条件および予め定めた撮影シーケンスに従って、静磁場の中に置かれた被検体に高周波磁場および傾斜磁場を印加して、被検体から発生するエコー信号を計測するものであり、静磁場を発生するマグネット１０１と、傾斜磁場を発生する傾斜磁場コイル１０２と、シーケンサ１０４と、傾斜磁場電源１０５と、高周波磁場発生器１０６と、高周波磁場を照射するとともに核磁気共鳴信号を検出する送受信コイル１０７と、受信器１０８とを備える。送受信コイル１０７として、図では単一のものを示しているが送信コイルと受信コイルとを別個に備えていてもよい。被検体（例えば、生体）１０３はマグネット１０１の発生する静磁場空間内の寝台（テーブル）に載置され、撮影が行われる。

シーケンサ１０４は、通常、予めプログラムされたタイミング、強度で各装置が動作するように制御を行う。プログラムのうち、特に、高周波磁場、傾斜磁場、信号受信のタイミングや強度を記述したものはパルスシーケンスと呼ばれる。シーケンサ１０４は、設定されたパルスシーケンスとユーザが任意に設定した撮影条件とに従って、撮影シーケンスを計算し、撮影シーケンスに従って、傾斜磁場電源１０５と高周波磁場発生器１０６に命令を送り、それぞれ傾斜磁場および高周波磁場を発生させる。

高周波磁場発生器１０６から発生された高周波磁場は、送受信コイル１０７を通じて被検体１０３に印加される。被検体１０３から発生した核磁気共鳴信号は送受信コイル１０７によって受波され、受信器１０８で検波が行われる。検波の基準とする核磁気共鳴周波数（検波基準周波数ｆ０）は、シーケンサ１０４によりセットされる。検波された信号は、計算機１１０に送られ、ここで画像再構成などの信号処理が行われる。その結果は、ディスプレイ１２０に表示される。必要に応じて、記憶媒体１４０に検波された信号や撮影条件を記憶させることもできる。

計算機１１０は、ＣＰＵとメモリとを備え、計測部１００の各部の動作を制御するとともに、得られたエコー信号に対し、各種の信号処理を施し、所望の画像を得る。これらの処理を実現するため、本実施形態の計算機１１０は、図２に示すように、制御部２１０と、演算部２２０とを備える。計算機１１０が実現する各機能は、例えば、記憶媒体１４０に格納されたプログラムを、計算機１１０のＣＰＵがメモリにロードして実行することにより実現される。

本実施形態のＭＲＩ装置は、撮影シーケンスによって決まる信号関数を用いて定量値を得て、得られた定量値から計算画像を得る。このため、演算部２２０は、計測したエコー信号から再構成画像を得る画像再構成部２２１と、数値シミュレーションによって、撮影シーケンス毎の信号関数を生成する信号関数生成部２２２と、撮影シーケンス毎の信号関数を用いて定量値（被検体パラメータや装置パラメータ）を推定し、計算画像（パラメータ画像）を生成するパラメータ推定部２２３と、パラメータ推定のための撮影条件を決定する撮影条件決定部２２４とを備える。撮影条件決定部２２４は、多数の撮影条件の組み合わせから除外すべき撮影条件を探索する撮影条件探索部２２５を含む。

被検体パラメータには、縦緩和時間（Ｔ１）、横緩和時間（Ｔ２）、スピン密度（ρ）、共鳴周波数差（Δｆ０）、拡散係数（ｂ）、高周波磁場の照射強度分布（Ｂ１）などがある。共鳴周波数差Δｆ０は、各ピクセルの共鳴周波数と基準周波数ｆ０との差である。装置パラメータには、静磁場強度（Ｂ０）、受信コイルの感度分布（Ｓｃ）などがある。ここでは、被検体パラメータと装置パラメータとを総称して、定量値という。また、本実施形態では、ピクセル毎に得た定量値の分布を、計算画像或いはマップと呼ぶ。

また、撮影シーケンス実行時に、ユーザが任意に設定可能なパラメータを撮影条件あるいは撮影パラメータと呼ぶ。撮影条件には、例えば、繰り返し時間（ＴＲ）、エコー時間（ＴＥ）、高周波磁場の設定強度（フリップ角（Ｆｌｉｐ Ａｎｇｌｅ：ＦＡ））、高周波磁場の位相（位相増分値θ）などがある。

なお、図２に示す計算機１１０の諸機能の一部は、ＰＬＤ（ｐｒｏｇｒａｍｍａｂｌｅ ｌｏｇｉｃ ｄｅｖｉｃｅ）等のハードウエアで実現してもよい。また、信号関数生成部２２２及びパラメータ推定部２２３は、ＭＲＩ装置１００とは独立に設けられた計算機であって、ＭＲＩ装置１００の計算機１１０とデータの送受信が可能な計算機上に構築することも可能である。

以上の構成を踏まえ、本実施形態のＭＲＩ装置１０の処理、主としてパラメータ推定に関わる処理の実施形態を説明する。

＜実施形態１＞
本実施形態では、撮影パラメータとして、ＦＡ、ＴＲ及びθの組み合わせをパラメータセットとして、パラメータセットを構成する１以上のパラメータの値を異ならせた複数のパラメータセットを用いて、定量値の推定（パラメータ推定）を行う場合を説明する。

本実施形態の処理の流れを図３に示す。図示するように、数値シミュレーションによる信号関数の生成（Ｓ３０１）、撮影条件決定部２２４による複数の撮影パラメータセットの決定（Ｓ３０２）、計測部１１０による撮影パラメータセットの条件での信号計測（Ｓ３０３）、画像再構成部２２１による計測データを用いた画像再構成（Ｓ３０４）、再構成された画像（元画像）を用いたパラメータ推定（Ｓ３０５）、及び推定されたパラメータを用いた画像生成（Ｓ３０６）、出力（Ｓ３０７）を行う。

以下、各ステップの詳細を説明する。

［Ｓ３０１：信号関数の生成］
まず、信号関数生成部２２２が数値シミュレーションによって信号関数３２を作成する。信号関数３２は、撮影に用いるパルスシーケンスに依存する。ここではパルスシーケンスとして、GEタイプのRF-spoiled GEシーケンスを用いる。図４（ａ）に、RF-spoiled GEシーケンスを表すタイミングチャートを示す。本図において、ＲＦ、Ｇｓ、Ｇｐ、Ｇｒはそれぞれ、高周波磁場、スライス傾斜磁場、位相エンコード傾斜磁場、リードアウト傾斜磁場を表す。

このパルスシーケンスでは、まず、スライス傾斜磁場パルス401の印加とともに高周波磁場（ＲＦ）パルス402を照射し、対象物体内のあるスライスの磁化を励起する。次いでスライスリフェーズ傾斜磁場パルス403と磁化の位相に位相エンコード方向の位置情報を付加するための位相エンコード傾斜磁場パルス404、ディフェーズ用リードアウト傾斜磁場405を印加した後、リードアウト方向の位置情報を付加するためのリードアウト傾斜磁場パルス406を印加しながら磁気共鳴信号（エコー）407を計測する。そして最後にディフェーズ用位相エンコード傾斜磁場パルス409を印加する。

以上の手順を位相エンコード傾斜磁場パルス404と409の強度（位相エンコード量kp）を変化させるとともに、ＲＦパルスの位相増分値θを117度としてRF位相を変化させながら (n番目のＲＦパルスの位相はθ(n) = θ(n-1) + 117nとなる)、繰り返し時間ＴＲで繰り返し、１枚の画像を得るのに必要なエコーを計測する。各エコーは図４（ｂ）に示すようにｋ空間上に配置され、２次元逆フーリエ変換によって画像が再構成される。この撮影シーケンスによって得られる画像のピクセル値は、被検体パラメータＴ１（縦緩和時間）、Ｔ２（横緩和時間）、ρ（スピン密度）、Ｂ１（ＲＦ照射強度）に依存し、特に、Ｔ１を強調した画像が得られる特徴をもつ。

このRF-spoiled GEで変更可能な撮影パラメータは、ＦＡ(フリップ角)、ＴＲ(繰り返し時間)、ＴＥ(エコー時間)、θ(RF位相増分値)である。このうち、ＲＦ位相増分値は、高速撮影法の一つであるFLASH（登録商標）と同等のＴ２依存性の少ない画像コントラストが得られるように一般には117度に固定されている。このθを変化させると、画像コントラストのＴ２依存性が大きく変化する。

以上の各パラメータを用いて、RF-spoiled GEの信号関数ｆｓは式（１）表すことができる。

ここで、θは、ＲＦパルスの位相増分値、Ｓｃは装置パラメータの受信コイル感度である。Ｂ１は撮影時にはＦＡの係数となるため、ＦＡとの積の形にしておく。ρとＳｃは信号強度に対して比例係数として作用するため、関数の外側に出しておく。

数値シミュレーションでは、ｆｓにおける被検体パラメータのＴ１、Ｔ２の任意の値に対して撮影パラメータＦＡ、ＴＲ、θを網羅的に変化させて信号を作成し、補間により信号関数を作成する。撮影対象のスピン密度ρとＢ１、Ｓｃは一定(例えば１)とする。

撮影パラメータと被検体パラメータは例えば以下のように変化させる。それぞれ、実際の撮影に用いる撮影パラメータの範囲と、被検体のＴ１、Ｔ２の範囲が含まれるようにする。
TR 4個 [ms]: 10, 20, 30, 40
FA 13個 [度]: 3, 5, 7, 10, 12, 15, 20, 25, 30, 35, 40, 45, 50
θ 21個 [度]: 160度から180度まで（１度きざみ）
T2 17個 [s]: 0.01, 0.02, 0.03, 0.04, 0.05, 0.07, 0.1, 0.14, 0.19, 0.27, 0.38, 0.53, 0.74, 1.0, 1.4, 2.0, 2.8
T1 15個 [s]: 0.05, 0.07, 0.1, 0.14, 0.19, 0.27, 0.38, 0.53, 0.74, 1.0, 1.5, 2.0, 2.8, 4.0, 5.6

以上の撮影パラメータと被検体パラメータのすべての組み合わせからなる102000個の撮影パラメータセット(310)を構成し、それぞれの信号値を計算機シミュレーションによって算出する。数値シミュレーションは、格子点上にスピンを配置した被検体モデルと撮影シーケンス、撮影パラメータ、装置パラメータを入力とし、磁気共鳴現象の基礎方程式であるBlochの式を解いて磁気共鳴信号を出力する。被検体モデルはスピンの空間分布(γ, M0, T1, T2, Sc)として与えられる。ここで、γは磁気回転比、M0は熱平衡磁化(スピン密度)、T1とT2はそれぞれ縦緩和時間と横緩和時間である。磁気共鳴信号を画像再構成することにより、与えられた条件での画像を得ることができる。

Blochの式は1階線形常微分方程式であり、次式（２）で表される。

ここで、(x, y, z)は３次元の直交座標系を表し、zは静磁場(強度がＢ０)の向きに等しい。また、(Mx, My, Mz)はスピン、Gx, Gy, Gzはそれぞれ添字方向の傾斜磁場強度、Ｂ０は静磁場強度、Ｈ１は高周波磁場強度、Δf0は回転座標系の周波数である。

数値シミュレーションによって得られた信号値Ｉから、補間により信号関数fsを構成する。補間には１次から３次程度の線形補間やスプライン補間を用いることが可能である。

上述のようにして作成した信号関数の強度の一部を図５に示す。ここでは、被検体パラメータをＴ１=900ms、Ｔ２=100msとし、位相増分値θ=176度の場合について、横軸と縦軸をそれぞれＦＡ、ＴＲとして表示した。なお、信号関数は、一度作成して保存しておけば計算画像撮影のたびに作成する必要はなく、繰り返し使用可能である。

［Ｓ３０２：撮影条件探索ステップ］
撮影条件探索部２４０は、誤差伝播の法則を用いて元画像の撮影に用いる複数の撮影条件(撮影パラメータセット)を探索する。

本実施形態では、一例として、ヒト頭部の主な組織(灰白質GM、白質WM、脂肪fat、脳脊髄液CSF)を対象としてＴ１、Ｔ２、Ｂ１、ａの各マップの推定誤差（元画像に含まれるノイズが推定結果の各画像に与えるノイズの大きさ）が最小となるよう誤差伝播の法則を用いて最適化するものとする。上述した４つの組織の被検体パラメータの値は、図６に示すとおりである。信号関数からこれらを未知数として算出する場合、未知数の数が４(T1、T2、B1、a)であることから、４以上の撮影パラメータセットを決定する。ここでは６組の撮影パラメータセットを求めるものとする。

各撮影パラメータ(ＦＡ、ＴＲ、θ)の値は、例えば、以下のとおりとする。
ＦＡ[度]: 10, 20, 30, 40
ＴＲ[ms]: 10, 20, 30, 40
θ: 160度から179度まで（１度きざみ）

これらのいずれかを組み合わせて１組とし、合計３２０組（４×４×２０）の組み合わせのうち、６組のＴＲの合計が所定の範囲となる組み合わせを探索範囲として、探索する。ＴＲの合計に制限を設けたのは、探索の範囲を広げすぎないためと、撮影時間を長時間化しないためである。具体的にはＴＲの合計を100 ms～120 ms程度とする。

撮影パラメータセットの探索は、上述したように基本的に誤差伝搬の法則に従って行うが、撮影パラメータセットの中に、得られる画像の画質の良くないものや安定しないものが含まれると定量値マップもその影響を受ける恐れがある。本実施形態では、画質が劣化する或いは不安定となる撮影パラメータを除く処理を行う。画質が劣化する或いは不安定となる大きな要因は、血液や脳脊髄液の流れの変動であり、流れの変動による画質の劣化はアーチファクトの発生や連続撮影した場合の画像間の輝度の変動として現れる。本実施形態では、撮影条件毎に、予め画像間の輝度変動をあらかじめ調べて、画像間に輝度変動を生じない撮影パラメータセットを決定する。輝度変動を用いることで、明確にアーチファクトとして明確には現れない画質の劣化或いは目視では確認できない画質の劣化も評価することができ、また画質の評価を容易に自動化できる。

撮影パラメータセット探索に先立って行う処理の流れを、図７に示すフローを参照して説明する。ここでは、３つの撮影パラメータ（ＦＡ、ＴＲ、θ）のうち、画質が安定するθをＦＡ、ＴＲごとに決定する。一般に、個々の撮影パラメータについては、ＴＲは短いほうが流れの影響を受けやすく、ＦＡは小さいほうが流れの影響を受けにくいことが知られており、またθが１８０に近いほど流れのある領域のコントラストが高くなるため、流れの影響を受けやすい。本実施形態では、これらの知見をもとに、まず、最も輝度変動の出にくいＦＡとＴＲを初期値として指定し、θについては１８０度に近い値に初期値に設定し、順に許容できるθを確定するように探索を行う(S701)。ここでは、初期値として、ＦＡを最小の10度、ＴＲとθをそれぞれ最大の40 msと179度とする。

この条件にて連続して被検者を撮影する(S702)。撮影回数は、輝度変動が判断できれば良いので最低２回、できれば３回以上であれば良い。そして、図８に示すように、連続撮影した画像で輝度変動がないか評価する(S703)。輝度変動の評価は、目視により行うことも可能であるが、自動的に評価を行う場合は、連続して撮影した画像のうち１回目の画像と２回目以降の画像の差分を計算し、その変動が数パーセント(例えば５％)以下の場合に実質的に輝度変動がないと判断してもよい。差分は、単純な画素値間の差分でもよいが、例えば二乗平均平方根誤差（RMSE: root mean squared error)を算出する。

図８は、ＦＡ：30度、ＴＲ：20 msのときに、θを174度、173度、172度にそれぞれ異ならせて、３回の連続撮影を行って得た再構成画像を示した図であり、θが174度、173度では、図中、白矢印で示すように、画像間で輝度が明るくなったり暗くなったりしている。この輝度変動の主な原因は、血液や脳脊髄液の流れの変動である。これに対し、θが172度のときには、このような変動がなく、流れの影響がでにくい撮影パラメータセットであることがわかる。

ステップS703において、輝度変動があると判断されれば、θを更新（S704）、同様に連続撮影と輝度変動の評価を繰り返す。ここではθの初期値を最大（179度）に設定しているので、更新毎にθを１度小さくする。これを輝度変動がなくなるまで繰り返し、輝度変動がなくなったθをＦＡ、ＴＲとともに保存する(S705)。

次に、ＦＡを20度に更新し(S706,S707)、上述と同様に連続撮影(S702)、輝度変動判定(S703)を繰り返し、輝度変動がなくなったθをＦＡ、ＴＲとともに保存する(S705)。ＦＡを更新する際、θは初期値（10度）に戻さずＦＡ更新直前のままでよい。これは、例えばＦＡが10度と20度では、20度の方が輝度変動が出やすいため、ＦＡを20度にしたときに輝度変動のないθは、ＦＡが10度の場合と同じかそれより小さくなるからである。

こうして全てのＦＡについて評価が完了すれば(S706)、ＴＲを一つ小さい値（30 ms）に更新し、ＦＡを最小の10度に戻す(S709)。また、θを、更新前のＴＲと最小のＦＡ 10度において輝度変動がないと判定されたθの値に更新する。これは、ＴＲが小さいほど輝度変動が発生しやすいため、ＴＲを一つ小さく更新した場合の輝度変動のないθは、ＴＲ更新前のθと同じか小さくなるからである。そして、上述と同様に輝度変動の評価を繰り返す。

以上を繰り返して全てのＴＲについて完了すれば処理を終了する(S708)。これにより、誤差伝播の法則を用いて撮影パラメータセットの組み合わせを最適化する際に、輝度変動が発生する撮影パラメータを除いておくことができる。

以上の処理によって、輝度変動がないと判定した撮影パラメータセットの組み合わせ例を図９に示す。図９において（ａ）～（ｄ）は、それぞれ、ＴＲを10 ms、20 ms、30 ms、40 msとしたときに、ステップS705で選択（保存）された撮影パラメータセットを示す図であり、「〇」が保存された撮影パラメータセットを示し、「〇」よりもθの値が小さいものは、全て輝度変動がない撮影パラメータセットとなる。すなわち上記処理では、所定のＦＡ、ＴＲで輝度変動がないことが判定された最大のθが求められれば、同じＦＡ、ＴＲでは輝度変動の評価を省略することができる。

次いで、合計ＴＲの制限(ここでは120 ms以下)を満たすＦＡ、ＴＲ、θの全ての組合せの中から、上述のように、予め輝度変動のある条件を除いた後（図９の「○」をつけたθより大きいθを除いた後）、残る撮影パラメータセットについて、予め求めた信号関数３２を用い、誤差伝搬法則を従って計算画像取得用の撮影パラメータセット３３を決定する。

撮影パラメータセットは、上述したように、未知数の数と同数かそれ以上である必要がある。一例として、被検体パラメータのＴ１、Ｔ２、Ｂ１、ρと装置パラメータＳｃの積であるａ（＝ρＳｃ）の４つのパラメータを推定するものとし、次のような６組の撮影パラメータセット（ＦＡ、ＴＲ、θ）が決定される。
(40, 40, 173)、(10, 30, 177)、(40, 20, 160)、(40, 10, 162)、(40, 10, 162)、(40, 10, 162)

上記６組の撮影パラメータセットのうち、４つはそれぞれパラメータの値の組み合わせが異なるが、残りの２つは、４つのうちの一つと同一である。これら同じ撮影パラメータセットで取得した画像は加算して用いられる。すなわち、撮影パラメータのうちＳＮが低いことが予想される撮影パラメータについては、異なる撮影パラメータセットが未知数以上という条件を満たせば、重複させてＳＮを改善することが可能である。

［Ｓ３０３：エコー計測］
計測部１１０は、予め定めたパルスシーケンスと撮影条件決定部２２４が決定した撮影パラメータセット３３に従って、撮影シーケンスをシーケンサ１０４に設定し、エコー信号の計測を行う。そして、計測されたエコー信号をｋ空間に配置し計測データ（ｋ空間データ）３４を得る。

［Ｓ３０４：画像再構成］
画像再構成部２２１は、計測部１１０が収集した計測データに対し、逆フーリエ変換を行い、画像を再構成する。これにより６組の撮影パラメータセットで撮影した６個の再構成画像（元画像）３５が得られる。

［Ｓ３０５、Ｓ３０６：パラメータ推定・マップ生成］
パレメータ推定部２２３は、６個の元画像３５と信号関数３２を用いて被検体パラメータと装置パラメータを推定する。具体的には、ピクセルごとの信号値Iを、式（１）を変形した式（３）の関数ｆに対してフィッティングすることにより、４つのパラメータＴ１、Ｔ２、Ｂ１、ａ（ρ×Ｓｃ）の推定を行う。

関数フィッティングは、式（４）で表される最小二乗法により行う。

ここで、χは信号関数とファントムのピクセル値の残差の総和、Ｉは上述した撮影パラメータ(ＦＡ、ＴＲ、θ)で撮影して得られた再構成画像（元画像）におけるピクセル値である。この計算を、各ピクセルについて行うことで、再構成画像と同じマトリクスサイズの計算画像（パラメータマップ）を生成することができる。計算画像は、Ｔ１、Ｔ２、Ｂ１及びρＳｃの４種が生成可能である。また計算画像間の計算によって、新たな計算画像を生成してもよい。

［Ｓ３０７：画像表示］
計算画像は、必要に応じて、ディスプレイ１２０に表示される。

こうして推定した健常人頭部のＴ１、Ｔ２マップの一例を図１０（Ａ）、（Ｂ）に示す。また、比較のために画質の良くない撮影条件を除外せずに探索した撮影パラメータセットを用いて撮影し、生成したＴ１、Ｔ２マップを図１１（Ａ）、（Ｂ）に示す。図１０及び図１１において、（Ａ）はＴ１マップ、（Ｂ）はＴ２マップであり、各マップの下段は上段の画像に示す直線(前頭部)におけるプロファイルである。

図１１（Ａ）に示すＴ１マップでは、前頭白質部分の値が左右で異なっているのに対して、図１０（Ａ）のＴ１マップでは左右が等しくなっており、輝度が不安定な元画像が含まれていないことを示している。また、Ｔ２マップを比較すると、図１０（Ｂ）の方がマップ中心付近のノイズが少なくなっており、プロファイルのギザギザも小さくなっている。これも、画質の良くない元画像を除外したことの効果である。

＜変形例＞
本実施形態では、画質の良くない画像として、脳画像を例にして、脳脊髄液ＣＳＦによる画質の劣化がある場合を説明したが、肩関節や膝関節などの関節部位における血流の影響や、呼吸動や拍動の影響を受けやすい部位（胸部など）についても、撮影条件によって、フローや体動の画質への影響は異なる。これら部位についても、上述した手法と同様の手法（例えば図７の手順）で、撮影パラメータセットから、画質が良くないものを除去することができ、同様の効果を得ることができる。

また本実施形態では、連続撮影を行い、得られる複数の画像間の輝度の変動から除外すべき撮影パラメータセットを自動的に判断する場合を説明したが、部位やフロー或いは体動によって、画像上に観察できるアーチファクトが生じる場合もあり、そのような場合には画像上のアーチファクトの有無を判断してもよい。アーチファクトの有無の判断は、ユーザの目視で行ってもよいし、フローによるアーチファクトがでないことがわかっている撮影条件で撮影した画像と比較することにより行ってもよい。

また本実施形態では、撮影パラメータセットとして、ＦＡ、ＴＲ、θの３つのパラメータの組み合わせを説明したが、さらに別の撮影パラメータ（例えばＴＥ）を追加したり、上記３つのパラメータのうち一つ（例えばＴＲ）を固定して、他のパラメータのみを変化させて除外すべき撮影パラメータセットを探索するなど種々の変更が可能である。

さらに上記実施形態では、信号関数は、GEタイプのRF-spoiled GEシーケンスの信号関数を数値シミュレーションによって求める場合を説明したが、撮影シーケンスが数値計算によって信号関数が求められるものの場合には、図３のステップＳ３０１は省略することができる。

以上、説明したように、本実施形態及びその変形例によれば、パラメータ推定のために複数の撮影を行って複数の元画像を取得する際に、複数の撮影パラメータを適切に設定することができ、これにより、フローの影響などにより画質が良くない画像が元画像に含まれることを防ぎ、安定して信頼性のあるパラメータ推定及びマップ生成を行うことができる。また本実施形態によれば、撮影パラメータセットを誤差伝搬法則に従って最適化する際に、探索する撮影パラメータセットの範囲を限定することで、最適化に係る演算量を軽減することができる。

１０１：静磁場を発生するマグネット、１０２：傾斜磁場コイル、１０３：被検体、１０４：シーケンサ、１０５：傾斜磁場電源、１０６：高周波磁場発生器、１０７：プローブ、１０８：受信器、１１０：計算機、１２０：ディスプレイ、１３０：入力デバイス、１４０：記憶媒体、２１０：制御部、２２０：演算部、２２１：画像再構成部、２２２：信号関数算出部、２２３：パラメータ推定部、２２４：撮影条件決定部、２２５：撮影条件探索部。

Claims (11)

1. 予め定めた撮影条件および予め定めた撮影シーケンスに従って、静磁場の中に置かれた被検体に高周波磁場および傾斜磁場を印加して、前記被検体から発生するエコー信号を計測する計測部と、
   計測した前記エコー信号から再構成画像を作成する画像再構成部と、
   前記撮影条件を変えて前記撮影シーケンスを実行してそれぞれ得られた再構成画像と前記撮影シーケンスによって定まる信号関数とを用いて前記被検体に依存する被検体パラメータあるいは装置に依存する装置パラメータの分布を推定するパラメータ推定部と、
   前記パラメータ推定部が用いる複数の再構成画像を取得するための複数の撮影条件の組み合わせを決定する撮影条件決定部と、を備え、
   前記撮影条件決定部は、再構成画像の画質を評価し、前記複数の撮影条件の組み合わせから除外すべき撮影条件の組み合わせを探索する撮影条件探索部を有し、当該除外すべき撮影条件の組み合わせを除外して、前記複数の撮影条件の組み合わせを決定することを特徴とする磁気共鳴イメージング装置。
2. 請求項１に記載の磁気共鳴イメージング装置であって、
   前記撮影条件決定部は、再構成画像に発生するアーチファクトの大きさに基づき、パラメータ推定に必要な前記撮影条件を決定することを特徴とする磁気共鳴イメージング装置。
3. 請求項１に記載の磁気共鳴イメージング装置であって、
   前記撮影条件決定部は、同じ撮影条件で複数回撮影して得られる再構成画像間の輝度変動を評価し、その結果に基づき、パラメータ推定に必要な撮影条件を決定することを特徴とする磁気共鳴イメージング装置。
4. 請求項３に記載の磁気共鳴イメージング装置であって、
   前記撮影条件決定部は、再構成画像間の差分または二乗平均平方根誤差をもとに輝度変動を評価することを特徴とする磁気共鳴イメージング装置。
5. 請求項１に記載の磁気共鳴イメージング装置であって、
   前記撮影条件決定部が決定する前記撮影条件は、前記撮影シーケンスの繰り返し時間（ＴＲ）、前記高周波磁場のフリップ角（ＦＡ）、及び前記高周波磁場の位相増分値（θ）の組合せであることを特徴とする磁気共鳴イメージング装置。
6. 請求項５に記載の磁気共鳴イメージング装置であって、
   前記撮影条件決定部は、所定の繰り返し時間（ＴＲ）及び所定のフリップ角（ＦＡ）において、同じ撮影条件で複数回撮影して得られる再構成画像間の輝度変動をもとに当該繰り返し時間及びフリップ角における高周波磁場の位相増分値（θ）を決定することを特徴とする磁気共鳴イメージング装置。
7. 請求項５に記載の磁気共鳴イメージング装置であって、
   前記撮影条件決定部は、前記撮影条件の組み合わせとして、前記パラメータ推定部が推定するパラメータの数と同数又はそれ以上の組み合わせを決定することを特徴とする磁気共鳴イメージング装置。
8. 請求項７に記載の磁気共鳴イメージング装置であって、
   前記撮影条件決定部が決定する組み合わせは、同一の組み合わせを含むことを特徴とする磁気共鳴イメージング装置。
9. 磁気共鳴イメージング装置において異なる撮影条件で取得した複数の再構成画像を処理し、被検体に依存する被検体パラメータあるいは装置に依存する装置パラメータの分布を推定する画像処理方法であって、
   撮影条件毎に再構成画像の画質を評価し、前記磁気共鳴イメージング装置による前記分布の推定に用いる撮影条件を探索するステップを含むことを特徴とする画像処理方法。
10. 請求項９に記載の画像処理方法であって、
    前記撮影条件を探索するステップは、再構成画像における流れ又は体動によるアーチファクトを評価し、再構成画像におけるアーチファクトが少ない撮影条件を、前記分布の推定に用いる撮影条件として決定することを特徴とする画像処理方法。
11. 請求項１０に記載の画像処理方法であって、
    前記撮影条件を探索するステップは、同じ撮影条件で複数回撮影して得られる再構成画像間の輝度変動から前記アーチファクトを評価することを特徴とする画像処理方法。

Images