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

Description

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

パラレルイメージング（ＰＩ（Parallel Imaging））を用いた磁気共鳴イメージングにおいては、例えば、複数のコイルを用いてｋ空間データが収集され、その後、画像再構成処理が行われて、医用画像（出力対象画像）が生成される。この画像再構成処理は、収集されたｋ空間データから、出力対象画像を復元する逆問題である。逆問題としての画像再構成処理は、最適化問題の解が入力データの小さな変化に敏感に依存する場合がある。この結果、例えば、最小二乗法により逆問題を解く場合、解が一意に定まらず発散する場合がある。

このような場合、解の安定性を高めるため、目的関数に対して所定の追加の項をペナルティ項として導入する操作、すなわち正則化（regularization）が、しばしば行われる。

しかし、出力対象の画像に対して単純に所定のノルムを用いて正則化を行った場合、再構成後の出力対象画像の画質が、かえって低下しまうことがある。

K.P.Pruessmann et al., "SENSE : Sensitivity Encoding for Fast MRI," MRM, Vol.42, pp.952-962, 1999.

本発明が解決しようとする課題は、パラレルイメージングにおける出力対象画像の画質を向上することができる磁気共鳴イメージング装置及び画像処理装置を提供することである。

実施形態に係る磁気共鳴イメージング装置は、シーケンス制御部と、収集部と、推定部と、生成部とを備える。シーケンス制御部は、アンダーサンプリングを行いながらパルスシーケンスを実行する。収集部は、前記パルスシーケンスを用いてｋ空間データを収集する。推定部は、前記ｋ空間データを用いて生成されたデータに基づいて、出力対象画像の各画素の画素位置が観察対象領域に含まれているか否かを推定する。生成部は、前記推定部が観察対象に含まれていると推定した画素位置に対する重みと、前記推定部が観察対象領域に含まれないと推定した画素位置に対する重みとが異なる値となるような正則化項を用いて、単時相において、または時系列データの各時刻ごとに正則化を行って、前記ｋ空間データに対して画像再構成処理を行い前記出力対象画像を生成する。

図１は、第１の実施形態に係る磁気共鳴イメージング装置を示すブロック図。 図２は、第１の実施形態に係る磁気共鳴イメージング装置が行う画像生成処理の概略を説明するための図（１）。 図３は、第１の実施形態に係る磁気共鳴イメージング装置が行う画像生成処理の概略を説明するための図（２）。 図４は、第１の実施形態に係る磁気共鳴イメージング装置の行う処理手順を示すフローチャート（１）。 図５は、第１の実施形態に係る磁気共鳴イメージング装置の行う処理手順を示すフローチャート（２）。 図６は、第１の実施形態に係る磁気共鳴イメージング装置の行う処理手順を示すフローチャート（３）。 図７は、第１の実施形態に係る磁気共鳴イメージング装置の行う処理手順を示すフローチャート（４）。 図８は、第１の実施形態に係る磁気共鳴イメージング装置が行う領域重み係数の算出処理について説明した図（１）。 図９は、第１の実施形態に係る磁気共鳴イメージング装置が行う領域重み係数の算出処理について説明した図（２）。 図１０は、第１の実施形態に係る磁気共鳴イメージング装置が行う領域重み係数の算出処理について説明した図（３）。 図１１は、第１の実施形態に係る磁気共鳴イメージング装置が行う領域重み係数の算出処理について説明した図（４）。 図１２は、第１の実施形態に係る磁気共鳴イメージング装置が行う領域重み係数の算出処理について説明した図（５）。 図１３は、第１の実施形態に係る磁気共鳴イメージング装置が行う画像生成処理について説明した図（１）。 図１４は、第１の実施形態に係る磁気共鳴イメージング装置が行う画像生成処理について説明した図（２）。 図１５は、第１の実施形態に係る磁気共鳴イメージング装置が行う画像生成処理を従来法と比較した図。 図１６は、第２の実施形態に係る磁気共鳴イメージング装置を示すブロック図。

以下、本発明の実施形態について図面を参照しながら説明する。ここで、互いに同じ構成には共通の符号を付して、重複する説明は省略する。

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

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

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

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

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

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

受信回路１０９は、受信コイル１０８から出力される磁気共鳴信号に基づいて磁気共鳴データを生成する。具体的には、受信回路１０９は、受信コイル１０８から出力される磁気共鳴信号をデジタル変換することによって磁気共鳴データを生成する。また、受信回路１０９は、生成した磁気共鳴データをシーケンス制御回路１１０へ送信する。なお、受信回路１０９は、静磁場磁石１０１や傾斜磁場コイル１０２等を備える架台装置側に備えられていてもよい。ここで、第１の実施形態において、受信コイル１０８の各コイルエレメントから出力される磁気共鳴信号は、適宜分配合成されることで、チャネル等と呼ばれる単位で受信回路１０９に出力される。このため、磁気共鳴データは、受信回路１０９以降の後段の処理においてチャネル毎に取り扱われる。コイルエレメントの総数とチャネルの総数との関係は、同一の場合もあれば、コイルエレメントの総数に対してチャネルの総数が少ない場合、あるいは反対に、コイルエレメントの総数に対してチャネルの総数が多い場合もある。以下において、「チャネル毎」のように表記する場合、その処理が、コイルエレメント毎に行われてもよいし、あるいは、コイルエレメントが分配合成されたチャネル毎に行われてもよいことを示す。なお、分配合成のタイミングは、上述したタイミングに限られるものではない。磁気共鳴信号若しくは磁気共鳴データは、後述する処理回路１５０の画像生成機能１２２ｅを用いた処理の前までに、チャネル単位に分配合成されればよい。

シーケンス制御回路１１０は、計算機システム１２０から送信されるシーケンス情報に基づいて、傾斜磁場電源１０３、送信回路１０７及び受信回路１０９を駆動することによって、被検体Ｐの撮像を行う。ここで、シーケンス情報は、撮像を行うための手順を定義した情報である。シーケンス情報には、傾斜磁場電源１０３が傾斜磁場コイル１０２に供給する電源の強さや電源を供給するタイミング、送信回路１０７が送信コイル１０６に送信するＲＦ（Ｒａｄｉｏ Ｆｒｅｑｕｅｎｃｙ）パルスの強さやＲＦパルスを印加するタイミング、受信回路１０９が磁気共鳴信号を検出するタイミング等が定義される。

なお、シーケンス制御回路１１０は、傾斜磁場電源１０３、送信回路１０７及び受信回路１０９を駆動して被検体Ｐを撮像した結果、受信回路１０９から磁気共鳴データを受信すると、受信した磁気共鳴データを計算機システム１２０へ転送する。

計算機システム１２０は、磁気共鳴イメージング装置１００の全体制御や、データ収集、画像生成等を行い、処理回路１５０、記録回路１２３、入力装置１２４、ディスプレイ１２５を有する。処理回路１５０は、インタフェース機能１２１、感度マップ生成機能１２２ａ、推定機能１２２ｂ、アンフォールド機能１２２ｃ、収集機能１２２ｄ、画像生成機能１２２ｅ、制御機能１２６を有する。推定機能１２２ｂ、アンフォールド機能１２２ｃの具体的処理、また画像生成機能１２２ｅの詳細については後述する。

第１の実施形態では、インタフェース機能１２１、感度マップ生成機能１２２ａ、推定機能１２２ｂ、アンフォールド機能１２２ｃ、収集機能１２２ｄ、画像生成機能１２２ｅ、制御機能１２６にて行われる各処理機能は、コンピュータによって実行可能なプログラムの形態で記憶回路１２３へ記憶されている。処理回路１５０はプログラムを記憶回路１２３から読み出し、実行することで各プログラムに対応する機能を実現するプロセッサである。換言すると、各プログラムを読みだした状態の処理回路１５０は、図１の処理回路１５０内に示された各機能を有することになる。なお、図１においては単一の処理回路１５０にて、インタフェース機能１２１、感度マップ生成機能１２２ａ、推定機能１２２ｂ、アンフォールド機能１２２ｃ、収集機能１２２ｄ、画像生成機能１２２ｅ、制御機能１２６にて行われる処理機能が実現されるものとして説明したが、複数の独立したプロセッサを組み合わせて処理回路１５０を構成し、各プロセッサがプログラムを実行することにより機能を実現するものとしても構わない。

換言すると、上述のそれぞれの機能がプログラムとして構成され、１つの処理回路が各プログラムを実行する場合であってもよいし、特定の機能が専用の独立したプログラム実行回路に実装される場合であってもよい。なお、シーケンス制御回路１１０、処理回路１５０の有する推定機能１２２ｂ、収集機能１２２ｄ、画像生成機能１２２ｅは、それぞれシーケンス制御部、推定部、収集部、生成部の一例である。

上記説明において用いた「プロセッサ」という文言は、例えば、ＣＰＵ（Ｃｅｎｔｒａｌ Ｐｒｏｃｅｓｓｉｎｇ Ｕｎｉｔ）、ＧＰＵ（Ｇｒａｐｈｉｃａｌ Ｐｒｏｃｅｓｓｉｎｇ Ｕｎｉｔ）或いは、特定用途向け集積回路（Ａｐｐｌｉｃａｔｉｏｎ Ｓｐｅｃｉｆｉｃ Ｉｎｔｅｇｒａｔｅｄ Ｃｉｒｃｕｉｔ：ＡＳＩＣ）、プログラマブル論理デバイス（例えば、単純プログラマブル論理デバイス（Ｓｉｍｐｌｅ Ｐｒｏｇｒａｍｍａｂｌｅ Ｌｏｇｉｃ Ｄｅｖｉｃｅ：ＳＰＬＤ）、複合プログラマブル論理デバイス（Ｃｏｍｐｌｅｘ Ｐｒｏｇｒａｍｍａｂｌｅ Ｌｏｇｉｃ Ｄｅｖｉｃｅ：ＣＰＬＤ）、及びフィールドプログラマブルゲートアレイ（Ｆｉｅｌｄ Ｐｒｏｇｒａｍｍａｂｌｅ Ｇａｔｅ Ａｒｒａｙ：ＦＰＧＡ））等の回路を意味する。プロセッサは記憶回路１２３に保存されたプログラムを読み出し実行することで機能を実現する。

なお、記憶回路１２３にプログラムを保存する代わりに、プロセッサの回路内にプログラムを直接組み込むよう構成しても構わない。この場合、プロセッサは回路内に組み込まれたプログラムを読み出し実行することで機能を実現する。なお、寝台制御回路１０５、送信回路１０７、受信回路１０９等も同様に、上記のプロセッサ等の電子回路により構成される。

処理回路１５０は、インタフェース機能１２１により、シーケンス情報をシーケンス制御回路１１０へ送信し、シーケンス制御回路１１０から磁気共鳴データを受信する。また、インタフェース機能１２１を有する処理回路１５０は、磁気共鳴データを受信すると、受信した磁気共鳴データを記憶回路１２３に格納する。記憶回路１２３は、単一または複数チャネル分の時系列のｋ空間データを記憶する。

画像生成機能１２２ｅを有する処理回路１５０は、収集機能１２２ｄにより収集されたデータや、記憶回路１２３で記憶されたデータを用いて、画像の生成を行う。画像生成機能１２２ｅを有する処理回路１５０には、異なる収集方法で収集された２種類の時系列ｋ空間データが入力される。２種類の時系列ｋ空間データとは、血管描出用シーケンスであり、例えばＡＳＬ−ＭＲＡ（Ａｒｔｅｒｉａｌ Ｓｐｉｎ Ｌａｂｅｌｉｎｇ−Ｍａｇｎｅｔｉｃ Ｒｅｓｏｎａｎｃｅ Ａｎｇｉｏｇｒａｐｈｙ）で収集されるような特定領域の血液を励起させて標識化した画像（ｔａｇ画像）とそうでない画像（ｃｏｎｔｒｏｌ画像）を指す。以下の実施形態の説明では、２種類の時系列ｋ空間データがＡＳＬ−ＭＲＡにおけるｔａｇ画像とｃｏｎｔｒｏｌ画像の場合について述べるが，本実施形態はこれに限定されるものではなく、例えば、造影ＭＲＡ（Ｍａｇｎｅｔｉｃ Ｒｅｓｏｎａｎｃｅ Ａｎｇｉｏｇｒａｐｈｙ）、ＦｌｏｗＳｐｏｉｌｅｄ−ＦＢＩ（Ｆｒｅｓｈ Ｂｌｏｏｄ Ｉｍａｇｉｎｇ）、Ｐｈａｓｅ ｃｏｎｔｒａｓｔなどの方法にも同様に適用可能である。これらの磁気共鳴イメージングシーケンスでは、両者の差分をとることで所望の像を描出する。画像生成機能１２２ｅを有する処理回路１５０は、感度マップ生成機能１２２ａ、推定機能１２２ｂ、アンフォールド機能１２２ｃを用いて、所定の画像を生成する。なお、画像生成機能１２２ｅによって処理回路１５０が得た画像は、必要に応じてディスプレイ１２５で表示されたり、記憶回路１２３に送られて記憶されたりする。

処理回路１５０は、感度マップ生成機能１２２ａにより、準備スキャンまたは本スキャンで収集機能１２２ｄによって収集された各チャネルの収集データから各コイルの感度分布情報である感度マップを生成する。処理回路１５０は、推定機能１２２ｂにより、収集機能１２２ｄによって収集された各チャネルの時系列のｋ空間データから、ｘ空間の各画素の血管領域らしさを推定する。感度マップ生成機能１２２ａを有する処理回路１５０は、感度マップとして、ｔａｇ画像再構成用とｃｏｎｔｒｏｌ画像再構成用の２種類の感度マップを生成してもよいし、１種類の感度マップを生成し共通で使用してもよい。

ここで、以下の実施形態において、「ｘ空間」とは、水素分布画像空間（便宜上、本実施形態では水素原子以外を撮像の対象とする場合も含めて水素分布画像と呼ぶ）のことであり、「ｘ空間データ」とは、ｘ空間内における信号点の集合のことである。ｘ空間上の異なる信号点は、ｘ空間上の異なる位置の信号点と対応付けられる。例えば、脳の３次元の水素分布画像は、３次元のｘ空間データであり、心臓のある断面を撮像した２次元のｘ空間画像は、２次元のｘ空間データである。

また、処理回路１５０は、収集機能１２２ｄにより、空間方向のサンプリング位置を変化させながら、複数のチャネルについて、単時相あるいは時系列のｋ空間データを収集する。また、処理回路１５０は、例えば、収集機能１２２ｄにより、ｔａｇ画像収集用シーケンスとｃｏｎｔｒｏｌ画像収集用シーケンスの２種類のシーケンスによりデータを収集する。

処理回路１５０は、制御機能１２６により、磁気共鳴イメージング装置１００の全体制御を行う。具体的には、制御機能１２６を有する処理回路１５０は、入力装置１２４を介して操作者から入力される撮像条件に基づいてシーケンス情報を生成し、生成したシーケンス情報をシーケンス制御回路１１０に送信することによって撮像を制御する。また、処理回路１５０は、制御機能１２６により、撮像の結果としてシーケンス制御回路１１０から送られる磁気共鳴データを用いて行われる画像再構成処理を制御したり、ディスプレイ１２５による表示を制御したりする。例えば、制御機能１２６は、ＡＳＩＣ（Application Specific Integrated Circuit）、ＦＰＧＡ（Field Programmable Gate Array）等の集積回路、ＣＰＵ（Central Processing Unit）、ＭＰＵ（Micro Processing Unit）等の電子回路により実装される。

記憶回路１２３は、インタフェース機能１２１を有する処理回路１５０によって受信された磁気共鳴データや、処理回路１５０が画像生成機能１２２ｅによって生成した画像データ等を記憶する。例えば、記憶回路１２３は、ＲＡＭ（Random Access Memory）、フラッシュメモリ等の半導体メモリ素子、ハードディスク、光ディスク等である。

入力装置１２４は、操作者からの各種指示や情報入力を受け付ける。入力装置１２４は、例えば、マウスやトラックボール等のポインティングデバイス、あるいはキーボード等の入力デバイスである。ディスプレイ１２５は、処理回路１５０の有する制御機能１２６による制御のもと、画像データ等の各種の情報を表示する。ディスプレイ１２５は、例えば、液晶表示器等の表示デバイスである。

まず、実施形態に係る磁気共鳴イメージング装置に係る背景として、パラレルイメージングによる画像再構成及び正則化を中心に、簡単に説明する。

磁気共鳴イメージング装置１００は、核磁気共鳴現象を利用して被検体内部の情報を画像化するための装置である。磁気共鳴イメージング装置は、コイルを利用して、被検体内部にある特定の原子核（例えば、水素原子の原子核）からの核磁気共鳴信号をサンプリングすることで、ｋ空間データと呼ばれるデータを収集し、ｋ空間データにフーリエ変換を適用することで、磁気共鳴画像を得る。

核磁気共鳴信号は、１次元データとしてサンプリングされる。そこで磁気共鳴イメージング装置１００では、２次元あるいは３次元の磁気共鳴画像を得るために、ｋ空間上で１次元のサンプリングを繰り返し行い、磁気共鳴画像の生成に必要なデータを収集する。出力対象の磁気共鳴画像と同じ解像度（フルサンプリング）でｋ空間データがサンプリングされれば、得られたｋ空間データにフーリエ変換を適用することで、磁気共鳴画像の生成が可能である。

磁気共鳴イメージングにおいてサンプリングには時間がかかることが知られている。特に、時系列データの撮像においては、サンプリングに時間がかかるため、患者負担軽減のためにもサンプリング時間は短縮されることが望ましい。そのため、高速な撮像を実現するための様々な技術が研究、開発されている。そのような技術の１つとして、パラレルイメージング（ＰＩ（Parallel Imaging））と呼ばれる技術がある。パラレルイメージングでは、フルサンプリングに対して間引いた（アンダーサンプリングを行った）サンプリングパターンを用い、複数のコイルを用いてｋ空間データを収集する。アンダーサンプリングされたｋ空間データにそのままフーリエ変換を適用するとエイリアシング（画像の折り返し）が生じるが、パラレルイメージングでは、コイルの物理的な配置によって生じる感度の違いを利用することで、アンダーサンプリングによるエイリアシングのない磁気共鳴イメージング画像を生成する。

パラレルイメージングにおいてよく用いられる手法に、ＳＥＮＳＥ（sensitivity encoding）がある。パラレルイメージングは、大別して、ＳＥＮＳＥ、ｍＳＥＮＳＥ、ｋ−ｔ ＳＥＮＳＥ等のＳＥＮＳＥ系と呼ばれる手法と、ＧＲＡＰＰＡ（ＧｅｎｅＲａｉｚｅｄ Ａｕｔｏｃａｌｉｂｒａｔｉｎｇ Ｐａｒｔｉａｌｌｙ Ｐａｒａｌｌｅｌ Ａｃｑｕｉｓｉｔｉｏｎ）等のＳＭＡＳＨ（Ｓｉｍｕｌｔａｎａｕｅｏｕｓ Ａｃｑｕｉｓｉｔｉｏｎ ｏｆ Ｓｐａｔｉａｌ Ｈａｒｍｏｎｉｃｓ）系の手法がある。前者は、アレイコイルを用いてリファレンススキャンを行った後に、ｋ空間上のデータを間引いて収集した本スキャンの画像をフーリエ変換し、準備スキャン等により求めたコイルの感度マップを用いてエイリアシングのある画像を展開（アンフォールド）する方法である。後者は、ｋ空間上のデータをアンダーサンプリングしながら収集し、フーリエ変換前に、ｋ空間上で未収集のサンプル値を収集ｋ空間データを用いて生成した後、フーリエ変換を行いエイリアシングのない画像を得る方法である。本実施形態では、以下、再構成手法として、ＳＥＮＳＥを用いた例で説明するが、実施形態はこれに限られず、再構成手法として、他の再構成手法、例えばＳＥＮＳＥ系の他の手法を用いてもよい。

ＳＥＮＳＥの手順は次のとおりである。まず、シーケンス制御回路１１０は、準備スキャン等によって各コイルの感度分布情報を収集する。処理回路１５０は、収集機能１２２ｄにより、各コイルの感度分布情報をシーケンス制御回路１１０を通じて取得する。その後、処理回路１５０は、感度マップ生成機能１２２ａにより、各コイルの感度マップを生成する。

続いて、処理回路１５０は、画像生成機能１２２ｅにより、本スキャン（イメージングスキャン）で得られた各コイルのアンダーサンリングされたｋ空間データを用いてフーリエ変換を行い、各コイルの磁気共鳴画像を得る。

その後、処理回路１５０は、アンフォールド機能１２２ｃにより、「各コイルの磁気共鳴画像は、真の磁気共鳴画像に対してコイルの感度分布情報とエイリアシングのある各コイルにおける信号との積和演算を施して得られた画像である」という関係式に基づいて、真の磁気共鳴画像（出力対象画像）を推定する。

ＳＥＮＳＥでは、処理回路１５０は、アンフォールド機能１２２ｃにより、式（１）のエネルギー関数（目的関数）Ｅ（ｘ）を最小化することで、エイリアシングのない磁気共鳴画像x（出力対象画像）を再構成する。

ここで、cはコイル番号を、ｙ_ｃは、ｃ番目のコイルの、エイリアシングを含む磁気共鳴画像の信号値を表す。また、演算子Ａ_ｃは、ｃ番目のコイル感度をＳ_ｃ、フーリエ変換の操作をＦＴ、サンプリング操作をΓとすると、Ａ_ｃ＝Γ・ＦＴ・Ｓ_ｃと表現することができる演算子である。

ところで、式（１）よりわかるように、エネルギー関数Ｅ（ｘ）は、ｙ_ｃ-Ａ_ｃｘのＬ_２ノルムを、各コイルについて和をとったものになる。ノルムは一般に負の値を取り得ないから、エネルギー関数Ｅ（ｘ）は、すべてのコイル番号ｃについて、ｙ_ｃ-Ａ_ｃｘ＝０が成り立つならば、最小値０を取る。Ａ_ｃの逆演算子が一意的に定まるならば、ｙ_ｃ-Ａ_ｃｘ＝０をｘについて解くと、ｘ＝Ａ_ｃ ^−１ｙ_ｃとなる。従って、ある場合には、処理回路１５０は、アンフォールド機能１２２ｃにより、エイリアシングを含む磁気共鳴画像データｙ_ｃに演算子Ａ_ｃ ^−１を作用することにより、出力対象画像ｘを生成することができる。例えば、コイルの総数ｎがＲｅｄｕｃｔｉｏｎ ｆａｃｔｏｒ Ｒと等しい場合、多くの場合このケースに該当する。ここで、Ｒｅｄｕｃｔｉｏｎ ｆａｃｔｏｒは、フルサンプリングデータに対するサンプリングデータの間引き率を示す値である。例えば、サンプリングデータのデータサイズがフルサンプリングデータのデータサイズの１／２である場合、Ｒｅｄｕｃｔｉｏｎ ｆａｃｔｏｒは２となる。

しかしながら、Ｒｅｄｕｃｔｉｏｎ ｆａｃｔｏｒ Ｒがコイルの総数ｎと等しくないときには、この状況は一般には当てはまらない。

例えばＲ＞ｎのとき、式（１）のエネルギー関数Ｅ（ｘ）が０となる解ｘは一般に複数存在するから、出力対象画像ｘは、式（１）の最適化条件だけでは一意に定まらない。

逆にＲ＜ｎのときは、全てのコイル番号ｃについてｙ_ｃ-Ａ_ｃx＝０となるような解ｘは、多くの場合一般には存在しない。従って、この場合、式（１）のエネルギー関数Ｅ（ｘ）の値を、全ての出力対象画像ｘの候補について調べて、エネルギー関数Ｅ（ｘ）が最小になるような解が、出力対象画像ｘとなる。例えば、出力対象画像ｘが２５６画素であった場合、それぞれの画素に対応する２５６個の変数の値を変化させてエネルギー関数Ｅ（ｘ）を計算し、それらの値が最小になるような出力対象画像ｘを探索する。

しかしながら、出力対象画像ｘの画素数が多い場合、総当たり方式で全ての出力対象画像ｘの候補について調べることは計算時間の面で難しいケースも多い。その場合には、何らかの方法で出力対象画像ｘの候補を絞った上で、絞られた候補に対してのみエネルギー関数Ｅ（ｘ）の値を算出する。このようにして、収集されたｋ空間データから、出力対象画像を復元する逆問題が、解かれることになる。

しかし、逆問題としての画像再構成処理は、最適化問題の解が入力データの小さな変化に敏感に依存するという意味で、しばしばいわゆる不良設定問題となる。この場合、式（１）のエネルギー関数Ｅ（ｘ）を直接法で最適化するのは、有効な解法になりにくい。例えば、最小二乗法により逆問題を解く場合、解が一意に定まらず発散する可能性がある。

このような場合、解の安定性を高めるため、目的関数に対して所定の追加の項をペナルティ項として導入する操作、すなわち正則化が、しばしば行われる。具体的には、例えば、処理回路１５０は、画像生成機能１２２ｅにより、解の安定性を高めるため、式（２）のような正則化項が導入する。

ここで、λは正則化項の強さを制御するための正則化パラメータ（未定乗数）である。式（２）の場合では、処理回路１５０は、画像生成機能１２２ｅにより、出力対象画像に対するＬ_２ノルムを正則化項として、目的関数（式（１）のエネルギー関数Ｅ（ｘ））に加える。処理回路１５０は、画像生成機能１２２ｅにより、目的関数と正則化項との和（式（２）のエネルギー関数Ｅ（ｘ））が最小となるように、出力対象画像ｘを定める。

ここで、正則化処理を行う理由としては、例えば以下のような理由が挙げられる。まず、当初の最適化問題、すなわち式（１）のエネルギー関数Ｅ（ｘ）の最適化問題を解くにあたっては、現実的な計算時間では全ての出力対象画像ｘの候補について調べることができないので、近似や付加的な過程をおくことより、少ない数の出力対象画像ｘについてのみエネルギー関数Ｅ（ｘ）を算出することになる。従って、得られた出力対象画像ｘは、真のエネルギー関数Ｅ（ｘ）の最小値にならず、最小値に近似する値になる。しかしながら、医用画像など多くの目的においては、出力対象画像ｘの大まかな傾向が正しければ、細部の微小な信号値の違いは問題にならないことも多い。この場合、目的関数（式（１）のエネルギーＥ（ｘ））の真の最小値を与える出力対象画像ｘが得られる必要はなく、近似的な最小値を与える出力対象画像ｘを得られれば十分である。従って、近似的な最小値を与える出力対象画像ｘの複数の候補、言い換えれば、おおむね妥当な値を持つ出力対象画像ｘの複数の候補中から、真の出力対象画像ｘが満たすべき条件に適合する出力対象画像ｘを選択してもよい。

すなわち、正則化項とは、真の出力対象画像が満たすべき条件として期待する条件に適合しているかどうかを判定するための項である。換言すると、処理回路１５０は、画像生成機能１２２ｅにより、真の出力対象画像らしい画像に対して低い値を取り、真の出力対象画像らしくない画像に対してペナルティとして高い値を取る項を正則化項として目的関数に加え、目的関数と正則化項との和を最小化する。このことにより、目的関数が０に近いが真の出力対象画像らしくない画像が排除され、目的関数が０に近く、かつ真の出力対象画像らしい画像を得ることができる。

正則化項としては、スパースな解が得られるＬ_０ノルム、Ｌ_１ノルムや、計算が簡便なＬ_２ノルムが用いられる。出力画像に対するＬ_２ノルムは、例えば出力対象画像全体の分散を抑制するためにしばしば用いられる。

しかし、出力対象の画像に対して単純に所定のノルムを用いて正則化を行った場合、ペナルティ項として導入された正則化項の影響で、再構成後の出力対象画像の画質が、かえって低下しまうことがある。すなわち、正則化項は、出力対象画像全体の分散を抑制する項であるため、正則化パラメータを強くすると、再構成後の出力対象画像の画質に影響が出てしまう。

実施形態に係る磁気共鳴イメージング装置１００は、かかる背景に基づくものである。すなわち、第１の実施形態では、次の式（３）で表されるように、領域重み係数Ｗ（観察対象領域らしさを表す重み係数）を出力対象画像ｘに掛け合わせた項に関するノルムを正則化項として導入する。

式（３）において、領域重み係数Ｗは、後ほど詳述するが、出力対象画像の各画素の画素位置に対して定義された係数である。領域重み係数Ｗは、例えば、観察対象に含まれそうな画素に対する重みが、観察対象に含まれなさそうな画素に対する重みよりも小さくなるように、定められる。これにより、処理回路１５０は、出力対象画像の安定性を向上さて、かつ背景領域のアーチファクトを抑制することができる。

かかる点について、図２及び図３を用いて説明する。図２及び図３は、第１の実施形態に係る磁気共鳴イメージング装置が行う画像生成処理の概略を説明するための図である。なお、図２及び図３では、観察対象として、血管の画像を描出する場合（ＡＳＬ−ＭＲＡ）を例にとって説明する。

図２の左の図（正方形１０）は、出力対象画像の信号値（式（２）及び式（３）のｘに対応する量）を表す模式図である。縦長の長方形１３及び横長の長方形１５は、血管領域を表す。また、それ以外の部分は、背景信号の領域を表す。ＡＳＬ−ＭＲＡにおいては、血管領域以外の領域は背景信号になっているので、画像内に占める背景信号の面積が、比較的大きくなっている。

図２の中央の図（正方形１１）は、領域重み係数Ｗの値を表す模式図である。領域重み係数Ｗの詳細な算出方法については、後述する。正方形１１において、例えば、図の「１」と書かれている部分は、領域重み係数Ｗが１である領域を示している。また、図の「０」と書かれている部分は、領域重み係数Ｗが０である領域を示している。領域重み係数Ｗが１である領域は、観察対象領域でないと判断されている領域である。また、領域重み係数Ｗが０である領域は、観察対象領域であると判断されている領域である。

図２の左の図と中央の図とを比較してみるとわかるように、実際の血管領域（縦長の長方形１３及び横長の長方形１５）と、観察対象領域であると判断されている領域（領域重み係数Ｗが０である領域）は、必ずしも一致しなくても良い。また、観察対象領域であるか否かの判定は、必ずしも厳密に正確なものである必要もない。すなわち、この例でいくと、領域重み係数Ｗは、実際の血管が存在すると期待されるおおまかな位置を示す程度のものであればよい。

図２の右の図（正方形１２）は、実施形態に係る正則化項（例えば、式（３）の第２項）の値を表す模式図である。式（３）のｘに対応する量は正方形１０に示されており、また、式（３）の領域重み係数Ｗに対応する量は正方形１１に示されている。この二つを掛け算すると、縦長の長方形１３及び横長の長方形１５で示される血管領域は、領域重み係数Ｗの値「０」にマスクされて、正則化項への寄与は０になる。これに対して、血管領域以外の領域は、領域重み係数Ｗの値「０」にマスクされて、正則化項への寄与が０になるか、あるいは、領域重み係数Ｗの値「１」が掛け合わされ、もともとの信号値がそのまま正則化項への寄与となる。図２の右図をみるとわかるように、実施形態に係る正則化項（式（３）の第２項）は、血管領域以外の領域の信号値に対するペナルティ関数となっている。

これに対して、領域重み係数Ｗを伴わない正則化項（例えば、式（２）の第２項）の場合は、図２の左の図を参照するとわかるように、血管領域以外の領域の信号値に起因する寄与の他に、血管領域の信号値に起因する寄与が混じってしまう。血管領域は、全体では小さな面積を占めるにすぎないものの、信号の強度は高いことから、血管領域の信号値に起因する寄与が、正則化項にペナルティ関数として生じてしまう。この結果、画質が劣化する。逆に、画像生成機能１２２ｅを有する処理回路１５０は、領域重み係数Ｗを伴う正則化項を用いることで、出力対象画像の画質を向上させることができる。

次に、図３を用いて、画像生成機能１２２ｅを有する処理回路１５０が領域重み係数Ｗを伴う正則化項を用いることで、出力対象画像の画質を向上させることができる点について別の角度から説明する。なお、図３においては、説明の便宜のため、正則化項として、Ｌ_１ノルムである場合について説明するが、他のノルムを取る場合でも、同様の説明が成り立つ。

図３の例では、画素数が４×４＝１６個の場合について説明する。今、正方形２０において、斜線で示されている画素（ピクセル）が血管領域であると仮定する。

先に説明したように、正則化を行わない場合（式（１）の場合）、出力対象画像ｘの探索空間は現実的な計算時間スケールでは十分広くとることが難しいため、出力対象画像ｘの近似解としては、適切な解の他に、不適切な解も現れる。図３上段の左の図（正方形２１）は、得られた適切な解、図３上段の中央の図及び右の図（正方形２２及び正方形２３）は、得られた不適切な解の例を表している。

次に、領域重み係数Ｗを伴わない正則化を行った場合について説明する。図３の上段（正方形２１、２２、２３）は、領域重み係数Ｗを伴わない正則化を行った場合（例えば、式（２）でＬ_２ノルムをＬ_１ノルムに変えた場合）について説明している。この場合、すべての画素値を足すと、正方形２１、正方形２２、正方形２３はそれぞれ、画素値合計（Ｌ_１ノルム）は４４となる。従って、正方形２１、正方形２２、正方形２３は、正則化項という観点では等価であり、これらの解を区別できない。従って、処理回路１５０は、不適切な解を排除することができない。

次に、領域重み係数Ｗを用いた正則化を行った場合について説明する。図３の下段（正方形２４、２５、２６）は、領域重み係数Ｗを伴う正則化を行った場合（例えば、式（３）でＬ_２ノルムをＬ_１ノルムに変えた場合）について説明している。この場合、例えば正方形２１、正方形２２、正方形２３に対して、正方形２０の斜線で示されている領域に対して領域重み係数Ｗによるマスクが行われて、正則化項への寄与は、それぞれ正方形２４、正方形２５、正方形２６のようになる。この場合、各場合において正則化項を計算すると、図３の下段の左の図（正方形２４）の場合は、「９」、図３の下段の中央の図（正方形２５）の場合は「２３」、図３の下段の右側の図（正方形２６）の場合は、「２９」となる。従って、領域重み係数Ｗを用いた正則化を行った場合、図３の中央の図及び右の図の場合はペナルティ関数が大きくなることから解の候補としては除外され、図３の左の図の場合が最終的な出力画像として選択される。これにより、画質が向上する。

以上まとめると、処理回路１５０は、推定機能１２２ｂにより、第１段階として、後段の最適化に比較して簡便な方法を用いて領域重み係数Ｗを算出する。次に、処理回路１５０は、画像生成機能１２２ｅにより、第２段階として、算出された領域重み係数Ｗを用いた正則化項により正則化を行って、出力対象画像を生成する。これにより、処理回路１５０は、出力対象画像の安定性を向上させ、かつ背景領域のアーチファクトを抑制することができる。

続いて、第１の実施形態に係る磁気共鳴イメージング装置の行う処理手順について、必要に応じて適宜図８〜１５を参照しながら、図４〜７を用いて説明する。

なお、以下のフローチャートでは、例えばＡＳＬ−ＭＲＡで収集されるような特定領域の血液を励起させて標識化した画像（ｔａｇ画像）とそうでない画像（ｃｏｎｔｒｏｌ画像）に係る２種類の時系列ｋ空間データを、アンダーサンプリングを行いながら収集するシーケンスを実行し、パラレルイメージングにより画像再構成を行い血管の描出をする場合で説明するが、実施形態はこれに限られない。

図４〜７は、第１の実施形態に係る磁気共鳴イメージング装置の行う処理手順を示すフローチャートである。図４において、まず、処理回路１５０が、制御機能１２６により、操作者から撮像条件の入力を受け付ける（ステップＳ１０１）。続いて、処理回路１５０は、収集機能１２２ｄにより、操作者から入力された撮像条件に基づいてシーケンス情報を生成し、生成したシーケンス情報をシーケンス制御回路１１０に送信することで、準備スキャンの実行を制御する。シーケンス制御回路１１０は、準備スキャンを行う。（ステップＳ１０２）。ここで、準備スキャンには、例えば、位置決め用の画像を収集するスキャンや、静磁場の不均一性を補正するシミングスキャン、感度分布情報を収集するスキャン等が含まれる。

準備スキャンの終了後、続いて、処理回路１５０は、収集機能１２２ｄにより、操作者から入力された撮像条件に基づいてシーケンス情報を生成し、生成したシーケンス情報をシーケンス制御回路１１０に送信することで、出力対象画像（例えば、診断用に出力される画像）を収集するイメージングスキャンの実行を制御する（ステップＳ１０３）。シーケンス制御回路１１０は、シーケンス情報を処理回路１５０から受信すると、アンダーサンプリングを行いながらパルスシーケンスを実行する。具体的には、シーケンス制御回路１１０は、第１のパルスシーケンスと、第２のパルスシーケンスを実行する。ここで、第１のパルスシーケンスは、撮像領域内の標識化領域に対して標識化パルスを印加するパルスシーケンス、すなわちｔａｇ画像に係るパルスシーケンスである。また、第２のパルスシーケンスは、標識化パルスを印加しないパルスシーケンスである。この時、シーケンス制御回路１１０は、第１のパルスシーケンスと第２のパルスシーケンスとのうち少なくとも一方をアンダーサンプリングを行いながら実行する。処理回路１５０は、収集機能１２２ｄにより、第１のパルスシーケンスを用いて第１のｋ空間データを収集し、第２のパルスシーケンスを用いて第２のｋ空間データを収集する。ここで、第１のｋ空間データはｔａｇ画像生成用のｋ空間データである。また、第２のｋ空間データは、ｃｏｎｔｒｏｌ画像生成用のｋ空間データである。また、この時収集されるデータは時系列データである。すなわち、処理回路１５０は、収集機能１２２ｄにより、ｋ空間データとして、複数チャネルで時系列ｋ空間データを収集する。

第１の実施形態においては、シーケンス制御回路１１０は、Ｒｅｄｕｃｔｉｏｎ ｆａｃｔｏｒ Ｒ（フルサンプリングデータに対するサンプリングデータの間引き率）に従って、データを間引きながら時系列ｋ空間データを収集する。すなわち、シーケンス制御回路１１０が収集する時系列ｋ空間データの個数は、フルサンプリングに比べて１／Ｒとなる。具体的には、シーケンス制御回路１１０は、ｋ空間データとして、ＲＯ（Ｒｅａｄ ｏｕｔ）方向の信号点の個数×ＰＥ（Ｐｈａｓｅ Ｅｎｃｏｄｅ）方向の信号点の個数×ＳＥ（ＳｌｉｃｅＥｎｃｏｄｅ）方向の信号点の個数×時間方向のフレーム数÷Ｒｅｄｕｃｔｉｏｎ ｆａｃｔｏｒ Ｒ個×チャネル数、で求められる個数の信号点を収集する。

そして、処理回路１５０は、画像生成機能１２２ｅにより、ステップＳ１０３において収集された時系列のｋ空間データ、または、ステップＳ１０３で収集され記憶回路１２３に格納された複数チャネル分の時系列のｋ空間データを用いて、出力対象の画像の生成を行う（ステップＳ１０４）。生成された出力対象の画像は、必要に応じて、記憶回路１２３に格納されたり、ディスプレイ１２５に表示されたりする。

図５を用いて、第１の実施形態におけるステップＳ１０４の処理を詳述する。

まず、処理回路１５０は、感度マップ生成機能１２２ａにより、各コイルの感度分布情報である感度マップを生成する（ステップＳ２０１）。感度マップ生成機能１２２ａを備える処理回路１５０は、準備スキャンまたは本スキャンにおいて収集したデータから各コイルの感度マップを生成する。処理回路１５０は、続いて、アンフォールド機能１２２ｃにより処理を行う。

処理回路１５０は、推定機能１２２ｂにより、ｔａｇ画像生成用データとｃｏｎｔｒｏｌ画像生成用データから領域重み係数Ｗを算出する（ステップＳ２０２）。ここで、ｔａｇ画像生成用データ及びｃｏｎｔｒｏｌ画像生成用データそれぞれは、４Ｄ（３Ｄ＋ｔ（時系列））のｋ空間データである。この４Ｄのｋ空間データは、各々がチャネルの個数分だけ用意されている。領域重み係数Ｗは、ＡＳＬ−ＭＲＡの場合、ｘ空間の各画素における血管領域らしさを表す量であり、例えば３Ｄのｘ空間の各画素の血管領域らしさを示す０から１の値である。すなわち、推定機能１２２ｂを有する処理回路１５０は、ｔａｇ画像生成用データとｃｏｎｔｒｏｌ画像生成用データの２種類の４Ｄ（３Ｄ＋ｔ（時系列））ｋ空間データから、３Ｄｘ空間の各画素の血管領域らしさを示す領域重み係数Ｗを算出する。換言すると、処理回路１５０は、推定機能１２２ｂにより、ｋ空間データを用いて生成されたデータに基づいて、出力対象画像の各画素の画素位置が観察対象領域に含まれているか否かを推定する。ここで、ｋ空間データを用いて生成されたデータとは、より具体的には、第１のｋ空間データ（ｔａｇ画像生成用データ）と第２のｋ空間データ（ｃｏｎｔｒｏｌ画像生成用データ）とを用いて生成されたデータである。

図６はステップＳ２０２の詳細な処理手順を示すフローチャートである。まず、処理回路１５０は、推定機能１２２ｂにより、４Ｄのｋ空間データであるｔａｇ画像用データ（第１のｋ空間データ）とｃｏｎｔｒｏｌ画像生成用データ（第２のｋ空間データ）との間で差分処理を行って、４Ｄのｋ空間の差分画像データ（第１の差分ｋ空間データ）を生成する（ステップＳ２０２ａ）。ここで、第１のｋ空間データ及び第２のｋ空間データはいずれも複素数データであることから、ステップＳ２０２ａの処理は複素差分処理となる。また、第１のｋ空間データ及び第２のｋ空間データは、例えば、ｋ空間のＲＯ方向、ｋ空間のＰＥ方向、ｋ空間のＳＥ方向、時間、がそれぞれ軸となる時系列データになっており、各々がチャネルの個数分だけ用意されている。

また、第１のｋ空間データ（ｔａｇ画像用データ）と第２のｋ空間データ（ｃｏｎｔｒｏｌ画像用データ）は、いずれも同一のサンプリングパターンで間引き収集されたデータ、でも良いし、異なるサンプリングパターンで間引き収集した後に同一のサンプリングパターンになるように変化したデータである。

続いて、処理回路１５０は、推定機能１２２ｂにより、生成された第１の差分ｋ空間データに対して時間平均処理を行って、平均差分ｋ空間データを生成する（ステップＳ２０２ｂ）。第１のｋ空間データと第２のｋ空間データのサンプリングパターンが、時間軸方向に射影した場合にｋ空間を全て埋めるパターンの場合、平均差分ｋ空間データはフルサンプリング相当のｋ空間データとなる。この時点で、平均差分ｋ空間データは時系列データではなくなり、ＲＯ×ＰＥ×ＳＥサイズのｋ空間データがチャネルの個数分あるデータとなる。時間平均処理は全時刻のｋ空間データを対象に行ってもよいし、ある範囲の時刻のｋ空間データのみを対象に行ってもよい。また、平均処理ではなく、加算処理であってもよい。

続いて、処理回路１５０は、推定機能１２２ｂにより、ステップＳ２０２ｂで算出した平均差分ｋ空間データに対して、３次元（逆）フーリエ変換を行って、平均差分ｘ空間データを生成する（ステップＳ２０２ｃ）。すなわち、この平均差分ｘ空間データは、時系列ｋ空間データを用いて生成された時系列ではないデータである。ステップＳ２０２ｂで算出した平均差分ｋ空間データがフルサンプリング相当のデータであるため、ステップＳ２０２ｃで得られる平均差分ｘ空間データには折り返しがわずかである。また、ステップＳ２０２ａでｔａｇ画像用データとｃｏｎｔｒｏｌ画像用データの差分を算出しているため、ステップＳ２０２ｃで描出される平均差分ｘ空間データは血液信号となる。

本実施形態で示したステップＳ２０２ａ〜ステップＳ２０２ｃの処理順序は一例であり、まずステップＳ２０２ｂでｔａｇ画像用データとｃｏｎｔｒｏｌ画像用データの時間平均処理を行い、その後ステップＳ２０２ａの複素差分処理を行った後に、ステップＳ２０２ｃの処理を行ってもよいし、まずステップＳ２０２ｂで時間平均処理を行い、その後ステップＳ２０２ｃで画像変換処理を行った後に、ステップＳ２０２ａの差分処理を行ってもよい。

続いて、処理回路１５０は、推定機能１２２ｂにより、ステップＳ２０２ｃで算出した全チャネルの平均差分ｘ空間データを用いてｘ空間の各画素の血管領域らしさを判定する（ステップＳ２０２ｄ）。すなわち、処理回路１５０は、推定機能１２２ｂにより、時系列ｋ空間データを用いて生成された時間軸を持たない空間データに基づいて、前記出力対象画像の各画素の前記画素位置が前記観察対象領域に含まれているか否かを推定する。

図７はステップＳ２０２ｄの詳細な処理手順を示すフローチャートである。推定機能１２２ｂを有する処理回路１５０は、ステップＳ４０１及びステップＳ４０２の処理を、各チャネルごとに行う。すなわち、推定機能１２２ｂを有する処理回路１５０は、n_{max_ch}をチャネル数として、それぞれのチャネルについて、ステップＳ４０１及びステップＳ４０２の処理を繰り返す。ステップＳ４０１及びステップＳ４０２の処理は、各チャネルについて、直列的に行われても良いし、並列的に行われても良い。なお、ここではステップＳ４０１及びステップＳ４０２をチャネルごとに行う場合について述べるが、予め複数チャネルのデータを統合したのちにステップＳ４０１及びステップＳ４０２の処理を行ってもよい。

まず、ｉ番目のチャネルについて、処理回路１５０は、推定機能１２２ｂにより、ステップＳ２０２ｃで算出した平均差分ｘ空間データからノイズ分布を推定する。具体的には、ｉ番目のチャネルについて、処理回路１５０は、推定機能１２２ｂにより、平均差分ｘ空間データで信号値が低い領域である空気領域を特定し、特定した平均差分ｘ空間データに基づいて、画素集合の信号値の平均値μ_ｉと標準偏差σ_ｉを求める（ステップＳ４０１）。空気領域で生じるノイズは複素空間でガウス分布となることが知られている。空気領域の特定は、予め設定した任意の領域であってもよいし、準備スキャンやイメージングスキャンで得られた収集データを利用してもよい。

続いて、ｉ番目のチャネルについて、処理回路１５０は、推定機能１２２ｂにより、各画素の画素位置が観察対象領域に含まれている確率Ｐ_ｉの推定値（血管領域らしさ）を算出する（ステップＳ４０２）。最も典型的な例として、処理回路１５０は推定機能１２２ｂにより、例えばステップＳ４０１で算出したノイズ分布を基準とし、ステップＳ２０２ｃで算出した平均差分ｘ空間データの各画素において、信号値が基準より小さい値を取る場合には、背景領域と判定し、基準より大きい値を取る場合には被写体領域（観察対象領域）と判定する。このようにして、平均差分ｘ空間データの各画素の信号値が背景領域か被写体領域（観察対象領域）かを判定する。なお、判定結果は、観察対象領域か否かといった２値の判定に限られない。また、別の例として、処理回路１５０は、推定機能１２２ｂにより、ステップＳ４０１において、画像全体から、画素集合の信号値の平均値μ_ｉと標準偏差σ_ｉを求め、ステップＳ４０２において、求めた平均値μ_ｉ及び標準偏差σ_ｉに基づいて、各画素の画素位置が観察対象領域に含まれている確率Ｐ_ｉの推定値を算出してもよい。

図８及び図９に、ステップＳ４０２の処理の典型的な具体例が示されている。図８及び図９は、第１の実施形態に係る磁気共鳴イメージング装置が行う領域重み係数Ｗ_i（及び確率Ｐi）の算出処理について説明した図である。なお、図８〜１１において、説明を簡略化するために、平均差分ｘ空間データ全体から、画素集合の信号値の平均値μ_ｉと標準偏差σ_ｉを求める例で説明している。しかしながら、この説明は、実施形態における処理方法である、平均差分ｘ空間データのうちの空気領域のデータから、画素集合の信号値の平均値μ_ｉと標準偏差σ_ｉを求め、それを基に領域重み係数Ｗ_ｉを算出する場合についても同様の考え方で適用できる。また、図８〜１１においては、説明を簡略化するために、画像データが実数の場合で説明しているが、実際の画像データは複素数のデータである。例えば、複素数の絶対値に対して同様の考え方を適用する、又は実部、虚部それぞれに対して同様の考え方を適用するなどして、複素数である実際の画像データに対しても、同様な考え方を適用できる。

図８では、まず、ｉ番目のチャネルの領域重み係数Ｗ_ｉ（または、観察対象領域に含まれている確率Ｐ_ｉ）の算出方法として、０か１かの２値の判定を行う場合について説明する。

図８に示す例では、推定機能１２２ｂを有する処理回路１５０は、各画素について、画素値が、画素集合の信号値の平均値μ_ｉと比較して、所定の許容標準偏差以内に収まっているかどうかで、観察対象領域か否かを判定する。

例えば、図８において、ｉ番目のチャネルの画素集合の平均値μ_ｉ＝２０、画素集合の標準偏差σ_ｉ＝５、許容標準偏差２σであるので、観察対象領域に含まれるための画素値ｙの基準は２０−５×２≦ｙ≦２０＋５×２となる。従って、処理回路１５０は、推定機能１２２ｂにより、１０≦ｙ≦３０の範囲に含まれる画素が、観察対象領域であると判定する。従って、処理回路１５０は、画素値が２３の画素１を、「○」、すなわち観察対象領域にあると判定する。また、処理回路１５０は、画素値が８の画素２、及び画素値が３２の画素３を、「×」、すなわち観察対象領域にないと判定する。今、各画素の画素位置が観察対象である確率Ｐ_ｉは０か１かの２値化されているので、推定機能１２２ｂを有する処理回路１５０は、画素１に対しては、Ｐ_ｉ＝１、画素２及び画素３に対しては、Ｐ_ｉ＝０と算出する。

続いて、推定機能１２２ｂを有する処理回路１５０は、各画素について、領域重み係数Ｗ_ｉを、各画素の画素位置が観察対象である確率Ｐ_ｉに基づいて定める。一例として、処理回路１５０は、推定機能１２２ｂにより、領域重み係数Ｗ_ｉを、Ｗ_ｉ＝１−Ｐ_ｉに基づいて定める。換言すると、正則化項に用いられる領域重み係数Ｗ_ｉは、各画素の画素位置が観察対象領域に含まれている場合には０であり、各画素の画素位置が観察対象領域に含まれていない場合には１である。例えば、処理回路１５０は、推定機能１２２ｂにより、画素１に対しては、Ｗ_ｉ＝０、画素２及び画素３に対しては、Ｗ_ｉ＝１と算出する。

図９の例では、図８と同様であるが、領域重み係数Ｗ_ｉ及び観察対象領域に含まれている確率Ｐ_ｉの算出方法として、２値の判定ではなく、連続量を用いて算出する場合について説明する。

図９に示す例では、推定機能１２２ｂを有する処理回路１５０は、各画素について、画素値が、画素集合の信号値の平均値μ_ｉと比較して、どれだけずれているかを表す量ｘを算出する。例えば、図９において、ｉ番目のチャネルの画素集合の平均値μ_ｉ＝２０、画素集合の標準偏差σ_ｉ＝５であるので、画素値が２３の画素１は、平均から、ｘ＝（２３−２０）÷５＝＋０．６σ_ｉだけずれていることになる。同様に、画素２、画素３の平均からのずれは、ｘ＝−２．４σ_ｉ、ｘ＝＋２．４σ_ｉとなる。続いて、推定機能１２２ｂを有する処理回路１５０は、各画素の画素位置が観察対象である確率Ｐ_ｉを、中心極限定理が成り立つと仮定し、平均値から離れるにしたがって指数関数的に減少すると仮定して、ｘ＝０のときとの相対値として、Ｐ_ｉ＝ｅｘｐ（−ａｂｓ（ｘ）／σ_ｉ）で算出する。

続いて、推定機能１２２ｂを有する処理回路１５０は、各画素について、領域重み係数Ｗ_ｉを、各画素の画素位置が観察対象である確率Ｐ_ｉに基づいて定める。換言すると、正則化項に用いられる重みは、各画素の画素位置が観察対象に含まれる確率の推定値に基づいて算出された値となる。一例として、処理回路１５０は、推定機能１２２ｂにより、領域重み係数Ｗ_ｉを、Ｗ_ｉ＝１−Ｐ_ｉに基づいて算出する。この場合、Ｗ_ｉ＝１−ｅｘｐ（−ａｂｓ（ｘ）／σ_ｉ）となる。すなわち、画素１については、Ｗ_ｉ＝０．４５１、画素２については、Ｗ_ｉ＝０．９０９、画素３については、Ｗ_ｉ＝０．９０９となる。

図７に戻り、処理回路１５０は、推定機能１２２ｂにより、ステップＳ４０２で得られた各チャネルの判定結果を統合する（ステップＳ４０３）。換言すると、正則化項に用いられる重みは、ｋ空間データの収集単位であるチャネルごとに算出された値を統合することにより算出された値となる。統合方法は画素ごとに全チャネルの単純加算でもよいし、重み付け加算でもよい。また、その他の統合方法を用いてもよい。図１０及び図１１に、各チャネルの判定結果の統合の具体例を示す。図１０及び図１１は、第１の実施形態に係る磁気共鳴イメージング装置が行う領域重み係数の算出処理について説明した図である。

図１０には、領域重み係数Ｗ_ｉが２値化されている場合についての領域重み係数の統合について示されている。図８で示されているのと同様に、処理回路１５０は、推定機能１２２ｂにより、各チャネルでの領域重み係数Ｗ_ｉを、それぞれのチャネルでの画素集合の平均値μ_ｉ、標準偏差σ_ｉ及び許容標準偏差に基づいて定める。処理回路１５０は、推定機能１２２ｂにより総合判定を行い、一例として、各チャネルの領域重み係数Ｗ_ｉの平均値を、領域重み係数Ｗの値として算出する。例えば、画素１に対しては、１番目のチャネルの領域重み係数Ｗ_１＝０であり、２番目のチャネルの領域重み係数Ｗ_２＝０であり、３番目のチャネルの領域重み係数Ｗ_３＝１であるから、処理回路１５０は、推定機能１２２ｂにより、画素１に対する領域重み係数ＷをＷ＝（Ｗ_１＋Ｗ_２＋Ｗ_３）／３＝１／３と算出する。同様に、処理回路１５０は、推定機能１２２ｂにより、画素２、画素３に対する領域重み係数Ｗを、それぞれＷ＝１、Ｗ＝２／３と算出する。

なお、実施形態はこれに限られない。例えば、処理回路１５０は、推定機能１２２ｂにより、領域重み係数Ｗの値と、各チャネルでの領域重み係数Ｗ_ｉの平均ではなく和、最大値、最小値などとしてもよく、また、各チャネルでの領域重み係数Ｗ_ｉを、論理演算「ＡＮＤ」や「ＯＲ」で結合して領域重み係数Ｗの値を算出してもよい。また、処理回路１５０は、推定機能１２２ｂにより、チャネルごとの感度の値を反映した重みで、各チャネルの領域重み係数Ｗ_ｉの加重線形和を取っても良い。

図１１には、領域重み係数Ｗ_ｉの算出方法として、連続量を用いて算出する場合が示されている。処理回路１５０は、推定機能１２２ｂにより、各チャネルの領域重み係数Ｗ_ｉを、図９で計算したのと同様の方法で算出する。続いて、処理回路１５０は、推定機能１２２ｂにより、一例として、各チャネルの領域重み係数Ｗ_ｉの二乗平均平方根（root mean square）を、領域重み係数Ｗの値として算出する。また、別の例として、続いて、処理回路１５０は、推定機能１２２ｂにより、一例として、各チャネルの領域重み係数Ｗ_ｉの最大値や最小値を、領域重み係数Ｗの値として算出しても良い。

図１２には、このようにして推定機能１２２ｂを有する処理回路１５０が算出した領域重み係数Ｗの一例が示されている。図１２は、第１の実施形態に係る磁気共鳴イメージング装置が行う領域重み係数の算出処理について説明した図である。パネル３０に、複数の断面での領域重み係数の値が示されている。画像３１、画像３２、画像３３は、それぞれ、Ｘ−Ｙ平面、Ｚ−Ｙ平面、Ｘ−Ｚ平面での領域重み係数Ｗを表している。領域７０、領域７１、領域７２、領域７３は、領域重み係数Ｗの値が低値になっており、観察対象領域（血管領域）になっている。また、その他の領域は、領域重み係数Ｗの値が高値になっており、観察対象領域以外の領域になっている。

図６に戻り、処理回路１５０は、推定機能１２２ｂにより、ステップＳ２０２ｄにおける算出結果に基づいて、アンフォールド処理において使用する領域重み係数を最終的に算出する（ステップＳ２０２ｅ）。具体的には、推定機能１２２ｂを有する処理回路１５０は、ステップＳ２０２ｅにおいて、例えばシグモイド関数や指数関数を用いて変換処理を行ってもよい。領域重み係数Ｗは背景領域ほど１に近い値をとり、血管領域ほど０に近い値をとる。

図５に戻り、推定機能１２２ｂを有する処理回路１５０がステップＳ２０２で、ｔａｇ画像及びｃｏｎｔｒｏｌ画像それぞれについて領域重み係数Ｗの算出を終えると、処理はアンフォールド機能１２２ｃ及び画像生成機能１２２ｅに受け継がれる。すなわち、処理回路１５０は、アンフォールド機能１２２ｃ及び画像生成機能１２２ｅにより以降の処理を行う。

続いて、処理回路１５０は、アンフォールド機能１２２ｃ及び画像生成機能１２２ｅにより、ｔａｇ画像とｃｏｎｔｒｏｌ画像それぞれについて、再構成後のｘ空間データを生成し、生成したｘ空間データを基に、出力対象の画像を生成する（ステップＳ２０３）。具体的には、アンフォールド機能１２２ｃ及び画像生成機能１２２ｅを有する処理回路１５０は、ステップＳ１０３において収集機能１２２ｄにより収集したｔａｇ画像用ｋ空間データと、ｃｏｎｔｒｏｌ画像用ｋ空間データの２種類の４Ｄ（３Ｄ＋ｔ）ｋ空間データと、ステップＳ２０１で感度マップ生成機能１２２ａにより得られた感度マップと、ステップＳ２０２で推定機能１２２ｂで得られた重み係数を用いて、エイリアシングの無い、１つの再構成後のｔａｇ画像ｘ空間データと１つの再構成後のｃｏｎｔｒｏｌ画像ｘ空間データを生成する。アンフォールド機能１２２ｃ及び画像生成機能１２２ｅを有する処理回路１５０は、式（４）のエネルギー関数を最小化することで再構成を行う。

ここで、添え字のＴは、第１の画像（ｔａｇ画像）を示す添え字である。また、添え字のＯは、第２の画像（ｃｏｎｔｒｏｌ画像）を示す添え字である。また、ｘ_Ｔは、第１の画像、ｘ_Ｏは、第２の画像を表す。また、ｙ_Ｔｃは、ｃ番目のコイルに関する第１のｋ空間データ（ｔａｇ画像用のｋ空間データ）、ｙ_Ｏｃは、ｃ番目のコイルに関する第２のｋ空間データ（ｃｏｎｔｒｏｌ画像用のｋ空間データ）である。また、Ｅ（ｘ）は、エネルギー関数を表す。ここでｘ＝（ｘ_Ｔ、ｘ_Ｏ）、すなわち第１の画像及び第２の画像を表す。換言すると、エネルギー関数は、第１の画像ｘ_Ｔ及び第２の画像ｘ_Ｃの両方を引数に持つ関数である。

また、演算子Ａ_Ｔｃは、第１の画像に関するｃ番目のコイル感度をＳ_Ｔｃ、フーリエ変換の操作をＦＴ、第１の画像に係るパルスシーケンスにおけるサンプリング操作をΓ_ｃとすると、Ａ_Ｔｃ＝Γ_ｃ・ＦＴ・Ｓ_Ｔｃと表現することができる演算子である。演算子Ａ_Ｏｃは、第２の画像に関するｃ番目のコイル感度をＳ_Ｔｏ、フーリエ変換の操作をＦＴ、第２の画像に係るパルスシーケンスにおけるサンプリング操作をΓ_Ｏとすると、Ａ_Ｏｃ＝Γ_Ｏ・ＦＴ・Ｓ_Ｏｃと表現することができる演算子である。

また、Ｗは、処理回路１５０が推定機能１２２ｂにより生成した領域重み係数である。また、式（４）の第１項は、第１の画像の再構成条件ｙ_Ｔｃ＝Ａ_Ｔｃｘ_Ｔに関係する目的関数であり、式（４）の第２項は、第２の画像の再構成条件ｙ_Ｏｃ＝Ａ_Ｏｃｘ_Ｏに関係する目的関数であり、式（４）の第３項は、正則化項である。すなわち、処理回路１５０は、画像生成機能１２２ｅにより、正則化項における重みと、出力対象画像の各画素における値とを画素ごとに乗じて得られるベクトルに対するノルムを正則化項として目的関数に加え、目的関数と正則化項との和が最小になるような出力対象画像を算出することで、出力対象画像を生成する。

この際、第１段階の処理として、処理回路１５０は、画像生成機能１２２ｅにより正則化を行って、第１のｋ空間データ（ｔａｇ画像用のｋ空間データ）に対して再構成処理を行い第１の画像（ｔａｇ画像）を生成し、第２のｋ空間データ（ｃｏｎｔｒｏｌ画像用のｋ空間データ）に対して再構成処理を行い第２の画像（ｃｏｎｔｒｏｌ画像）を生成する。この時、第１の画像と第２の画像とは、式（４）により同時最適化されて生成される。続いて、第２段階の処理として、処理回路１５０は、画像生成機能１２２ｅにより、第１の画像と第２の画像とに基づいて出力対象画像（血管画像）を生成する。具体的には、処理回路１５０は、画像生成機能１２２ｅにより、第１の画像と第２の画像との間で差分処理を行い、差分画像である出力対象画像（血管画像）を生成する。

図１３及び図１４に、このようにして生成された出力対象画像（血管画像）の例が示されている。図１３及び図１４は、第１の実施形態に係る磁気共鳴イメージング装置が行う画像生成処理について説明した図である。

図１３は、実施形態に係る再構成処理により生成された画像の例を模式的に示したものである。図１３の左の図は、第１の画像（ｔａｇ画像）である。また、図１３の中央の図は、第２の画像（ｃｏｎｔｒｏｌ画像）である。処理回路１５０は、画像生成機能１２２ｅにより、第１の画像（図１３の左の図）から第２の画像（図１３の中央の図）を差し引くことにより、図１３の右の図で示される出力対象画像（血管画像（Ａｎｇｉｏｇｒａｐｈｙ ｉｍａｇｅ））が得られる。ここで、出力対象画像は、観察対象である血管のみが、高コントラストで描出された画像になっている。

図１４は、このようにして生成された出力対象画像の例を示している。パネル４１は３次元データの各断面を表している。例えば、画像４２は、Ｘ−Ｙ断面での画像である。また、画像４３は、Ｚ−Ｙ断面での画像である。また、画像４４は、Ｘ−Ｚ断面での画像である。

図１５は、第１の実施形態に係る磁気共鳴イメージング装置が行う画像生成処理を従来法と比較した図である。領域重み係数Ｗを用いない従来法と、領域重み係数Ｗを用いる実施形態に係る方法、いずれもＲｅｄｕｃｔｉｏｎ Ｆａｃｔｏｒ Ｒ＝４で比較を行っている。画像５０は、従来法のｔａｇ画像である。画像５1は、実施形態に係る方法のｔａｇ画像である。画像５２は、従来法のｃｏｎｔｒｏｌ画像である。画像５３は、実施形態に係る方法のｃｏｎｔｒｏｌ画像である。図１５の左側のｔａｇ画像から、図１５の中央のｃｏｎｔｒｏｌ画像を差し引くことで、図１５の右側の血管画像（Ａｎｇｉｏｇｒａｐｈｙ ｉｍａｇｅ（ｔａｇ−ｃｏｎｔｒｏｌ画像））を得ることができる。画像５４は、従来法の血管画像である。画像５５は、実施形態に係る方法の血管画像である。従来法の画像５４と実施形態に係る方法の画像５５を比較すると、従来法の画像５４は、折り返しの影響で画質が低下している（例えば領域６０）。実施形態に係る方法の画像５５は、背景のアーチファクトが領域重み係数Ｗを用いた正則化により抑制されているので、高画質な画像を得ることが可能になっている。

なお、式（４）の第３項の正則化項は、第１の画像と第２の画像との差分画像ｘ_Ｔ−ｘ_Ｏに対する領域重み係数ＷつきのＬ_２ノルムになっているので、第１のｋ空間データと第２のｋ空間データ（ｃｏｎｔｒｏｌ画像用ｋ空間データ）とが等しいｋ空間データである場合には０となるような正則化項になっている。このように、処理回路１５０は、画像生成機能１２２ｅにより、差分画像ｘ_Ｔ−ｘ_Ｏ（すなわち、標識化パルスが印加されたデータと標識化パルスが印加されていないデータの差分を取ることにより血管が描出されたデータ）に領域重み係数Ｗを乗じた正則化項を用いて正則化をすることで、血管画像を効果的に描出することができる。

また、本実施形態では、正則化は時系列データの各時刻ごとに行われる。換言すると、処理回路１５０は、画像生成機能１２２ｅにより、時系列ｋ空間データのうち各時刻ごとに画像再構成処理を行い出力対象画像を生成する。

以上をまとめると、第１の実施形態では、処理回路１５０は、画像生成機能１２２ｅにより、領域重み係数Ｗを用いて正則化を行って、ｋ空間データに対して画像再構成処理を行い出力対象画像を生成する。この際、処理回路１５０は、推定機能１２２ｂにより観察対象領域に含まれていると推定した画素位置に対する重みが、観察対象領域に含まれていないと推定した画素位置に対する重みよりも、小さな重みとなるような領域重み係数Ｗを用いる。

第１の実施形態では、画像生成機能１２２ｅを有する処理回路１５０が正則化項としてＬ_２ノルムを用いて正則化を行う場合について説明したが、実施形態はこれに限られない。例えば、画像生成機能１２２ｅを有する処理回路１５０は、正則化項として、その他のノルム、例えばＬ_１ノルムを用いて正則化を行っても良い。また、画像生成機能１２２ｅを有する処理回路１５０は、その他の正則化をさらに追加して正則化を行っても良い。例えば、画像生成機能１２２ｅを有する処理回路１５０は、ＴＶ（Ｔｏｔａｌ Ｖａｒｉａｔｉｏｎ）正則化項を追加して正則化を行っても良い。

第１の実施形態では、シーケンス制御回路１１０が、撮影領域内の標識化領域に対して標識化パルスを印加する第１のパルスシーケンスと、標識化パルスを印加する第２のパルスシーケンスとを実行する場合について説明したが、実施形態はこれに限られない。例えば、シーケンス制御回路１１０は、被検体に造影剤を用いた後に実行する第１のパルスシーケンスと、被検体に造影剤を用いない状態で実行する第２のパルスシーケンスとを実行する場合であってもよい。

第１の実施形態では、収集機能１２２ｄを有する処理回路１５０が、時系列ｋ空間データを収集する場合について説明したが、実施形態はこれに限られない。例えば、収集機能１２２ｄを有する処理回路１５０は、時系列ｋ空間データではなく、単時相でｋ空間データを収集してもよい。

第１の実施形態では、シーケンス制御回路１１０が、ｔａｇ画像生成用の第１のパルスシーケンスと、ｃｏｎｔｒｏｌ画像生成用の第２のパルスシーケンスとを実行し、収集機能１２２ｄを有する処理回路１５０が、ｔａｇ画像生成用とｃｏｎｔｒｏｌ画像生成用の２つのｋ空間データを収集し、処理回路１５０が、画像生成機能１２２ｅにより、ｔａｇ画像及びｃｏｎｔｒｏｌ画像を再構成する場合について説明したが、実施形態はこれに限られない。例えば、シーケンス制御回路１１０が、一つのパルスシーケンスを実行し、収集機能１２２ｄを有する処理回路１５０が、一つのｋ空間データを収集し、処理回路１５０が、画像生成機能１２２ｅにより、当該ｋ空間データから出力対象の画像を生成してもよい。この場合、処理回路１５０は、画像生成機能１２２ｅにより、例えば式（３）に基づいて正則化を行って最適化を行い出力対象の画像を生成する。

第１の実施形態によれば、再構成前のｔａｇ／ｃｏｎｔｒｏｌのｋ空間データから血管領域画素を推定し、正則化項に導入することで、解の安定性を向上し、かつ、背景領域のアーチファクトを抑制することができる。

（第１の実施形態の第１の変形例）
第１の実施形態では、推定機能１２２ｂを有する処理回路１５０が算出する観察対象領域の領域重み係数Ｗが、観察対象領域ではない領域の領域重み係数Ｗより小さい場合について説明したが、実施形態はこれに限られない。例えば、推定部１２２ｂを有する処理回路１５０は、観察対象領域の領域重み係数Ｗが、観察対象領域ではない領域の領域重み係数より大きくなるように算出してもよい。この場合、例えば、推定機能１２２ｂを有する処理回路１５０は、観察対象領域の重み係数Ｗが１、観察対象ではない領域の重み係数Ｗが０となるように、領域重み係数Ｗを算出する。これと連動して、画像生成機能１２２ｅを有する処理回路１５０は、例えば、以下の式（５）のような正則化項を用いて画像生成
を行う。

また、実施形態は、画像生成機能１２２ｅを有する処理回路１５０が、領域重み係数Ｗと、出力対象画像ｘの積の形の正則化項を用いて画像生成を行う場合について説明したが、実施形態はこれに限らない。例えば、画像生成機能１２２ｅを有する処理回路１５０は、その他の形の正則化項を用いて画像生成を行っても良い。例えば、画像生成機能１２２ｅを有する処理回路１５０は、式（４）や式（５）の正則化項に対して、ノルムをとったものに対して連続関数を作用させたものを正則化項として使用してもよいし、逆に連続関数を作用させてからノルムをとったものを正則化項として使用してもよい。また、画像生成機能１２２ｅを有する処理回路１５０は、直交関数等所定の基底で展開してからノルムをとったものを、正則化項として使用してもよい。

（第２の実施形態）
第１の実施形態では、画像生成機能１２２ｅを有する処理回路１５０が、式（４）を用いて、ｔａｇ画像とｃｏｎｔｒｏｌ画像を同時に最適化、再構成を行い、再構成されたｔａｇ画像と、再構成されたｃｏｎｔｒｏｌ画像との間に差分処理を行って、出力対象の画像を生成する場合について説明した。実施形態はこれに限られない。例えば、ｔａｇ画像用のｋ空間データとｃｏｎｔｒｏｌ画像用のｋ空間データの間で最初に差分処理を行い、差分処理を行ったデータに対して画像再構成が行われても良い。

図１６は、第２の実施形態に係る磁気共鳴イメージング装置１００の有する計算機システム１２０を表すブロック図である。処理回路１５０は、第１の実施形態に加えて、複素差分機能１２２ｆを備える。複素差分機能１２２ｆ以外の機能や、処理回路１５０以外の構成要素は、第１の実施形態に係る磁気共鳴イメージング装置と同様である。

第２の実施形態では、画像生成機能１２２ｅを有する処理回路１５０は、ｔａｇ画像とｃｏｎｔｒｏｌ画像との差分画像に対して画像再構成処理を行う。具体的には、処理回路１５０は、複素差分機能１１２ｆにより、図５のステップＳ２０３の処理の開始までに、ｔａｇとｃｏｎｔｒｏｌの２種類の４Ｄ ｋ空間データの複素差分を算出し、算出された複素差分データを、アンフォールド機能２２２ｃのために用意する。これは処理回路１５０が領域推定機能１２２ｂによりステップＳ２０２ａにおいて行う処理と同様の処理であるため、処理回路１５０は、複素差分機能２２２ｄによる処理を、ステップＳ２０２ａに変えても良い。

アンフォールド機能１２２ｃ及び画像生成機能１２２ｅを有する処理回路１５０は、複素差分機能２２２ｄにより生成された差分４Ｄｋ空間データと、感度マップ生成機能１２２ａにより得られた感度マップと、推定機能１２２ｂにより得られた重み係数を用いて、エイリアシングの無い、１つの差分ｘ空間データを生成する。具体的には、画像生成機能１２２ｅを有する処理回路１５０は、式（６）のエネルギー関数を最小化することで画像の再構成を行う。

ここで、添え字のＡは、ｔａｇ画像とｃｏｎｔｒｏｌ画像との差分画像である血管画像を表す。第２の実施形態に係る磁気共鳴イメージング装置によれば、差分画像のみを再構成するため、再構成処理にかかるメモリ量及び計算時間を削減することができる。

（画像処理装置）
上述した実施形態においては、医用画像診断装置である磁気共鳴イメージング装置１００が各種処理を実行する場合を説明したが、実施形態はこれに限られるものではない。例えば、磁気共鳴イメージング装置１００に替わり、画像処理装置や、磁気共鳴イメージング装置１００と画像処理装置とを含む画像処理システムが、上述した各種処理を実行してもよい。ここで、画像処理装置とは、例えば、ワークステーション、ＰＡＣＳ（Picture Archiving and Communication System）の画像保管装置（画像サーバ）やビューワ、電子カルテシステムの各種装置等である。この場合、例えば、画像処理装置は、磁気共鳴イメージング１００によって収集されたｋ空間データを、磁気共鳴イメージング装置１００から、若しくは、画像サーバからネットワーク経由で受信することで、あるいは、記録媒体を介して操作者から入力されること等で、受け付けて、記憶回路１２３に記憶する。そして、画像処理装置は、記憶回路１２３に記憶したこのｋ空間データを対象として、上述した各種処理を実行すれば良い。

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

また、記憶媒体からコンピュータや組み込みシステムにインストールされたプログラムの指示に基づきコンピュータ上で稼働しているＯＳ（オペレーティングシステム）や、データベース管理ソフト、ネットワーク等のＭＷ（ミドルウェア）等が、上述した実施形態を実現するための各処理の一部を実行してもよい。

更に、記憶媒体は、コンピュータあるいは組み込みシステムと独立した媒体に限らず、ＬＡＮ（Local Area Network）やインターネット等により伝達されたプログラムをダウンロードして記憶又は一時記憶した記憶媒体も含まれる。

また、記憶媒体は１つに限られず、複数の媒体から、上述した実施形態における処理が実行される場合も、実施形態における記憶媒体に含まれ、媒体の構成は何れの構成であってもよい。

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

また、実施形態におけるコンピュータとは、パソコンに限らず、情報処理機器に含まれる演算処理装置、マイコン等も含み、プログラムによって実施形態における機能を実現することが可能な機器、装置を総称している。

以上述べた少なくとも一つの実施形態によれば、パラレルイメージングにおける出力対象画像の画質を向上することができる。

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

１２２ｂ 推定機能
１２２ｄ 収集機能
１２２ｅ 画像生成機能

Claims (13)

1. アンダーサンプリングを行いながらパルスシーケンスを実行するシーケンス制御部と、
   前記パルスシーケンスを用いてｋ空間データを収集する収集部と、
   前記ｋ空間データを用いて生成されたデータに基づいて、出力対象画像の各画素の画素位置が観察対象領域に含まれているか否かを推定する推定部と、
   前記推定部が観察対象領域に含まれていると推定した画素位置に対する重みと、前記推定部が観察対象領域に含まれないと推定した画素位置に対する重みとが異なる値となるような正則化項を用いて、単時相において、または時系列データの各時刻ごとに正則化を行って、前記ｋ空間データに対して画像再構成処理を行い前記出力対象画像を生成する生成部と
   を備える磁気共鳴イメージング装置。
2. 前記生成部は、前記推定部が観察対象領域に含まれていると推定した画素位置に対する重みが、前記推定部が観察対象領域に含まれないと推定した画素位置に対する重みよりも、小さな重みとなる前記正則化項を用いて前記出力対象画像を生成する、請求項１に記載の磁気共鳴イメージング装置。
3. 前記正則化項に用いられる重みは、各画素の画素位置が前記観察対象領域に含まれている場合には０であり、各画素の画素位置が前記観察対象領域に含まれていない場合には１である、請求項１又は２に記載の磁気共鳴イメージング装置。
4. アンダーサンプリングを行いながらパルスシーケンスを実行するシーケンス制御部と、
   前記パルスシーケンスを用いてｋ空間データを収集する収集部と、
   前記ｋ空間データを用いて生成されたデータに基づいて、出力対象画像の各画素の画素位置が観察対象領域に含まれているか否かを推定する推定部と、
   前記推定部が観察対象領域に含まれていると推定した画素位置に対する重みと、前記推定部が観察対象領域に含まれないと推定した画素位置に対する重みとが異なる値となるような正則化項を用いて正則化を行って、前記ｋ空間データに対して画像再構成処理を行い前記出力対象画像を生成する生成部と
   を備え、
   前記生成部は、前記正則化項における重みと、前記出力対象画像の各画素における値とを画素ごとに乗じて得られるベクトルに対するノルムを前記正則化項として目的関数に加え、前記目的関数と前記正則化項との和が最小になるような前記出力対象画像を算出することで、前記出力対象画像を生成する、磁気共鳴イメージング装置。
5. アンダーサンプリングを行いながらパルスシーケンスを実行するシーケンス制御部と、
   前記パルスシーケンスを用いてｋ空間データを収集する収集部と、
   前記ｋ空間データを用いて生成されたデータに基づいて、出力対象画像の各画素の画素位置が観察対象領域に含まれているか否かを推定する推定部と、
   前記推定部が観察対象領域に含まれていると推定した画素位置に対する重みと、前記推定部が観察対象領域に含まれないと推定した画素位置に対する重みとが異なる値となるような正則化項を用いて正則化を行って、前記ｋ空間データに対して画像再構成処理を行い前記出力対象画像を生成する生成部と
   を備え、
   前記シーケンス制御部は、第１のパルスシーケンスと、第２のパルスシーケンスとを、
   前記第１のパルスシーケンスと前記第２のパルスシーケンスとのうち少なくとも一方をアンダーサンプリングを行いながら実行し、
   前記収集部は、前記第１のパルスシーケンスを用いて第１のｋ空間データを収集し、前記第２のパルスシーケンスを用いて第２のｋ空間データを収集し、
   前記推定部は、前記第１のｋ空間データと前記第２のｋ空間データとを用いて生成されたデータに基づいて、前記出力対象画像の各画素の画素位置が観察対象領域に含まれているか否かを推定し、
   前記生成部は、前記正則化を行って、前記第１のｋ空間データに対して画像再構成処理を行い第１の画像を生成し、前記第２のｋ空間データに対して画像再構成処理を行い第２の画像を生成し、前記第１の画像と前記第２の画像とに基づいて前記出力対象画像を生成する、磁気共鳴イメージング装置。
6. 前記正則化項は、前記第１のｋ空間データと前記第２のｋ空間データとが等しいｋ空間データである場合には０になるような正則化項である、請求項５に記載の磁気共鳴イメージング装置。
7. 前記第１のパルスシーケンスは、撮像領域内の標識化領域に対して標識化パルスを印加するパルスシーケンスであり、前記第２のパルスシーケンスは、標識化パルスを印加しないパルスシーケンスである、請求項５に記載の磁気共鳴イメージング装置。
8. 前記第１のパルスシーケンスは、被検体に造影剤を用いた後に実行されたパルスシーケンスであり、前記第２のパルスシーケンスは、前記被検体に造影剤を用いない状態で実行されたパルスシーケンスである、請求項５に記載の磁気共鳴イメージング装置。
9. アンダーサンプリングを行いながらパルスシーケンスを実行するシーケンス制御部と、
   前記パルスシーケンスを用いてｋ空間データを収集する収集部と、
   前記ｋ空間データを用いて生成されたデータに基づいて、出力対象画像の各画素の画素位置が観察対象領域に含まれているか否かを推定する推定部と、
   前記推定部が観察対象領域に含まれていると推定した画素位置に対する重みと、前記推定部が観察対象領域に含まれないと推定した画素位置に対する重みとが異なる値となるような正則化項を用いて正則化を行って、前記ｋ空間データに対して画像再構成処理を行い前記出力対象画像を生成する生成部と
   を備え、
   前記正則化項に用いられる重みは、各画素の画素位置が前記観察対象領域に含まれている確率の推定値に基づいて算出された値である、磁気共鳴イメージング装置。
10. アンダーサンプリングを行いながらパルスシーケンスを実行するシーケンス制御部と、
    前記パルスシーケンスを用いてｋ空間データを収集する収集部と、
    前記ｋ空間データを用いて生成されたデータに基づいて、出力対象画像の各画素の画素位置が観察対象領域に含まれているか否かを推定する推定部と、
    前記推定部が観察対象領域に含まれていると推定した画素位置に対する重みと、前記推定部が観察対象領域に含まれないと推定した画素位置に対する重みとが異なる値となるような正則化項を用いて正則化を行って、前記ｋ空間データに対して画像再構成処理を行い前記出力対象画像を生成する生成部と
    を備え、
    前記正則化項に用いられる重みは、前記ｋ空間データの収集単位であるチャネルごとに算出された値を統合することにより算出された値である、磁気共鳴イメージング装置。
11. アンダーサンプリングを行いながらパルスシーケンスを実行するシーケンス制御部と、
    前記パルスシーケンスを用いてｋ空間データを収集する収集部と、
    前記ｋ空間データを用いて生成されたデータに基づいて、出力対象画像の各画素の画素位置が観察対象領域に含まれているか否かを推定する推定部と、
    前記推定部が観察対象領域に含まれていると推定した画素位置に対する重みと、前記推定部が観察対象領域に含まれないと推定した画素位置に対する重みとが異なる値となるような正則化項を用いて正則化を行って、前記ｋ空間データに対して画像再構成処理を行い前記出力対象画像を生成する生成部と
    を備え、
    前記収集部は、前記ｋ空間データとして、時系列ｋ空間データを収集し、
    前記推定部は、前記時系列ｋ空間データを用いて生成された時間軸を持たない空間データに基づいて、前記出力対象画像の各画素の前記画素位置が前記観察対象領域に含まれているか否かを推定する、磁気共鳴イメージング装置。
12. アンダーサンプリングを行いながらパルスシーケンスを実行するシーケンス制御部と、
    前記パルスシーケンスを用いてｋ空間データを収集する収集部と、
    前記ｋ空間データを用いて生成されたデータに基づいて、出力対象画像の各画素の画素位置が観察対象領域に含まれているか否かを推定する推定部と、
    前記推定部が観察対象領域に含まれていると推定した画素位置に対する重みと、前記推定部が観察対象領域に含まれないと推定した画素位置に対する重みとが異なる値となるような正則化項を用いて正則化を行って、前記ｋ空間データに対して画像再構成処理を行い前記出力対象画像を生成する生成部と
    を備え、
    前記収集部は、前記ｋ空間データとして、時系列ｋ空間データを収集し、
    前記生成部は、前記時系列ｋ空間データのうち各時刻ごとに、前記正則化を行って、前記時系列ｋ空間データに対して前記各時刻ごとに前記画像再構成処理を行い前記出力対象画像を生成する、磁気共鳴イメージング装置。
13. アンダーサンプリングを行いながら実行されたパルスシーケンスを用いて収集されたｋ空間データを用いて生成されたデータに基づいて、出力対象画像の各画素の画素位置が観察対象領域に含まれているか否かを推定する推定部と、
    前記推定部が観察対象に含まれていると推定した画素位置に対する重みと、前記推定部が観察対象領域に含まれないと推定した画素位置に対する重みとが異なる値となるような正則化項を用いて、単時相において、または時系列データの各時刻ごとに正則化を行って、前記ｋ空間データに対して画像再構成処理を行い前記出力対象画像を生成する生成部と
    を備える画像処理装置。