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

Description

本発明の実施形態は、画像処理装置及び磁気共鳴イメージング装置に関する。

磁気共鳴イメージング装置（以下、適宜「ＭＲＩ(Magnetic Resonance Imaging)装置」）は、核磁気共鳴現象を利用して、被検体内部の情報を画像化するための装置である。ＭＲＩ装置は、コイルを利用して物体内部にある特定の原子（例えば、水素原子）からの核磁気共鳴信号をサンプリングすることで、ｋ空間データと呼ばれるデータを取得する。また、ＭＲＩ装置は、ｋ空間データにフーリエ変換を適用することで、磁気共鳴画像（以下、適宜「ＭＲ（Magnetic Resonance）画像）」）を再構成する。

核磁気共鳴信号（以下、適宜「ＭＲ（Magnetic Resonance）信号」）は、１次元のデータとしてサンプリングされる。そこで、ＭＲＩ装置は、２次元あるいは３次元の再構成画像を得るために、ｋ空間上で１次元のサンプリングを繰り返し行い、再構成に必要なｋ空間データを取得する。ＭＲ画像と同じ解像度（フルサンプリング）でｋ空間データをサンプリングすれば、得られたｋ空間データにフーリエ変換を適用することで再構成が可能である。しかし、ＭＲＩ装置でフルサンプリングを行うと、その撮像時間は非常に長くなる。そこで、従来より、撮像時間短縮に向けて、撮像手法及び再構成手法の両面において技術開発がなされている。

例えば、撮像時間短縮のための再構成手法として、複数のコイルを利用し、フルサンプリングよりも少ないサンプル数で被検体の撮像を行い、コイル毎のコイル感度の違いを利用してＭＲ画像を再構成する、パラレルイメージング（以下、適宜「ＰＩ（Parallel Imaging）」）と呼ばれる技術がある。ＰＩの技術として、例えば、ＳＥＮＳＥ（sensitivity encoding）及びＧＲＡＰＰＡ（Generalized Auto calibrating Partially Parallel Acquisition）が広く知られている。

米国特許第８０２６７２０号明細書 米国特許第７２８２９１７号明細書 米国特許第７６８８６８号明細書 K.P.Pruessmann et al., "SENSE: Sensitivity Encoding for Fast MRI," MRM 42: 952-962, 1999. Mark A. Grisworld et al., "Generalized Autocalibrating Partially Parallel Acquisitions（GRAPPA）," MRM, vol.47，pp.1202-1210, 2002. A. C. S. Brau et al., Comparison of Reconstruction Accuracy and Efficiency Among Autocalibrating Data-driven Parallel Imaging Methods, MRM, vol.59, pp.382-395, 2008.

本発明が解決しようとする課題は、画質を向上することができる磁気共鳴イメージング装置及び画像処理装置を提供することである。

実施形態に係る磁気共鳴イメージング装置は、ｋ空間データ生成部と、画像生成部と、導出部と、算出部とを備える。ｋ空間データ生成部は、ｋ空間上のフルサンプリング位置の一部にサンプル値を有する第１のｋ空間データのうち、前記サンプル値を有しないサンプリング位置の少なくとも一部にサンプル値を与えて第２のｋ空間データを生成する。画像生成部は、前記第１のｋ空間データから第１画像を生成し、前記第２のｋ空間データから第２画像を生成する。導出部は、前記第１画像及び前記第２画像に対する重み係数を導出する。算出部は、前記第１画像及び前記第２画像に対して、前記重み係数を用いた重み付け加算を行い、磁気共鳴画像を算出する。

図１は、第１の実施形態に係るＭＲＩ装置１００を示す機能ブロック図。 図２は、ＧＲＡＰＰＡにおけるサンプル値の関係を説明するための図。 図３は、第１の実施形態におけるＭＲ画像の重み付け線形和を説明するための図。 図４は、第１の実施形態における処理手順の概要を示すフローチャート。 図５は、第１の実施形態における出力ＭＲ画像を得る処理を説明するための図。 図６は、第１の実施形態におけるｋ空間データを説明するための図。 図７は、第１の実施形態における処理手順を示すフローチャート。 図８は、第１の実施形態における処理手順を示すフローチャート。 図９は、第１の実施形態における処理手順を示すフローチャート。 図１０は、第１の実施形態の変形例を説明するための図。 図１１は、その他の実施形態におけるｋ空間データを説明するための図。 図１２は、その他の実施形態におけるｋ空間データを説明するための図。 図１３は、その他の実施形態におけるｋ空間データを説明するための図。 図１４は、その他の実施形態におけるｋ空間データを説明するための図。

各実施形態に係るＭＲＩ装置及び画像処理装置は、ＳＥＮＳＥあるいはＧＲＡＰＰＡと同様、フルサンプリングよりも少ないサンプル数で被検体の撮像を行い、コイル毎のコイル感度の違いを利用して、ＭＲ画像を再構成する。既存のＰＩ技術の大半は、ＳＥＮＳＥあるいはＧＲＡＰＰＡの改良技術であるが、以下では、既存の枠組みとは異なる定式化によって、ＰＩ技術を提供する。以下、図面を参照しながら、各実施形態に係るＭＲＩ装置及び画像処理装置を説明する。なお、実施形態は、以下の実施形態に限られるものではない。また、各実施形態において説明する内容は、原則として、他の実施形態においても同様に適用することができる。

（第１の実施形態）
図１は、第１の実施形態に係るＭＲＩ装置１００を示す機能ブロック図である。図１に示すように、ＭＲＩ装置１００は、静磁場磁石１０１と、静磁場電源１０２と、傾斜磁場コイル１０３と、傾斜磁場電源１０４と、寝台１０５と、寝台制御部１０６と、送信コイル１０７と、送信部１０８と、受信コイルアレイ１０９と、受信部１１０と、シーケンス制御部１２０と、計算機１３０とを備える。なお、ＭＲＩ装置１００に被検体Ｐ（例えば、人体）は含まれない。

静磁場磁石１０１は、中空の円筒形状に形成された磁石であり、内部の空間に静磁場を発生する。静磁場磁石１０１は、例えば、超伝導磁石等であり、静磁場電源１０２から電流の供給を受けて励磁する。静磁場電源１０２は、静磁場磁石１０１に電流を供給する。なお、静磁場磁石１０１は、永久磁石でもよく、この場合、ＭＲＩ装置１００は、静磁場電源１０２を備えなくてもよい。また、静磁場電源１０２は、ＭＲＩ装置１００とは別に備えられてもよい。

傾斜磁場コイル１０３は、中空の円筒形状に形成されたコイルであり、静磁場磁石１０１の内側に配置される。傾斜磁場コイル１０３は、互いに直交するＸ，Ｙ，Ｚの各軸に対応する３つのコイルが組み合わされて形成されており、これら３つのコイルは、傾斜磁場電源１０４から個別に電流を受けて、Ｘ，Ｙ，Ｚの各軸に沿って磁場強度が変化する傾斜磁場を発生させる。なお、Ｚ軸方向は、静磁場と同方向とする。

傾斜磁場電源１０４は、傾斜磁場コイル１０３に電流を供給する。ここで、傾斜磁場コイル１０３によって発生するＸ，Ｙ，Ｚ各軸の傾斜磁場は、例えば、スライス選択用傾斜磁場Ｇｓ、位相エンコード用傾斜磁場Ｇｅ、及びリードアウト用傾斜磁場Ｇｒにそれぞれ対応する。スライス選択用傾斜磁場Ｇｓは、任意に撮像断面を決めるために利用される。位相エンコード用傾斜磁場Ｇｅは、空間的位置に応じてＭＲ信号の位相を変化させるために利用される。リードアウト用傾斜磁場Ｇｒは、空間的位置に応じてＭＲ信号の周波数を変化させるために利用される。

寝台１０５は、被検体Ｐが載置される天板１０５ａを備え、寝台制御部１０６による制御の下、天板１０５ａを、被検体Ｐが載置された状態で傾斜磁場コイル１０３の空洞（撮像口）内へ挿入する。通常、寝台１０５は、長手方向が静磁場磁石１０１の中心軸と平行になるように設置される。寝台制御部１０６は、計算機１３０による制御の下、寝台１０５を駆動して天板１０５ａを長手方向及び上下方向へ移動する。

送信コイル１０７は、傾斜磁場コイル１０３の内側に配置され、送信部１０８からＲＦ（Radio Frequency）パルスの供給を受けて、高周波磁場を発生する。送信部１０８は、対象とする原子の種類及び磁場の強度で決まるラーモア周波数に対応するＲＦパルスを送信コイル１０７に供給する。

受信コイルアレイ１０９は、傾斜磁場コイル１０３の内側に配置され、高周波磁場の影響によって被検体Ｐから発せられるＭＲ信号を受信する。受信コイルアレイ１０９は、ＭＲ信号を受信すると、受信したＭＲ信号を受信部１１０へ出力する。なお、第１の実施形態において、受信コイルアレイ１０９は、複数の受信コイルを有するコイルアレイである。

受信部１１０は、受信コイルアレイ１０９から出力されるＭＲ信号に基づいてＭＲデータを生成する。具体的には、受信部１１０は、受信コイルアレイ１０９から出力されるＭＲ信号をデジタル変換することによってＭＲデータを生成する。また、受信部１１０は、生成したＭＲデータをシーケンス制御部１２０へ送信する。なお、受信部１１０は、静磁場磁石１０１や傾斜磁場コイル１０３等を備える架台装置側に備えられていてもよい。

ここで、第１の実施形態において、受信コイルアレイ１０９の各コイルエレメント（各受信コイル）から出力されるＭＲ信号は、適宜分配合成されることで、チャネル等と呼ばれる単位で受信部１１０に出力される。このため、ＭＲデータは、受信部１１０以降の後段の処理においてチャネル毎に取り扱われる。コイルエレメントの総数とチャネルの総数との関係は、同一の場合もあれば、コイルエレメントの総数に対してチャネルの総数が少ない場合、あるいは反対に、コイルエレメントの総数に対してチャネルの総数が多い場合もある。以下において、「チャネル」と表記する場合、その処理が、コイルエレメント単位で行われてもよいし、あるいは、コイルエレメントが分配合成されたチャネル単位で行われてもよいことを示す。なお、分配合成のタイミングは、上述したタイミングに限られるものではない。ＭＲ信号若しくはＭＲデータは、後述する画像の再構成処理の前までに、チャネル単位に分配合成されればよい。

シーケンス制御部１２０は、計算機１３０から送信されるシーケンス情報に基づいて、傾斜磁場電源１０４、送信部１０８及び受信部１１０を駆動することによって、被検体Ｐの撮像を行う。ここで、シーケンス情報は、撮像を行うための手順を定義した情報である。シーケンス情報には、傾斜磁場電源１０４が傾斜磁場コイル１０３に供給する電源の強さや電源を供給するタイミング、送信部１０８が送信コイル１０７に送信するＲＦパルスの強さやＲＦパルスを印加するタイミング、受信部１１０がＭＲ信号を検出するタイミング等が定義される。

なお、シーケンス制御部１２０は、傾斜磁場電源１０４、送信部１０８及び受信部１１０を駆動して被検体Ｐを撮像した結果、受信部１１０からＭＲデータを受信すると、受信したＭＲデータを計算機１３０へ転送する。

計算機１３０は、ＭＲＩ装置１００の全体制御や、データ収集、画像再構成等を行い、インタフェース部１３１、記憶部１３２、制御部１３３、入力部１３４、及び表示部１３５を有する。

インタフェース部１３１は、シーケンス情報をシーケンス制御部１２０へ送信し、シーケンス制御部１２０からＭＲデータを受信する。また、インタフェース部１３１は、ＭＲデータを受信すると、受信したＭＲデータを記憶部１３２に格納する。記憶部１３２に格納されたＭＲデータは、後述するｋ空間データ収集部１３３ａによってｋ空間に配置される。この結果、記憶部１３２は、複数チャネル分のｋ空間データを記憶する。

記憶部１３２は、インタフェース部１３１によって受信されたＭＲデータや、ｋ空間データ収集部１３３ａによってｋ空間に配置されたｋ空間データ、制御部１３３によって生成された画像データ等を記憶する。例えば、記憶部１３２は、ＲＡＭ（Random Access Memory）、フラッシュメモリ等の半導体メモリ素子、ハードディスク、光ディスク等である。

入力部１３４は、操作者からの各種指示や情報入力を受け付ける。入力部１３４は、例えば、マウスやトラックボール等のポインティングデバイス、モード切替スイッチ等の選択デバイス、あるいはキーボード等の入力デバイスである。表示部１３５は、制御部１３３による制御の下、スペクトラムデータや画像データ等の各種の情報を表示する。表示部１３５は、例えば、液晶ディスプレイ等の表示デバイスである。

制御部１３３は、ＭＲＩ装置１００の全体制御を行う。具体的には、制御部１３３は、入力部１３４を介して操作者から入力される撮像条件に基づいてシーケンス情報を生成し、生成したシーケンス情報をシーケンス制御部１２０に送信することによって撮像を制御する。また、制御部１３６は、撮像の結果としてシーケンス制御部１２０から送られるＭＲデータに基づいて行われる画像の再構成を制御したり、表示部１３５による表示を制御したりする。例えば、制御部１３３は、記憶部１３２によって記憶されたｋ空間データに対して、フーリエ変換処理等の再構成処理を施すことによって、スペクトラムデータや、画像データを生成する。例えば、制御部１３３は、ＡＳＩＣ（Application Specific Integrated Circuit）、ＦＰＧＡ（Field Programmable Gate Array）等の集積回路、ＣＰＵ（Central Processing Unit）、ＭＰＵ（Micro Processing Unit）等の電子回路である。

ここで、ＰＩの技術の１つであるＳＥＮＳＥ及びＧＲＡＰＰＡを説明する。ＳＥＮＳＥでは、コイル毎の感度を表す感度マップが予め収集される。また、ｋ空間データは、フルサンプリングに対して等間隔に間引かれたサンプリングにより収集される。間引かれたｋ空間データにフーリエ変換が適用され、再構成が行われると、コイル毎に、折り返し画像が得られる。なお、ＳＥＮＳＥにおいては、フルサンプリングによるサンプリング量を、間引きサンプリングによるサンプリング量で除した値を「Reduction Factor」と呼び、「Ｒ」で表す。

再構成された折り返し画像は、Ｒ個の折り返し信号が重ねられた画像である。このとき、コイルを「ｃ＝１，・・・，Ｃ」で表し、各コイルの折り返し画像の画素値を「ρ（ｃ）」で表すと、折り返し画像内で注目する１画素の画素値は、次の（１）式で表される。また、（１）式をコイルの本数分並べた線形システムは、行列表記を用いて、（２）式で表される。なお、「出力ＭＲ画像」は、出力対象のＭＲ画像のことである。

このとき、コイルのノイズ共分散行列Ψによって重み付けされた最小二乗法による解は、（３）式で表される。ＳＥＮＳＥでは、各画素について（３）式を適用することで、出力ＭＲ画像を再構成する。

一方、ＧＲＡＰＰＡでは、複数のコイルにより収集されたｋ空間データについて、未知のサンプル値が周辺のサンプル値の重み付け線形和に等しいと考える。また、ＧＲＡＰＰＡでは、ｋ空間の中心から離れた領域（周辺領域）については、ＳＥＮＳＥと同様、フルサンプリングに対して等間隔に間引かれたサンプリングが行われ、ｋ空間の中心に近い領域（中心領域）については、ＳＥＮＳＥと異なり、フルサンプリングでサンプリングが行われる。

ＧＲＡＰＰＡでは、「重み係数の算出」及び「未知のサンプル値の生成」の２つのステップが行われる。「重み係数の算出」のステップでは、フルサンプリングされた領域を利用して重み係数が推定される。具体的には、ＧＲＡＰＰＡでは、サンプル値を得ようとするｋ空間上のサンプリング位置を「ｋ」で表し、サンプル値を有するｋ空間上のサンプリング位置における「ｋ」からの移動量を「ｂ」で表すときに、（４）式が成り立つと考える。なお、ここでは、ｋ空間データが２次元の場合を説明するが、３次元の場合にも同様の考えが成り立つ。

「Ｆ（ｋ，ｃout）」は、コイルｃoutにおけるサンプリング位置ｋのサンプル値である。また、「ｗ（ｂ，ｃ，ｃout）」は、移動量ｂ、コイルｃに対応する、コイルｃoutにおける重み係数である。また、「Ｆ（ｋ＋ｂ，ｃ）」は、コイルｃにおいて、サンプリング位置「ｋ」から「ｂ」だけ移動したサンプリング位置でサンプリングされたサンプル値である。

図２は、ＧＲＡＰＰＡにおけるサンプル値の関係を説明するための図である。図２は、（４）式を具体的な数値例で説明するものであり、一例に過ぎない。図２に示すように、サンプル値Ｆ（ｋ，１）は、コイル１におけるサンプリング位置ｋのサンプル値である。また、サンプル値Ｆ（ｋ＋１，３）は、コイル３において、サンプリング位置ｋから「１」だけ移動したサンプリング位置のサンプル値である。図２に示すように、サンプル値Ｆ（ｋ，１）は、コイル１〜３それぞれにおいて、サンプリング位置ｋから「１」だけ移動したサンプリング位置のサンプル値、コイル１〜３それぞれにおいて、サンプリング位置ｋから「−１」だけ移動したサンプリング位置のサンプル値の重み付け線形和で表される。重み係数は、それぞれ、ｗ（１，１，１）、ｗ（１，２，１）、ｗ（１，３，１）、ｗ（−１，１，１）、ｗ（−１，２，１）、ｗ（−１，３，１）である。

フルサンプリングされた領域では全てのサンプル値（即ち、Ｆ（ｋ，ｃout）及びＦ（ｋ＋ｂ，ｃ））が既知であるから、例えば、最小二乗法により重み係数ｗ（ｂ，ｃ，ｃout）を推定することができる。そして、「未知のサンプル値の生成」のステップでは、「重み係数の算出」のステップで得られた重み係数を（４）式にあてはめて、未知のサンプル値を生成する。

ここで、第１の実施形態においては、出力ＭＲ画像を、この出力ＭＲ画像と同じ解像度を有する複数のＭＲ画像の重み付け加算で表現できると考える。また、第１の実施形態においては、重み係数の推定において、重み係数に制約を与える条件の一例として、（４）式と同等のデータ拘束を与える手法を説明する。なお、「（４）式と同等のデータ拘束」とは、フルサンプリングされた領域における、全てのサンプリング位置「ｋ」に対する、（４）式の左辺と右辺との差を２乗した値（誤差）の総和のことである。もっとも、仮に重み係数の推定において（４）式と同等のデータ拘束を用いたとしても、第１の実施形態において、出力ＭＲ画像の再構成処理は、ｋ空間におけるサンプル値の補間ではなく、ＭＲ画像間の重み付け線形和によって行われるので、ＧＲＡＰＰＡとは全く異なるものである。

出力ＭＲ画像は、予め定められた出力解像度で配置された格子点上で、画素値を持つ。もっとも、ＭＲＩ装置１００は、ｋ空間データのサンプリングを行う際に、出力解像度に対応した全ての格子点でサンプリングを行うとは限らない。以下では、説明の便宜上、出力解像度に対応した格子点の全てについてサンプリングを行うことを「フルサンプリング」と呼ぶ。また、出力解像度に対応した各格子点の全部あるいは一部でサンプリングを行う手法を「Ｃａｒｔｅｓｉａｎサンプリング」と呼び、格子点でない点も含めてサンプリングを行う手法を「非Ｃａｒｔｅｓｉａｎサンプリング」と呼ぶ。

Ｃａｒｔｅｓｉａｎサンプリング、非Ｃａｒｔｅｓｉａｎサンプリングにかかわらず、サンプリング間隔が大きい場合、即ち、ＭＲ画像の再構成に必要なサンプリング密度よりもサンプル値が少ない場合、再構成されたＭＲ画像は、エイリアシングにより複数の信号が重ねられた画素値を持つ。ＰＩ技術では、コイルの感度を利用して、重ねられた複数の信号から、エイリアシングされていない信号値を復元する。

第１の実施形態においては、線形の重みを考え、エイリアシングされていない（復元の目的である）出力ＭＲ画像は、エイリアシングされたＭＲ画像を重み付けした形で表せると考え、その重み係数を推定した上で、その重み係数を利用して、出力ＭＲ画像を復元する。

図３は、第１の実施形態におけるＭＲ画像の重み付け線形和を説明するための図である。ＰＩの考え方から、複数のコイル（チャネル）で収集されたデータに基づいてＭＲ画像を生成し、生成したＭＲ画像に適切な重み係数を乗じた上で合成すれば、エイリアシングのない再構成画像が得られるはずである。言い換えれば、図３に示すように、複数のコイル（チャネル）で収集されたデータから生成された複数のＭＲ画像の重み付け線形和によって表現可能な出力解像度の空間における部分空間上に、出力ＭＲ画像は存在すると考えることができる。また、線形和の重み係数を推定する方法の１例として、コイルの感度とリンクした制約条件を用いた重み係数の推定を利用することができる。なお、目的は、所望の再構成画像を得ることであって、線形和の重み係数を推定する方法は、必ずしもコイルの感度と密接にリンクしている必要はない。また、図３に示すように、出力ＭＲ画像の解像度を「２５６×２５６ピクセル」とした場合、その変数の数は「６５５３６」である。これに対し、仮に、チャネル数を「１６」、各チャネルから生成された補正用ＭＲ画像の数を「４」とした場合、その変数の数は、（折り返し画像の分を含めても）「６５」である。すべての画素を推定する手法と比べ、出力ＭＲ画像を複数のＭＲ画像の重み付け線形和として表現しその重みを推定する手法を用いた方が、推定すべき変数の数が少なくてすむことがわかる。

より具体的には、第１の実施形態に係るＭＲＩ装置１００は、フルサンプリングよりも少ないサンプル数でｋ空間データを収集すると、このｋ空間データに対して何等かの処理を施すことで別のｋ空間データを得る。そして、ＭＲＩ装置１００は、オリジナルのｋ空間データから得られた「折り返しのある」ＭＲ画像と、別のｋ空間データから得られたＭＲ画像との重み付け加算によって、「折り返しのない」出力ＭＲ画像を得る。

なお、以下では、重み付け加算の対象となるＭＲ画像を「補正用ＭＲ画像」と呼び、補正用ＭＲ画像の重み付け加算で得られるＭＲ画像を「出力ＭＲ画像」と呼ぶが、これらは単なる説明の便宜上の区別に過ぎない。例えば、補正用ＭＲ画像が出力されてもよい。また、例えば、出力ＭＲ画像に対して後処理が施されたＭＲ画像が、最終的に表示部１３５等に出力されるものであってもよい。

制御部１３３は、ｋ空間データ収集部１３３ａと、補正用ＭＲ画像生成部１３３ｂと、重み係数導出部１３３ｃと、出力ＭＲ画像算出部１３３ｄとを備える。

ｋ空間データ収集部１３３ａは、操作者から入力された撮像条件に従ってｋ空間データを収集する。具体的には、ｋ空間データ収集部１３３ａは、各種撮像条件の入力をＧＵＩ（Graphical User Interface）上で受け付け、受け付けた撮像条件に従ってシーケンス情報を生成し、生成したシーケンス情報をシーケンス制御部１２０へ送信する。また、ｋ空間データ収集部１３３ａは、シーケンス制御部１２０から受信したＭＲデータを、ｋ空間に配置する。

ここで、第１の実施形態に係るｋ空間データ収集部１３３ａは、フルサンプリングではなく、フルサンプリングよりも少ないサンプル数のｋ空間データを収集する。例えば、ｋ空間データ収集部１３３ａは、ｋ空間の中心領域においては、フルサンプリングと同一の解像度でサンプリングを行い、ｋ空間の周辺領域においては、間引きサンプリングを行う。間引きサンプリングにより得られたｋ空間データは、フルサンプリングに対応するｋ空間上のサンプリング位置のうち、一部のサンプリング位置にサンプル値を有するものとなる。なお、以下では、ｋ空間データ収集部１３３ａによって収集されたオリジナルのｋ空間データを、適宜「撮像ｋ空間データ」と呼び、オリジナルのｋ空間データに対して何等かの処理を施すことで得られた別のｋ空間データと区別する。

補正用ＭＲ画像生成部１３３ｂは、出力ＭＲ画像を算出するための補正用ＭＲ画像を生成する。具体的には、補正用ＭＲ画像生成部１３３ｂは、フルサンプリングよりも少ないサンプル数の撮像ｋ空間データに対して、サンプル値を有しないサンプリング位置の少なくとも一部にサンプル値を与える処理を施すことで、別のｋ空間データを生成する。また、補正用ＭＲ画像生成部１３３ｂは、撮像ｋ空間データから、補正用ＭＲ画像を生成し、撮像ｋ空間データから得られた別のｋ空間データからも、補正用ＭＲ画像を生成する。補正用ＭＲ画像生成部１３３ｂによる処理の詳細は、後述する。

重み係数導出部１３３ｃは、補正用ＭＲ画像生成部１３３ｂによって生成された補正用ＭＲ画像に関する重み係数を導出する。重み係数導出部１３３ｃによる処理の詳細は、後述する。

出力ＭＲ画像算出部１３３ｄは、補正用ＭＲ画像生成部１３３ｂによって生成された補正用ＭＲ画像に対して、重み係数導出部１３３ｃによって導出された重み係数を用いた重み付け加算を行い、出力ＭＲ画像を算出する。出力ＭＲ画像算出部１３３ｄによる処理の詳細は、後述する。

さて、具体的な処理の説明に入る前に、まず、第１の実施形態における処理の定式化を行う。以下では、出力ＭＲ画像を縦ベクトルＸで表し、Ｘと同じ大きさを持つＮ個（Ｎ≧２）の補正用ＭＲ画像を、縦ベクトルＸiで表し、重み係数（スカラー値）をｗiで表す。補正用ＭＲ画像Ｘiは、出力ＭＲ画像Ｘと同一画素数を持つ画像である。以下の実施形態において、出力ＭＲ画像Ｘは、（５）式で与えられる。補正用ＭＲ画像Ｘiが与えられると、出力ＭＲ画像Ｘの自由度は、未知の重み係数の数Ｎを超えることはない。なお、任意のｍ、ｎ（ｍ≠ｎ）に対し、ＸmとＸnとは、独立であっても、従属であってもよい。

このように、出力ＭＲ画像Ｘは、補正用ＭＲ画像Ｘiの線形和として表現可能な空間上の点であり、補正用ＭＲ画像Ｘiに対する重みは、ｗiで表現される。言い換えると、出力ＭＲ画像Ｘが所望の再構成画像であるためには、所望の再構成画像が補正用ＭＲ画像Ｘiの線形和として表現されなければならない。

そこで、以下の実施形態においては、出力ＭＲ画像Ｘを得る問題を、「補正用ＭＲ画像Ｘiを与える過程」と、「重み係数ｗiを導出する過程」とを含む問題として捉える。具体的には、補正用ＭＲ画像生成部１３３ｂが、補正用ＭＲ画像Ｘiを生成し、重み係数導出部１３３ｃが、重み係数ｗiを導出する。そして、出力ＭＲ画像算出部１３３ｄが、補正用ＭＲ画像Ｘiと、重み係数ｗiとを用いて、（５）式の関係から、出力ＭＲ画像Ｘを算出する。

なお、出力チャネルが複数の場合、チャネルｋ＝１，・・・，Ｃに対する出力ＭＲ画像Ｘkは（６）式で与えられる。この場合、補正用ＭＲ画像生成部１３３ｂが、補正用ＭＲ画像Ｘi,kを生成し、重み係数導出部１３３ｃが、重み係数ｗi,kを導出する。そして、出力ＭＲ画像算出部１３３ｄが、補正用ＭＲ画像Ｘi,kと、重み係数ｗi,kとを用いて、（６）式の関係から、出力ＭＲ画像Ｘkを算出する。出力チャネルが複数の場合には、各チャネルｋ＝１，・・・，Ｃに対して、同一の処理が適用される。以下では、説明の便宜上、チャネルｋを省略して表現するが、出力チャネルが複数の場合にも、以下の実施形態を同様に適用することができる。

続いて、第１の実施形態における処理手順の概要を説明する。図４は、第１の実施形態における処理手順の概要を示すフローチャートであり、図５は、第１の実施形態における出力ＭＲ画像を得る処理を説明するための図である。なお、図４に示す処理手順の前には、ｋ空間データ収集部１３３ａによるｋ空間データの収集が既に行われ、処理対象のｋ空間データが、記憶部１３２に記憶されているものとする。

ここで、以下で用いる用語を定義する。まず、出力解像度に対応した全ての格子点についてサンプリングを行う操作を「フルサンプリング」と呼び、出力解像度に対応した格子点のうちの一部についてサンプリングを行う操作を「間引きサンプリング」と呼ぶ。また、各格子点を「サンプリング位置」と呼び、サンプリングが行われた「サンプリング位置」のデータを「既知サンプル」（図５において黒丸）と呼び、間引きサンプリングによってサンプリングが行われなかった「サンプリング位置」のデータを「未知サンプル」（図５において白丸）と呼ぶ。即ち、「既知サンプル」は、サンプル値を有し、「未知サンプル」は、サンプル値を有しない。

ステップＳ１：補正用ＭＲ画像生成部１３３ｂが、記憶部１３２から、処理対象の撮像ｋ空間データを読み出す。第１の実施形態において、この撮像ｋ空間データは、図５に示すように、ｋ空間の中心領域においては、フルサンプリングと同一の解像度でサンプリングが行われ、ｋ空間の周辺領域においては、間引きサンプリングが行われたものである。既知サンプルは複数のチャネルのそれぞれについてサンプリングされたデータであってもよい。

ステップＳ２：補正用ＭＲ画像生成部１３３ｂは、記憶部１３２から読み出した撮像ｋ空間データを用いて、補正用ＭＲ画像を生成する。例えば、補正用ＭＲ画像生成部１３３ｂは、撮像ｋ空間データのうち、既知サンプルが存在しないサンプリング位置のサンプル値を「ゼロ」で埋めることで、フルサンプリングと同一の解像度のｋ空間データを生成する。そして、補正用ＭＲ画像生成部１３３ｂは、生成したｋ空間データを再構成することで、第１補正用ＭＲ画像を生成する。この第１補正用ＭＲ画像は、撮像ｋ空間データの既知サンプルを加工せずに用いて生成された画像である。

また、補正用ＭＲ画像生成部１３３ｂは、撮像ｋ空間データを用いて、別途定めた規則に従った既知サンプルの平行移動（シフト）を行い、フルサンプリングと同一の解像度の新たなｋ空間データ（以下、適宜「シフトｋ空間データ」）を生成する。そして、補正用ＭＲ画像生成部１３３ｂは、生成したシフトｋ空間データを再構成することで、第２補正用ＭＲ画像を生成する。この第２補正用ＭＲ画像は、撮像ｋ空間データの既知サンプルを加工して生成された画像である。既知サンプルが複数のチャネルに対するデータであれば、補正用ＭＲ画像生成部１３３ｂは、第２補正用ＭＲ画像を全チャネル分だけ生成する。

例えば、補正用ＭＲ画像生成部１３３ｂは、別途定めた平行移動の規則に従って、既知サンプルをｋ空間上で移動させる。次に、補正用ＭＲ画像生成部１３３ｂは、「平行移動されたサンプル値が存在し、且つ、既知サンプルが存在しない」という条件を満たすサンプリング位置のサンプル値を「平行移動されたサンプル値」とし、そうでないサンプリング位置のサンプル値を「ゼロ」として、フルサンプリングと同じ解像度にしたシフトｋ空間データを生成する。そして、補正用ＭＲ画像生成部１３３ｂは、そのｋ空間データを再構成した第２補正用ＭＲ画像を生成する。

ステップＳ３：重み係数導出部１３３ｃは、重み係数を導出する。例えば、重み係数導出部１３３ｃは、撮像ｋ空間データのうちフルサンプリングされた領域を利用し、重み係数を推定する。例えば、重み係数導出部１３３ｃは、「ｋ空間のあるサンプリング位置におけるサンプル値は周辺のサンプル値の重み付け線形和で表現される」という仮定の下、フルサンプリングされた領域で、（４）式と同等のデータ拘束を用いて、あるサンプリング位置におけるサンプル値を、周辺のサンプル値の重み付け線形和で表現する場合の重み係数を推定する。これを言い換えると、重み係数導出部１３３ｃは、注目チャネルのｋ空間データにおける各サンプリング位置のサンプル値と、注目チャネル以外のチャネルを含む、全チャネルのｋ空間データにおける周辺のサンプリング位置のサンプル値の重み付き線形和との差を小さくする項をデータ拘束として、重み係数を推定する。フルサンプリングされた領域では全ての値が既知であるので、例えば、最小二乗法を利用することで、重み係数を推定することができる。最小二乗法を最小化問題のデータ項として表すと、重み係数導出部１３３ｃは、フルサンプリングされた領域における、全てのサンプリング位置「ｋ」に対する、（４）式の左辺と右辺との差を２乗した値（誤差）の総和として表せる。この総和に対する最小化問題を解くことにより、重み係数ｗ（ｂ，ｃ，ｃout）を推定することができる。

ここで推定された重み係数は、ステップＳ４において、補正用ＭＲ画像の重み付け線形和で用いられる重み係数となる。フーリエ変換は線形変換であるので、（４）式と同等のデータ拘束を用いて推定された重み係数は、フーリエ変換後の画像空間においても（４）式の関係を維持し、補正用ＭＲ画像の重み付け線形和の重み係数としてそのまま利用することができる。なお、（４）式において、重み係数ｗが、移動量ｂとコイルｃとの組合せ毎に求められるということは、補正用ＭＲ画像の重み付け線形和で用いられる重み係数が、平行移動の規則とチャネルとの組合せ毎に求められることと対応する。

ステップＳ４：出力ＭＲ画像算出部１３３ｄは、出力ＭＲ画像を算出する。具体的には、出力ＭＲ画像算出部１３３ｄは、ステップＳ２で生成した各補正用ＭＲ画像に対して、ステップＳ３で導出した重み係数をかけることで、重み付け加算を行い、図５に示すように、出力ＭＲ画像を算出する。なお、複数のチャネルが入力され、複数のチャネルについてそれぞれ出力ＭＲ画像を算出する場合には、ステップＳ２〜Ｓ４の処理を全チャネルについて適用すれば良い。

次に、第１の実施形態における処理手順の詳細を説明する。図６は、第１の実施形態におけるｋ空間データを説明するための図であり、図７〜９は、第１の実施形態における処理手順を示すフローチャートである。以下では、ｋ空間データが２次元の場合の例を説明するが、実施形態はこれに限られるものではなく、ｋ空間データが３次元の場合でも、以下の実施形態を同様に適用することができる。

ステップＳ１０１：補正用ＭＲ画像生成部１３３ｂが、記憶部１３２から処理対象の撮像ｋ空間データを、各チャネル（チャネルｋ＝１，・・・，Ｃ）について読み出す。第１の実施形態において、この撮像ｋ空間データは、図６に示すように、位相エンコード（Phase Encode）方向では、ｋ空間の中心領域において、フルサンプリングと同一の解像度でサンプリングが行われ、周辺領域において、１ラインずつ規則的に間引きサンプリングが行われたものである。また、図６では図示を省略しているが、この撮像ｋ空間データは、リードアウト（Readout）方向では、フルサンプリングと同一の解像度でサンプリングが行われたものである。そして、補正用ＭＲ画像生成部１３３ｂは、撮像ｋ空間データのうち既知サンプルが存在しないサンプリング位置のサンプル値を「ゼロ」で埋めることで、フルサンプリングと同一の解像度のｋ空間データを生成する。

ステップＳ１０２：ステップＳ１０１で生成された各チャネルの撮像ｋ空間データについて、ステップＳ１０２−１〜Ｓ１０２−４に示す処理が行われる。ステップＳ１０２−１〜Ｓ１０２−４に示す処理は、例えば、チャネル数が「８」であれば、８回行われる。

ステップＳ１０２−１：補正用ＭＲ画像生成部１３３ｂは、注目するチャネルの撮像ｋ空間データにフーリエ変換を適用し、第１補正用ＭＲ画像を生成する。

ステップＳ１０２−２：補正用ＭＲ画像生成部１３３ｂは、注目するチャネルの撮像ｋ空間データを用いて別途定めた規則に従った平行移動を行い、平行移動後のシフトｋ空間データを再構成することで、第２補正用ＭＲ画像を生成する。また、補正用ＭＲ画像生成部１３３ｂは、各チャネルに対する各平行移動量について、ステップＳ１０２−３１〜Ｓ１０２−３３の処理を行う。例えば、平行移動の規則が、（Readout、Phase Encode）＝（０，１）（０，−１）（１，１）（１，−１）（−１，１）（−１，−１）の６通りの平行移動量を定めているとする。これは、位相エンコード（Phase Encode）方向で２種類（位相を＋１及び−１とする）、各位相エンコード（Phase Encode）方向の移動に対してリードアウト（Read Out）方向で３種類（位相を＋０、＋１、−１とする）とした、２×３の６通りの平行移動である。

この場合、補正用ＭＲ画像生成部１３３ｂは、１つのチャネルについて、６通りの平行移動量それぞれについて、ステップＳ１０２−３１〜Ｓ１０２−３３の処理を行う。また、チャネル数が８であれば、補正用ＭＲ画像生成部１３３ｂは、４８回、ステップＳ１０２−２１〜Ｓ１０２−２３の処理を行う。なお、第１の実施形態においては、６通りの平行移動を行う例を説明したが、これは一例に過ぎず、例えば、１通りの平行移動を行うものであってもよい。

ステップＳ１０２−２１：補正用ＭＲ画像生成部１３３ｂは、撮像ｋ空間データを、注目する平行移動量だけずらす。例えば、補正用ＭＲ画像生成部１３３ｂは、平行移動の規則が（Readout、Phase Encode）＝（０，１）の場合、図６に示すように、まず、撮像ｋ空間データ全体を、位相エンコード方向に「１」平行移動する。

ステップＳ１０２−２２：次に、補正用ＭＲ画像生成部１３３ｂは、平行移動後のシフトｋ空間データのうち、撮像ｋ空間データにおける既知サンプルのサンプリング位置に「ゼロ」を設定して、マスク済シフトｋ空間データを生成する。このマスク済シフトｋ空間データは、撮像ｋ空間データにおける未知サンプルのサンプリング位置にのみサンプル値が存在するｋ空間データである。

ステップＳ１０２−２３：そして、補正用ＭＲ画像生成部１３３ｂは、ステップＳ１０２−２２で生成されたマスク済シフトｋ空間データにフーリエ変換を適用し、第２補正用ＭＲ画像を生成する。

ステップＳＳ１０２−３：こうして、補正用ＭＲ画像生成部１３３ｂは、全てのチャネルについて、それぞれ６通りの平行移動量に従った第２補正用ＭＲ画像を生成する。例えば、チャネル数が８であれば、補正用ＭＲ画像生成部１３３ｂは、４８の第２補正用ＭＲ画像を生成する。ここで、撮像ｋ空間データから生成された第１補正用ＭＲ画像をＸ0、各チャネルｃ＝１，・・・，Ｃのマスク済シフトｋ空間データから生成された第２補正用ＭＲ画像をＸi（ｉ＝１，・・・，Ｃ）とすると、出力ＭＲ画像Ｘは、（６）式において重み係数ｗ0＝１とした、（７）式で表現することができる。

この例において、補正用ＭＲ画像は、第１補正用ＭＲ画像Ｘ0及び第２補正用ＭＲ画像Ｘi（ｉ＝１，・・・，Ｃ）である。フーリエ変換は線形変換であるから、出力ＭＲ画像Ｘを、第１補正用ＭＲ画像Ｘ0及び第２補正用ＭＲ画像Ｘiの重み付き加算として算出することは、撮像ｋ空間データ及びマスク済シフトｋ空間データに対して重み係数ｗiの重み付き加算を算出し、それをフーリエ変換によって再構成することと等価である。また、マスク済シフトｋ空間データは、撮像ｋ空間データにおける未知サンプルのサンプリング位置にのみ「ゼロ」でない値を持ち、また、撮像ｋ空間データ全体を位相エンコード方向に平行移動する操作によって生成されたものである。この操作は、ｋ空間の線形補間と数学的には同等の操作である。したがって、この場合、重み係数は、（４）式と同等のデータ拘束を用いた最小化問題を解くことによって推定することができる。なお、後述するように、（４）式のデータ拘束以外の拘束式と組み合わせることで、ＰＩの性能を向上することが可能である。

こうして、重み係数導出部１３３ｃが、（４）式と同等のデータ拘束を用いて重み係数を推定する。例えば、重み係数導出部１３３ｃは、データ拘束の誤差が小さくなるように、第２補正用ＭＲ画像のそれぞれに対する重み係数を推定する。また、第１の実施形態においては、第１補正用ＭＲ画像に対する重み係数は「１」とする。

ステップＳ１０２−４：そして、出力ＭＲ画像算出部１３３ｄは、各補正用ＭＲ画像について対応する重み係数をかけ、その総和を、注目するチャネルの出力ＭＲ画像とする。

ステップＳ１０３：また、出力ＭＲ画像算出部１３３ｄは、全チャネルの出力ＭＲ画像を統合して、最終的な出力ＭＲ画像を得る。例えば、出力ＭＲ画像算出部１３３ｄは、全チャネルの出力ＭＲ画像の２乗和の平均を算出して、最終的な出力ＭＲ画像とする。なお、ここでいう「最終出力ＭＲ画像」とは、説明の便宜上、「出力ＭＲ画像」と区別するための表現に過ぎない。例えば、最終出力ＭＲ画像に対して後処理が施されたＭＲ画像が、最終的に表示部１３５等に出力されるものであってもよい。

（第１の実施形態の変形例）
上述した第１の実施形態においては、ｋ空間データ全体に対して同じ種類の平行移動の規則を適用する手法を説明したが、実施形態はこれに限られるものではない。各サンプリング位置に対して平行移動の規則を定め、同一の平行移動の規則が適用されるサンプリング位置をグループ化した上で、第２補正用画像Ｘiを生成する手法を適用することもできる。例えば、これまで説明した処理手順のうち、一部の処理手順を以下に示す処理手順に置き換えて実現することができる。なお、以下に説明する手法は、各サンプリング位置に対して平行移動の規則を定める手法であるので、サンプリング密度が可変の場合にも対応することができる。

ステップＳ２０１：補正用ＭＲ画像生成部１３３ｂは、撮像ｋ空間データにおける未知サンプルのサンプリング位置それぞれについて、平行移動の規則を定める。なお、平行移動の規則は、事前に定めておいてもよい。

ステップＳ２０２：補正用ＭＲ画像生成部１３３ｂは、シフトｋ空間データにおいて、同一の平行移動の規則が適用されるサンプリング位置を、同一の重み係数ｗiを有するサンプリング位置の集合としてグループ化する。図１０は、第１の実施形態の変形例を説明するための図である。例えば、図１０に示すように、２次元Ｃａｒｔｅｓｉａｎサンプリングによりｋ空間データを収集する際に、中心部をフルサンプリング、周辺部のうち中心部に近い一部を１／２ライン毎のサンプリング、それ以外の部分を１／３ライン毎のサンプリングで収集する例を考える。このとき、１／２ライン毎のサンプリングを行う領域を第１のグループ（Group１）、１／３ライン毎のサンプリングを第２のグループ（Group２）として２つにグループ化したうえで、グループ１については（Readout、Phase Encode）＝（０，１）（０，−１）（１，１）（１，−１）（−１，１）（−１，−１）という平行移動の規則を利用し、グループ２については（Readout、Phase Encode）＝（０，１）（０，−１）（１，１）（１，−１）（−１，１）（−１，−１）（０，２）（０，−２）（１，２）（１，−２）（−１，２）（−１，−２）という平行移動の規則を利用する。この際、グループ外に平行移動されたサンプルはマスクの対象（即ち、サンプル値を「ゼロ」）とする。また、平行移動量が同一であっても別のグループである場合には、別々の平行移動として扱う。なお、補正用ＭＲ画像生成部１３３ｂは、必ずしも、同一の平行移動の規則が適用されるサンプリング位置全てを１つのグループに所属させる必要はない。

ステップＳ２０３：補正用ＭＲ画像生成部１３３ｂは、各グループ、各平行移動の規則、各チャネルについて、第２補正用ＭＲ画像Ｘiを生成する。上述したように、グループ毎に異なる平行移動の規則が適用される場合、第２補正用ＭＲ画像Ｘiは、平行移動の規則の数の全てのグループ分の総和と、チャネル数との積の数分、生成される。

ステップＳ２０４：重み係数導出部１３３ｃは、各平行移動の規則に対応した（なお、上述したように、平行移動量が同一であっても別のグループである場合には、別々の平行移動として扱う）、（４）式と同様のデータ拘束式を生成し、生成したデータ拘束式を用いて、重み係数ｗiを推定する。具体的には、重み係数導出部１３３ｃは、各チャネル、各グループの各平行移動の規則に対して、それぞれ、フルサンプリング領域を用いて（４）式の拘束式を生成し、誤差の総和である（４）式に基づくデータ拘束を最小化するような重み係数を推定する。例えば、重み係数導出部１３３ｃは、各チャネルにおいて、グループ１の領域に適用される平行移動の規則、及び、グループ２の領域に適用される平行移動の規則、それぞれに対して、フルサンプリング領域を用いて（４）式の拘束式を生成し、誤差の総和である（４）式に基づくデータ拘束を最小化するような重み係数を推定する。

ステップＳ２０５：出力ＭＲ画像算出部１３３ｄは、各補正用ＭＲ画像について対応する重み係数をかけ、その総和を、注目するチャネルの出力ＭＲ画像とする。

（第１の実施形態の効果）
上述してきたように、第１の実施形態によれば、フルサンプリングよりも少ないサンプル数で収集されたｋ空間データから、高画質の出力ＭＲ画像を得ることができる。

（第２の実施形態）
上述した第１の実施形態においては、サンプリングされたデータから得られる拘束式を用いて（データ拘束）重み係数を導出する手法を説明したが、実施形態はこれに限られるものではない。データ拘束の手法は、第１の実施形態で説明した手法に限られるものではない。例えば、重み係数は、データ拘束とプライア（事前知識）による拘束との組合せで導出されてもよい。また、例えば、重み係数は、プライアによる拘束のみで導出されてもよい。

データ拘束を改めて説明する。データ拘束は、出力ＭＲ画像に対応するｋ空間データをＹで表現し、各第２補正用ＭＲ画像Ｘiに対応するｋ空間データをＹiで表現し、ｋ空間データ全体からフルサンプリングに利用する領域のみを抽出するための行列をＭで表現し、重み係数を並べたベクトルをｗで表現した場合に、次の（８）式で表現される最適化問題を解くことに相当する。なお、この（８）式では、データ拘束の各項について等しい重みを考えているが、重みを考慮してもよい。

また、重み係数導出部１３３ｃは、重み係数の導出にあたり、プライア（事前知識）による拘束を加えてもよい。上述したように、重み係数ｗは、第２補正用ＭＲ画像Ｘiの係数であり、（５）式から、出力ＭＲ画像Ｘは、重み係数ｗの関数である。そこで、出力ＭＲ画像をＸ（ｗ）で表し、またデータ拘束からの誤りを評価する関数をＥdata（ｗ）で表現し、プライア拘束からの誤りを評価する関数をＥprior（ｗ）で表現し、データ拘束とプライア拘束とのバランスを調整するためのパラメータをλで表現すると、重み係数は、（９）式で導出される。

重み係数導出部１３３ｃは、この（９）式の最適化問題を解いて、重み係数ｗを推定してもよい。そして、出力ＭＲ画像算出部１３３ｄは、出力ＭＲ画像（ｗ）を算出して出力する。重みｗの推定問題は、（１０）式のエネルギー関数の最小化問題であるので、既存のエネルギー最小化の枠組みがそのまま利用できる。

重み係数導出部１３３ｃは、上述したように、Ｅprior（ｗ）＝０であり、Ｅdata（ｗ）が２次の拘束式であれば、最小二乗解を求めればよい。そうでないエネルギー関数については、エネルギー関数の形に依存するが、Gradient descent法、Conjugate gradient法やNonlinear conjugate gradient法のいずれかが利用できることが多い。なお、ここで例示した方法以外のエネルギー最小化方法を用いてもよい。

続いて、プライア拘束を説明する。プライア拘束からの誤りを評価する関数Ｅprior（ｗ）を、出力ＭＲ画像Ｘ（ｗ）の関数として与えることを考える。数式で表現すると、ｆをＭＲ画像の関数として、Ｅprior（ｗ）＝ｆ（Ｘ（ｗ））で表現することを考える。このとき、関数ｆは、画像に対する関数であるから、画像再構成やノイズ除去で用いられるプライア拘束と同様のプライアを利用できる。例えば、以下のプライア拘束を適用できる。なお、以下に説明する各手法は、適宜組み合わせることができる。

（手法１）重み係数導出部１３３ｃは、テンプレート画像（基準画像）によるレギュラライズのプライア拘束を適用することができる。例えば、重み係数導出部１３３ｃは、関数ｆを（１１）式のように定めたプライアを適用する。なお、（１１）式ではＬ２ノルムを用いているが、圧縮センシングのように、別のノルム（例えば（１２）式のＬ１ノルム）を用いてもよい。

テンプレート画像Ｘ0としては次のものを利用することができる。例えば、ＭＲＩ装置１００において時系列に沿って連続的に収集されたデータが存在する場合、重み係数導出部１３３ｃは、例えば、時系列順で１つ前のフレームの画像を、テンプレート画像Ｘ0として利用することができる。なお、１つ前のフレームの画像は一例に過ぎない。時系列順で一定程度の相関が認められるフレームの画像であれば、１つ前のフレームのみならず、１つ後のフレーム、数フレーム前のフレーム、数フレーム後のフレームについても、テンプレート画像Ｘ0として利用することができる。

また、例えば、ＭＲＩ装置１００において、複数スライスを含む領域を撮像領域とするスラブ撮像が行われる場合、重み係数導出部１３３ｃは、例えば、同一スラブ内であって、且つ、他のスライスの画像を、テンプレート画像Ｘ0として利用することができる。

また、例えば、ＭＲＩ装置１００において、複数のプロトコルが連続して実行される場合、重み係数導出部１３３ｃは、他のプロトコルで収集され、生成された画像を、テンプレート画像Ｘ0として利用することができる。例えば、ＭＲＩ装置１００による撮像では、位置決め画像を収集するプロトコル、感度マップを収集するプロトコル、といった準備スキャンや、Ｔ１強調画像を収集するプロトコル、Ｔ２強調画像を収集するプロトコル、といったイメージングスキャン等、複数のプロトコルが連続して実行される場合がある。このような場合に、重み係数導出部１３３ｃは、例えば、位置決め画像を収集するプロトコルで収集された位置決め画像を、テンプレート画像Ｘ0として利用することができる。また、例えば、重み係数導出部１３３ｃは、例えば、イメージングスキャン内のプロトコル間で、得られた画像をテンプレート画像Ｘ0として利用することができる。

その他、テンプレート画像Ｘ0は、必ずしも被検体本人の画像でなくてもよく、何等かの学習等によって得られた画像であってもよい。

（手法２）重み係数導出部１３３ｃは、特定領域の信号値を「ゼロ」であると仮定するプライア拘束を適用することができる。例えば、重み係数導出部１３３ｃは、ＭＲ画像Ｘに対し、信号値をゼロとする領域を抽出するための行列をＭzeroで表すとき、関数ｆを（１３）式のように定めたプライアを適用する。

なお、ここでいう「特定領域」とは、例えば、被検体の外側や、脳の画像を撮像する場合の四隅等のように、信号値が「ゼロ」であることが明らかである領域である。

（手法３）重み係数導出部１３３ｃは、スパース変換（例えば、wavelet変換やtotal variation）を利用したプライア拘束を適用することができる。このプライア拘束は、「圧縮センシング技術」という名称で知られている。例えば、重み係数導出部１３３ｃは、Ψをスパース変換とするとき、関数ｆを（１４）式のように定めたプライアを適用する。

また、重み係数導出部１３３ｃは、２種以上のスパース変換を利用してもよい。Ψ1、Ψ2をスパース変換とし、λ1、λ2をパラメータとするとき、重み係数導出部１３３ｃは、関数ｆを（１５）式のように定めたプライアを適用する。

（手法４）重み係数導出部１３３ｃは、上記手法１から手法３までのプライア拘束を、２つ以上の組み合わせたプライア拘束を適用することができる。この場合、重み係数導出部１３３ｃは、関数ｆを、上記各式で求められる関数ｆの重み付き加算で定める。

これらの例に示されるように、画像空間に特有の知識を、プライア拘束として、データ拘束と組み合わせる結果、より高精度な画像再構成を実現することができる。なお、この処理は、出力ＭＲ画像を複数チャネルデータに基づく部分空間で表現したうえで、部分空間内でデータ拘束とプライア拘束を利用して出力ＭＲ画像を推定していると解釈できる。

（その他の実施形態）
実施形態は、上述した第１及び第２の実施形態に限られるものではない。

（定式化 時系列撮像に対する平行移動パターン）
上述した実施形態は、時系列のｋ空間データを再構成する場合にも、同様に適用することができる。時系列のｋ空間データを再構成する場合、「ｋ−ｔ ＧＲＡＰＰＡ」と呼ばれる手法が用いられることもある。図１１は、その他の実施形態におけるｋ空間データを説明するための図である。

図１１に示すように、例えば、ｋ空間データ収集部１３３ａは、各時刻において１／４の間引きを行いながらサンプリングを行うが、そのサンプリング位置を、１単位時刻毎にずらしながらサンプリングを行っている。図１１において、サンプリング位置８０２は、サンプル値を有しないサンプリング位置を示し、サンプリング位置８０３は、サンプリング位置８０２の推定に利用する、サンプル値を有するサンプリング位置を示す。また、キャリブレーション用データとして、領域８０１は、フルサンプリングされているものとする。

例えば、サンプリング位置８０２は、サンプリング位置８０３の各点を、位相エンコード（Phase Encode）方向及び時間方向に平行移動することで生成できるものと考える。そして、補正用ＭＲ画像生成部１３３ｂは、各平行移動パターンに対してフーリエ変換を適用して、第２補正用ＭＲ画像Ｘiを生成する。重みｗiの推定において、重み係数導出部１３３ｃは、データ拘束として、フルサンプリングされた領域８０１において、時空間で相対位置が等しい、サンプリング位置８０４及びサンプリング位置８０５を利用することができる。

（分割領域を用いたデータ拘束）
上述した実施形態は、ｋ空間の領域が分割されていても、適用することができる。Ｒａｄｉａｌ ＧＲＡＰＰＡとは、角度を変化させながらｋ空間の中心を通る線分を順にスキャンすることでｋ空間データを（非Ｃａｒｔｅｓｉａｎで）取得するＲａｄｉａｌスキャンを用いて得られたｋ空間データに対するＧＲＡＰＰＡの１実装である。この手法では、スキャンを別途定めた基準によってグループ（セグメント）に分割し、各セグメント内のデータは近似的にＣａｒｔｅｓｉａｎで取得されたものであると考えて、個々のセグメントについてＧＲＡＰＰＡが適用される。

Ｒａｄｉａｌ ＧＲＡＰＰＡのように領域分割されている場合には、例えば、補正用ＭＲ画像生成部１３３ｂは、入力として与えられた撮像ｋ空間データに、Non−uniformフーリエ変換（ＮｕＤＦＴ）による再構成を適用し、得られたＭＲ画像を、第１補正用ＭＲ画像Ｘ0とする。また、補正用ＭＲ画像生成部１３３ｂは、個々の分割領域内で生成した補正用のｋ空間データを再構成し、得えられたＭＲ画像を、第２補正用ＭＲ画像Ｘiとする。また、重み係数導出部１３３ｃは、上述した実施形態と同様に、個々の分割領域に対して、データ拘束を与えればよい。

（任意のＴｒａｊｅｃｔｏｒｙに対する補正用ＭＲ画像の生成手法）
非Ｃａｒｔｅｓｉａｎサンプリングの種類によっては、上記のＣａｒｔｅｓｉａｎやＲａｄｉａｌのために設計されたＧＲＡＰＰＡが適用できない場合がある。図１２〜１４は、その他の実施形態におけるｋ空間データを説明するための図である。

例えば、図１２に示すＳｐｉｒａｌスキャンでは、先に述べたＣａｒｔｅｓｉａｎやＲａｄｉａｌスキャンに対する手法を適用することはできない。そのような場合には、予め、ｋ空間上での平行移動パターンを１グループ以上定めておく。補正用ＭＲ画像生成部１３３ｂは、各グループに対し、サンプル値を有する点の集合に平行移動パターンを適用し、ｋ空間上の点の集合を生成する。そして、補正用ＭＲ画像生成部１３３ｂは、生成された点の集合に対してNon−uniformフーリエ変換（ＮｕＤＦＴ）による再構成を適用し、第２補正用画像Ｘiを得る（例えば、図１４を参照）。

図１２においてサンプル値を有する点の集合９０１（図１３において、黒丸の渦巻き状の点の集合）を、ｋx方向に単位量だけ平行移動させた例を、図１３に示す。図１３において、パターン模様の丸の渦巻き状の点の集合１００１は、点の集合９０１を平行移動させることによって生成した点であり、これらの点に対して各チャネルの重み係数を総和するものと考える。なお、点の集合１００１は点の集合９０１全てを含む必要はなく、いくつかの点を除外してもよい。

なお、ＮｕＤＦＴには何らかの形でサンプリング密度の不均一を補正する操作が含まれる（例えば、グリッディングと呼ばれる手法を介してＣａｒｔｅｓｉａｎのＦＦＴに帰着させる手法を利用するのであれば、密度の不均一さを位置依存の重み値の除算によって補正する操作を行う）が、その部分については、全第２補正用ＭＲ画像Ｘiを生成するために利用した、平行移動されたｋ空間データも含むｋ空間全てをサンプリングデータとして、そのサンプリング密度を算出し、それに基づいて不均一の補正を行う。

（ｋ空間データの生成）
上述した実施形態においては、撮像ｋ空間データから別のｋ空間データを生成する際に平行移動の処理を施す例を説明したが、実施形態はこれに限られるものではない。例えば、ＭＲＩ装置１００において時系列に沿って連続的に収集されたデータが存在する場合、補正用ＭＲ画像生成部１３３ｂは、撮像ｋ空間データとは別のｋ空間データとして、他の時点のｋ空間データを用いてもよい。例えば、補正用ＭＲ画像生成部１３３ｂは、時系列順で隣接している時点のｋ空間データを用いて、撮像ｋ空間データのうちの未知サンプルにサンプル値を与え、別のｋ空間データを生成してもよい。

また、例えば、補正用ＭＲ画像生成部１３３ｂは、撮像ｋ空間データに対してＧＲＡＰＰＡを施すことで、撮像ｋ空間データのうちの未知サンプルにサンプル値を与え、別のｋ空間データを生成してもよい。

また、上述した実施形態においては、撮像ｋ空間データに対して平行移動の処理を施す手法を説明したが、更に、これに主成分分析（ＰＣＡ；Principal Component Analysis）の考え方を導入してもよい。例えば、上述した実施形態においては、全てのチャネルで得られた撮像ｋ空間データに対して平行移動の処理を施し、これらを全て、重み付け加算の対象としていた。しかしながら、実施形態はこれに限られるものではなく、主成分分析の考え方に基づいて、重み付け加算の対象とするチャネルを削減してもよい。この結果、計算量を減らすことができ、また、ノイズの影響を軽減する効果も期待できる。

また、上述した実施形態においては、ｋ空間で撮像ｋ空間データに対して平行移動の処理を施す手法を説明したが、ｋ空間ではなく、フーリエ変換後の画像空間で平行移動の処理を施してもよい。例えば、画像空間が１次元である場合を例に挙げて説明する。この場合、ｋ空間における元のデータをｋ（ｘ）で表すと、ｋ空間でａだけ平行移動されたデータはｋ（ｘ−ａ）で表される。元のデータｋ（ｘ）の逆離散フーリエ変換をＴ（ｋ（ｘ））で表すと、平行移動されたデータＴ（ｋ（ｘ−ａ））は、（１６）式で表される。ここで、（１６）式におけるｅｘｐ（２πｉａｘ／Ｎ）は、自然対数の底（ネイピア数）ｅの（２πｉａｘ／Ｎ）乗を指す。

したがって、各位置ｘに対してｅｘｐ（２πｉａｘ／Ｎ）を乗じることで、画像空間で平行移動を行うことができる。なお、上述したような画像空間での平行移動の手法を２以上の座標軸のそれぞれに対して行うことで、画像空間が２次元以上である場合にも、画像空間上での平行移動を行うことができる。

（撮像ｋ空間データ）
上述した実施形態においては、撮像ｋ空間データとして、ｋ空間の中心領域においてはフルサンプリングと同一の解像度でサンプリングが行われ、ｋ空間の周辺領域においては間引きサンプリングが行われたｋ空間データを想定したが、実施形態はこれに限られるものではない。撮像ｋ空間データは、フルサンプリングに対応するｋ空間上のサンプリング位置の一部のサンプリング位置にサンプル値を有するｋ空間データであればよい。また、間引きサンプリングを行う場合のサンプリング間隔は、可変（サンプリング位置によってサンプリング間隔が異なるもの）であってもよい。

（画像処理装置）
また、上述した実施形態においては、ＭＲＩ装置１００のｋ空間データ収集部１３３ａがｋ空間データを収集し、収集したｋ空間データを対象に、補正用ＭＲ画像生成部１３３ｂ、重み係数導出部１３３ｃ、及び出力ＭＲ画像算出部１３３ｄによる処理が行われる例を説明したが、実施形態はこれに限られるものではない。例えば、ＭＲＩ装置１００とは異なる画像処理装置が、フルサンプリングよりも少ないサンプル数のｋ空間データを対象に、補正用ＭＲ画像生成部１３３ｂ、重み係数導出部１３３ｃ、及び出力ＭＲ画像算出部１３３ｄによる処理と同等の処理を実行してもよい。この場合、画像処理装置は、ＭＲＩ装置１００やその他の外部装置からネットワーク経由でｋ空間データを受け取る、あるいは記憶媒体経由でｋ空間データを受け取る等する。また、フルサンプリングよりも少ないサンプル数のｋ空間データがどのようにＭＲＩ装置１００において収集されたかは、画像処理装置における処理とは直接関係するものではない。即ち、フルサンプリングよりも少ないサンプル数のｋ空間データは、必ずしも間引きサンプリングによって収集された場合に限られず、フルサンプリングで収集されたｋ空間データが、事後的に間引かれたものであってもよい。結局、処理の対象として、フルサンプリングよりも少ないサンプル数のｋ空間データが存在すればよい。

（プログラム）
上述した実施形態の中で示した処理手順に示された指示は、ソフトウェアであるプログラムに基づいて実行されることが可能である。汎用の計算機が、このプログラムを予め記憶しておき、このプログラムを読み込むことにより、上述した実施形態のＭＲＩ装置や画像処理装置による効果と同様な効果を得ることも可能である。上述した実施形態で記述された指示は、コンピュータに実行させることのできるプログラムとして、磁気ディスク（フレキシブルディスク、ハードディスク等）、光ディスク（ＣＤ−ＲＯＭ、ＣＤ−Ｒ、ＣＤ−ＲＷ、ＤＶＤ−ＲＯＭ、ＤＶＤ±Ｒ、ＤＶＤ±ＲＷ等）、半導体メモリ、又はこれに類する記録媒体に記録される。コンピュータ又は組み込みシステムが読み取り可能な記憶媒体であれば、その記憶形式は何れの形態であってもよい。コンピュータは、この記録媒体からプログラムを読み込み、このプログラムに基づいてプログラムに記述されている指示をＣＰＵで実行させれば、上述した実施形態のＭＲＩ装置や画像処理装置と同様な動作を実現することができる。もちろん、コンピュータがプログラムを取得する場合又は読み込む場合はネットワークを通じて取得又は読み込んでもよい。

また、記憶媒体からコンピュータや組み込みシステムにインストールされたプログラムの指示に基づきコンピュータ上で稼働しているＯＳ（オペレーティングシステム）や、データベース管理ソフト、ネットワーク等のＭＷ（ミドルウェア）等が、上述した実施形態を実現するための各処理の一部を実行してもよい。更に、記憶媒体は、コンピュータあるいは組み込みシステムと独立した媒体に限らず、ＬＡＮ（Local Area Network）やインターネット等により伝達されたプログラムをダウンロードして記憶又は一時記憶した記憶媒体も含まれる。また、記憶媒体は１つに限られず、複数の媒体から、上述した実施形態における処理が実行される場合も、実施形態における記憶媒体に含まれ、媒体の構成は何れの構成であってもよい。

なお、実施形態におけるコンピュータ又は組み込みシステムは、記憶媒体に記憶されたプログラムに基づき、上述した実施形態における各処理を実行するためのものであって、パソコン、マイコン等の１つからなる装置、複数の装置がネットワーク接続されたシステム等の何れの構成であってもよい。また、実施形態におけるコンピュータとは、パソコンに限らず、情報処理機器に含まれる演算処理装置、マイコン等も含み、プログラムによって実施形態における機能を実現することが可能な機器、装置を総称している。

以上述べた少なくとも一つの実施形態の磁気共鳴イメージング装置及び画像処理装置によれば、画質を向上することができる。

本発明のいくつかの実施形態を説明したが、これらの実施形態は、例として提示したものであり、発明の範囲を限定することは意図していない。これら実施形態は、その他の様々な形態で実施されることが可能であり、発明の要旨を逸脱しない範囲で、種々の省略、置き換え、変更を行うことができる。これら実施形態やその変形は、発明の範囲や要旨に含まれると同様に、特許請求の範囲に記載された発明とその均等の範囲に含まれるものである。

１００ ＭＲＩ装置
１２０ シーケンス制御部
１３３ 制御部
１３３ａ ｋ空間データ収集部
１３３ｂ 補正用ＭＲ画像生成部
１３３ｃ 重み係数導出部
１３３ｄ 出力ＭＲ画像算出部

Claims (10)

1. ｋ空間上のフルサンプリング位置の一部にサンプル値を有する第１のｋ空間データのうち、前記サンプル値を有しないサンプリング位置の少なくとも一部にサンプル値を与えて第２のｋ空間データを生成するｋ空間データ生成部と、
   前記第１のｋ空間データから第１画像を生成し、前記第２のｋ空間データから第２画像を生成する画像生成部と、
   前記第１画像及び前記第２画像に対する重み係数を導出する導出部と、
   前記第１画像及び前記第２画像に対して、前記重み係数を用いた重み付け加算を行い、磁気共鳴画像を算出する算出部と
   を備える、磁気共鳴イメージング装置。
2. 前記ｋ空間データ生成部は、前記サンプル値を有しないサンプリング位置に、前記第１のｋ空間データのうちサンプリング位置に既に与えられているサンプル値に基づいて定まる値を与えることで、前記第１のｋ空間データから前記第２のｋ空間データを生成する、請求項１に記載の磁気共鳴イメージング装置。
3. 前記ｋ空間データ生成部は、前記サンプル値を有しないサンプリング位置に、特定の位置関係にあるサンプリング位置における、前記第１のｋ空間データのサンプル値を与えることで、前記第１のｋ空間データから第２のｋ空間データを生成する、請求項１に記載の磁気共鳴イメージング装置。
4. 前記ｋ空間データ生成部は、前記第２のｋ空間データのうち、前記第１のｋ空間データではサンプル値を有するサンプリング位置に、特定の値を与えることで、前記第１のｋ空間データから第２のｋ空間データを生成する、請求項１〜３のいずれか１つに記載の磁気共鳴イメージング装置。
5. 前記ｋ空間データ生成部は、前記第１のｋ空間データのうち、前記サンプル値を有しないサンプリング位置に、前記サンプル値を有するサンプリング位置のサンプル値を平行移動により与えることで、前記第１のｋ空間データから第２のｋ空間データを生成する、請求項１〜３のいずれか１つに記載の磁気共鳴イメージング装置。
6. 前記導出部は、注目チャネルのｋ空間データにおける各サンプリング位置のサンプル値が、前記注目チャネル以外のチャネルを含む、全チャネルのｋ空間データにおける周辺のサンプリング位置のサンプル値の重み付き線形和と等しいことを拘束式として、前記重み係数を導出する、請求項１〜３のいずれか１つに記載の磁気共鳴イメージング装置。
7. 前記導出部は、前記磁気共鳴画像と同一の解像度を有する基準画像における画素値と、前記磁気共鳴画像における画素値との差を小さくするように、前記重み係数を決定する、請求項１〜３のいずれか１つに記載の磁気共鳴イメージング装置。
8. 前記導出部は、前記磁気共鳴画像内の特定の領域の信号値がゼロであると仮定して、前記重み係数を決定する、請求項１〜３のいずれか１つに記載の磁気共鳴イメージング装置。
9. ｋ空間上のフルサンプリング位置の一部にサンプル値を有する第１のｋ空間データのうち、前記サンプル値を有しないサンプリング位置の少なくとも一部にサンプル値を与えて第２のｋ空間データを生成するｋ空間データ生成部と、
   前記第１のｋ空間データから第１画像を生成し、前記第２のｋ空間データから第２画像を生成する画像生成部と、
   前記第１画像及び前記第２画像に対する重み係数を導出する導出部と、
   前記第１画像及び前記第２画像に対して、前記重み係数を用いた重み付け加算を行い、磁気共鳴画像を算出する算出部と
   を備える、画像処理装置。
10. フルサンプリングに対応するｋ空間上のサンプリング位置の一部のサンプリング位置にサンプル値を有する第１のｋ空間データから、第２のｋ空間データを生成するｋ空間データ生成部と、
    前記第１のｋ空間データから、折り返しのある第１画像を生成し、前記第２のｋ空間データから第２画像を生成する画像生成部と、
    前記第１画像及び前記第２画像を用いて、折り返しのない磁気共鳴画像を算出する算出部と
    を備える、磁気共鳴イメージング装置。