JP3928106B2

JP3928106B2 - Magnetic resonance imaging device

Info

Publication number: JP3928106B2
Application number: JP18927599A
Authority: JP
Inventors: 智嗣平田; 良孝尾藤; 啓二塚田
Original assignee: Hitachi Medical Corp
Current assignee: Hitachi Healthcare Manufacturing Ltd
Priority date: 1999-07-02
Filing date: 1999-07-02
Publication date: 2007-06-13
Anticipated expiration: 2019-07-02
Also published as: JP2001008919A

Description

【０００１】
【発明の属する技術分野】
本発明は、磁気共鳴撮影装置に係り、特にケミカルシフトに関する情報の空間分布を測定するケミカルシフト情報処理方法に関する。
【０００２】
【従来の技術】
磁気共鳴撮影装置は、静磁場中に置かれた被検体に対し、特定周波数の高周波磁場を照射して磁気共鳴現象を引き起こさせ、これにより被検体から発生する磁気共鳴信号を利用して、被検体を構成する物質の化学的・物理的な情報を計測し、これに基づいて例えば画像を作成する。特に、磁気共鳴スペクトロスコピックイメージング(Magnetic Resonance Spectroscopic Imaging 以下、MRSIと略す)は、様々な分子の化学結合の違いによる磁気共鳴周波数の差異を、ケミカルシフト情報として測定し、各分子ごとの空間分布像をケミカルシフト画像として得る方法として知られている。
【０００３】
磁気共鳴信号からケミカルシフト画像を取得するための情報処理方法は、一般に、周波数及び位相にエンコードした空間情報を付加して計測した一連の磁気共鳴信号（エコートレイン信号）を、ケミカルシフト軸を含めて計測軸毎に並べかえた後、各軸毎に逆フーリエ変換を施す。つまり、例えば計測する空間がケミカルシフト軸（δ＝時間軸)、X軸、 Y軸及びZ軸の４次元実空間である場合、フーリエ変換で結ばれる４次元フーリエ空間の各軸(kδ軸、kx軸、ky軸及びkz軸)が計測軸に相当する。そして、ケミカルシフト画像の画質及びスペクトルの精度を向上させるために、逆フーリエ変換する前に、窓関数フィルタリング、ゼロフィリング等の前処理を施したり、逆フーリエ変換後に位相補正等の後処理を行うことが適宜なされている。
【０００４】
【発明が解決しようとする課題】
ところで、プロトン（¹H）を含む化合物を測定対象とする¹H−MRSIでは、微量な代謝物質からの信号を検出するために、代謝物質の信号レベルに対して極めて大きな信号レベルを有する水信号や皮下脂肪信号の発生を抑える必要がある。例えば、水信号は代謝物質信号の数千倍程度、皮下脂肪信号は代謝物質信号の数百倍程度に達する。特に、水信号の抑圧は重要であり、通常、励起帯域を制限した3個のガウス型高周波磁場パルス (帯域を1〜2 ppm程度に制限し、水信号のみを選択励起する作用をもつ)と、各軸のスポイラー傾斜磁場(励起した磁化の位相をバラバラにさせ信号の発生を抑える作用をもつ)からなる水信号抑圧シーケンスが用いられる。この処理により水信号強度を非抑圧時信号強度の百分の一程度まで抑圧することができるが、抑圧しきれなかった水信号が残留水信号として残り、トランケーションアーチファクトを生じる場合がある。
【０００５】
ここで、残留水信号によるトランケーションアーチファクトの発生が問題となる理由について説明する。核磁化の励起により生じた全ての磁気共鳴信号は時間とともに減衰する。この過程を緩和と呼び、緩和には縦磁化の緩和(T1緩和)と横磁化の緩和がある。横磁化緩和は、さらに静磁場不均一に起因するT2＊緩和と、核磁化の内的状態に起因するT2緩和とに分けられる。静磁場不均一が大きく180°パルスの照射が有る場合は、T2＊緩和が支配的となって、エコータイム（TE）に最大値を有する釣り鐘型の包絡線のエコートレイン信号が発生する。一方、静磁場不均一の無い理想的な状態の場合は1サンプリング期間中の信号減衰はT2緩和にのみ支配され、180°パルスの有無にかかわらず常に指数関数的(exp(−t／T2)に減衰する。ここで、tはサンプリング時刻、T2は横緩和時間である。すなわち、包絡線が指数関数型のエコートレイン信号となる。したがって、通常は、信号は釣り鐘型と指数関数型を組み合わせた形状となる。ただ、一般にMRSI実行時には、スペクトル分解能を高め代謝物質画像の画質を良くするために、オートシミング等で静磁場均一度を向上させることが多く、この場合はT2緩和が支配的となって、指数関数型の信号形状に近づく。
【０００６】
指数関数型の信号に対して逆FFTを行うと、信号の始点と終点の信号レベルが異なるので、不連続性(ギャップ)に起因するギブス現象(不連続点近傍での周期的増減)、すなわちトランケーションアーチファクトが発生してしまう。また、このトランケーションアーチファクトの信号強度はギャップの大きさに比例するため、高信号強度でT2の短い残留水信号(T2≒100 ms)によるトランケーションアーチファクトが、低信号強度でT2の長い代謝物質の信号(T2≒250〜400 ms)に与える影響が非常に大きな問題となる。
【０００７】
このトランケーションアーチファクトを低減させるためには、前述したように窓関数フィルタリング処理を施し、信号の始点と終点を滑らかに接続する方法が有効である。しかし、このフィルタリング処理には、残留水信号だけではなく、代謝物質信号のピーク幅も増大させてしまうブロードニング作用があるため、スペクトル分解能を低下させるという問題がある。つまり、窓関数フィルタリング処理には、指数関数型の信号形状を釣り鐘型に変形させる作用があり、この作用は静磁場均一度を低下させた場合と同じ影響をスペクトルに与えるためである。
【０００８】
本発明が解決しようとする課題は、ケミカルシフト情報のスペクトル分解能を向上させることにある。
【０００９】
【課題を解決するための手段】
上記課題は、磁気共鳴撮影により計測されるエンコード情報が付加されたエコートレイン信号を取り込み、該エコートレイン信号をケミカルシフト軸に沿って配列した後、前記エンコード情報に係る軸方向の逆フーリエ変換を施し、該逆フーリエ変換して得られる時系列データに基づいて、該時系列データを時間軸方向の所定の範囲にわたって外挿により補完して拡張時系列データを生成し、該拡張時系列データを当該拡張時系列データの近似関数で除して正規化し、該正規化された拡張時系列データに基づいて、前記ケミカルシフト軸方向の逆フーリエ変換を施してケミカルシフト情報を得るケミカルシフト情報処理方法とすることにより、解決することができる。
【００１０】
すなわち、ケミカルシフト軸方向の逆フーリエ変換によりケミカルシフト情報を得る前に、空間情報に係る軸方向の逆フーリエ変換を施して得られる時系列データを外挿法により、計測していない時間帯のデータを補完する。これにより得られた拡張時系列データを正規化処理することにより、拡張時系列データの始点と終点の信号レベルを同等にして、つまりＴ２緩和による減衰を補償して、不連続性（ギャップ）を大幅に緩和ないし低減できる。その結果、ケミカルシフト軸方向の逆フーリエ変換により発生するトランケーションアーチファクトを低減して、ケミカルシフト情報のスペクトル分解能を向上できる。
【００１１】
また、上記処理に加えて、正規化された拡張時系列データを時系列データの時間範囲を含む窓関数フィルタによりフィルタ処理し、このフィルタ処理された拡張時系列データに基づいて、ケミカルシフト軸方向の逆フーリエ変換を施してケミカルシフト情報を得ることが好ましい。これによれば、一層、ケミカルシフト情報のスペクトル分解能を向上できる。
【００１２】
本発明に係る磁気共鳴撮影装置は、静磁場発生手段と、互いに異なる3方向の傾斜磁場を発生する傾斜磁場発生手段と、高周波磁場発生手段と、磁気共鳴信号を検出する信号検出手段と、該信号検出手段により得られた磁気共鳴信号に基づいて処理演算を行う演算手段と、記憶手段と、前記各手段の動作を制御するシーケンス制御手段とを含んで構成され、前記演算手段に、上記のケミカルシフト情報処理機能を組み込むことにより実現できる。
【００１３】
また、本発明に係るケミカルシフト情報処理は、磁気共鳴撮影装置により計測された一連の磁気共鳴信号（エコートレイン信号）を、信号伝送路又は記録媒体を介して別途設けられたコンピュータに取り込み、そのコンピュータで処理することができる。
【００１４】
なお、計測していない時系列データの外挿を行う際、線型予測法を用いることができる。また、時系列データを近似する関数を指数関数とし、この近似指数関数を求める際に、最小二乗法を用いることができる。また、制御手段は、互いに異なる3方向の傾斜磁場のうち少なくとも1方向の傾斜磁場の勾配極性を周期的に反転させて印加しながら磁気共鳴信号を検出させる制御をするのが好ましい。
【００１５】
【実施の形態】
以下、本発明の実施の形態を図を用いて説明する。図１に本発明のケミカルシフト情報処理方法の手順を示し、図２に本発明のケミカルシフト情報処理法を適用してなる一実施形態の磁気共鳴撮影装置の全体構成を示す。図３に本発明を適用可能な計測シーケンスの一例を示す。
【００１６】
図２の磁気共鳴撮影装置において、被検体1は、静磁場発生マグネット2により生成される静磁場及び傾斜磁場発生コイル3により生成される異なる3方向の傾斜磁場が印加される空間に置かれる。シムコイル１１は複数のコイルを備え、各コイルに流す電流を変化させることにより、該静磁場の均一度を調整することができる。しかし、シムコイル１１は必ずしも備える必要はない。プロ−ブ4は被検体１に照射する高周波磁場を発生して被検体１に磁気共鳴現象を起こさせるとともに、被検体1から発生する磁気共鳴信号をプロ−ブ4により検出する。計算機5は、磁気共鳴信号を処理して必要な情報、例えばケミカルシフトのスペクトル、ケミカルシフト画像を生成し、ディスプレイ6に表示させる。シムコイル１１の駆動用電源部１２、傾斜磁場発生コイル３の駆動用電源部７、送信器８及び受信器９は、シーケンス制御装置１０により制御される。また必要に応じて、記憶媒体１３に計測信号及び測定条件等を記憶させるようになっている。
【００１７】
図３のパルスシーケンスは、高速に磁気共鳴スペクトル画像が得られる高速MRSIの一手法であるエコープラナースペクトロスコピックイメージングであり、これは特開昭61-13143等に記載された松井等による提案の一例である。同図において、上から順に、励起用高周波磁場パルスRF、Z軸方向の傾斜磁場Gz、X軸方向の傾斜磁場Gｘ、Y軸方向の傾斜磁場Gy、磁気共鳴信号Sigのタイミングを示している。このパルスシーケンスでは、初めに励起用高周波磁場パルスRF1 (90°パルス) とZ軸方向のスライス (ボリューム) 選択用傾斜磁場Gs1を同時に印加する。これにより、測定対象とする領域内に含まれる核スピンを選択的に励起する。次に、励起により生じる磁化に対して、Y軸及びZ軸方向の位相エンコード用傾斜磁場Ge1及びGe2を印加することにより、各磁化の位相にY軸及びZ軸方向の空間情報を付与する。さらに、励起パルスRF1の印加からTE / 2後に (TE : エコータイム) 、反転用高周波磁場パルスRF2 (180°パルス) とZ軸方向のスライス(ボリューム)選択用傾斜磁場Gs2を同時に印加する。これにより、静磁場不均一により位相が互いにずれていた核スピンを再び収束させ、RF2の印加からTE / 2後にスピンエコー信号を発生させる。そして、このスピンエコー信号を計測する際に、X軸方向のリードアウト用傾斜磁場Grの勾配を周期的に反転させて印加することにより、ケミカルシフトとX軸方向の空間情報を含んだエコートレイン信号を発生させることができる。Grの波形としては、台形波もしくは正弦波（ｓｉｎ波）を用いる。図示例は、Grに台形波を用いた例であり、極性切り替え部が滑らかになっていないが、連続的に移行するように切り替えてもよい。
【００１８】
この一連のパルスシーケンスを、Ge1及びGe2の印加強度を図中に破線で示すように段階的に変化させ、Y軸及びZ軸方向のピクセル数に相当する回数分、TRの間隔で繰り返す。例えば、空間マトリクス数mx×my×mz の画像を得るためには、現在、最も広範に用いられている従来法(4D-CSI)では励起・計測をmx×my×mz 回繰り返す必要がある。しかし、図３の高速MRSIではmy×mz回の励起・計測で済むため、計測時間を1 / mxに短縮することが可能となる。
【００１９】
次に、この高速MRSIにおけるデータ処理方法について述べる。上記計測で得られるエコートレイン信号は、リードアウト用傾斜磁場Grによる外的磁気相互作用と、ケミカルシフトに相当する内的磁気相互作用の双方により周波数変調されている。すなわち、各エコー信号1つ1つに注目した場合、現在MRIにおいて広範に用いられている超高速撮影法EPI(Echo Planar Imaging)で得られる各エコー信号と同様に、リードアウト軸方向の空間情報が周波数エンコードされていることがわかる。また、各エコー信号の中心点に注目した場合、これらの点ではGrの傾斜磁場の勾配極性を周期的に反転させて印加していることから、Grの影響がキャンセルされている。そのため、各エコー信号の中心点を結んだエコートレイン信号の包絡線（envelope）は、ケミカルシフトにのみ変調されていることになる。したがって、以下に述べるデータ処理を行うことにより、エコートレイン信号からケミカルシフト情報とX軸方向の空間情報を抽出することができる。
【００２０】
すなわち、図４に示すように、位相変調された1つのエコートレイン信号E1、E2、E３、E4、E５、E6、…をエコー毎に分割し、奇数番目のエコーと偶数番目のエコーにグループ分けした後、時間軸kδ軸に沿って2次元的な再配列を行う。奇数番目と偶数番目ごとにグループ分けする理由は、Grの勾配極性が同じエコーごとに分けて処理することにより、アーチファクトを低減するためである。再配列した２次元データ配列において、新しく定義した時間軸kδ軸に沿った信号の時間変化は、純粋にケミカルシフトに相当する内的磁気相互作用によって決まる。これはGrの反転が周期的であることにより、Grによる外的磁気相互作用がkδ軸に沿って消失してしまうためである。一方、kx軸に沿った個々のエコー信号の形状はGrによる外的磁気相互作用により支配される。Y軸及びZ軸方向の空間情報は、Ge₁及びGe₂により、kx軸に垂直なky及びkz軸に沿って位相エンコードされる。このようにして一連のエコー信号を取得することにより、4次元実空間(ケミカルシフトδ軸、X軸、 Y軸及びZ軸)とフーリエ変換で結ばれる4次元フーリエ空間(kδ軸、kx軸、ky軸及びkz軸)を計測することができる。
【００２１】
ここで、各軸のマトリクス数は、リードアウト軸（X）方向では、EPIと同様に1エコー当たりのサンプリング点数で決まり、ケミカルシフト軸方向でのスペクトル点数は、1つのエコートレインに含まれるエコー数の1 / 2で定まり、位相エンコード軸方向のピクセル数は、Ge₁及びGe₂のステップ数で決まる。そして、これらの4次元データ配列に4次元逆フーリエ変換を施すことにより、ケミカルシフト画像を得ることが可能となる。奇数番目及び偶数番目のエコーグループから作成される2枚のケミカルシフト画像は、位相補正後、足し合わされる。なお、Grにサイン波駆動を用いた場合には、フーリエ空間を走査する軌跡がサイン波状となるので、サイン補正を含んだ4次元逆フーリエ変換を施す必要がある。
【００２２】
次に、本発明の特徴に係るケミカルシフト情報の処理法について詳しく説明する。計算機５は、磁気共鳴撮影により計測される空間情報がエンコ−ドされたエコートレイン信号Sigを取り込み、エコートレイン信号をケミカルシフト軸δに沿って配列した後、空間情報に係る軸方向（Ｘ,Ｙ,Ｚ）とケミカルシフト軸方向に逆フーリエ変換を施して、ケミカルシフト情報を得る。この際に、逆フーリエ変換して得られる時系列データ（元データ）に、窓関数フィルタリング処理を直接施した場合、前述したように、代謝物質信号に与えるブロードニング効果が大きすぎる。
【００２３】
そこで、本発明は、図１に示すように、ケミカルシフト軸方向の逆フーリエ変換を行う前に、周波数と位相にエンコードした空間情報に係る軸方向（Ｘ,Ｙ,Ｚ）に逆フーリエ変換して得られる時系列データを、時間軸方向の所定の範囲にわたって外挿により補完して拡張時系列データを生成し、これにより得られた拡張時系列データをその近似関数で除して正規化し、正規化された拡張時系列データにケミカルシフト軸方向の逆フーリエ変換を施してケミカルシフト情報を得ることを特徴とする。
【００２４】
（ステップＳ１）
図３のパルスシーケンスを実施して４次元ＭＲＳＩの計測を行い、周波数エンコードおよび位相エンコードした空間情報が付加されたエコートレイン信号Sig（磁気共鳴信号）を収集する。
【００２５】
（ステップＳ２）
エコートレイン信号をケミカルシフト軸δに沿って配列し、４次元データ配列に並べかえる。
【００２６】
（ステップＳ３）
それぞれの空間情報に係る軸方向（Ｘ,Ｙ,Ｚ）の逆フーリエ変換を施して時系列データを得る。
【００２７】
（ステップＳ４）
ステップＳ３により得られる各ボクセル（voxel）の時系列データに対して、線形予測を用いた外挿処理を行う。例えば、図５（a）に示す１次元データ（１voxel）の時系列データ（元データ）に対し、元データ取得時刻よりも早い時刻のデータと元データ取得時刻よりも遅い時刻のデータを、元データから外挿して予測することにより、図５(b)に示す拡張した時系列データを生成する。
【００２８】
（ステップＳ５）ステップＳ４により得られた各ボクセルの拡張時系列データに対して、指数関数へのフィッティング処理を行う。すなわち、拡張時系列データの減衰を近似する近似関数を求める。これにより得られる近似関数（フィティング関数）は、例えば、図５(d)に示す指数関数となる。
【００２９】
（ステップＳ６）
拡張時系列データと、近似関数である指数関数の逆数を乗算して正規化処理をする。つまり、拡張時系列データを近似関数で除算して、Ｔ２緩和の減衰を補償(正規化)する処理を行う。これにより、例えば、図５(c)に示す正規化された拡張時系列データが得られる。
【００３０】
（ステップＳ７）
正規化された各ボクセルの拡張時系列データに対して、窓関数フィルタリング処理を施す。窓関数フィルタリング処理に用いるフィルタリング関数は、図６に一例を示したように、サンプリング時間軸の始点及び終点の近傍の信号強度を選択的に減衰するパターンのものである。
【００３１】
（ステップＳ８）
フィルタリング処理された各ボクセルの拡張時系列データに対して、δ軸方向の逆フーリエ変換（ＦＦＴ）処理を施す。これにより、ケミカルシフト情報のスペクトルデータが得られる。
【００３２】
図１のケミカルシフト情報処理を施して得られるケミカルシフト情報のスペクトルデータを、従来法により得られるものと比較しながら、図７を参照して説明する。図７は、図５と同様に、1次元データ(1 voxelのスペクトル)に対する適用結果を示している。
【００３３】
図７(a)は、従来例を示すもので、人体頭部のある部分を計測してなる1 ボクセルの時系列データに対して、ゼロフィリング、1次元逆ＦＦＴ、位相補正を施して得られた1次元スペクトルを、リアル（real）成分表示したものであり、窓関数フィルタリング処理は施していない。ゼロフィリングは、図５(a)のサンプリング時間帯域がエコータイムＴＥを中心として６４点の時系列データであることから、これを２５６点の時系列データに拡張するための処理であり、６４点の時系列データの時間軸の両側にゼロを付加して、つまりゼロ詰めを行うことである。この従来例により得られたスペクトルによれば、図７(a)に示すように、残留水信号によるトランケーションアーチファクトが、コリン(Cho)、クレアチン(Cr)の信号付近まで及んでいることがわかる。
【００３４】
図７(b)も従来例を示すもので、同図(a)と同じ1 voxelの時系列データに対して、窓関数フィルタリング処理、ゼロフィリング(64点 → 256点)、1次元逆ＦＦＴ、位相補正を施して得られた1次元スペクトルをリアル成分表示したものである。フィルタリング処理の窓関数は、図６に示した関数であり、時間軸方向に余弦フィルタ関数及び線形フィルタ関数を乗算して得られる関数である。この窓関数フィルタリング処理により、トランケーションアーチファクトが低減されているが、ピークの幅が広がりChoとCrのピークが重畳してしまうことが解る。
【００３５】
図７(c)は、図１の本発明の処理を施して得られるスペクトルであり、同図(a)と同じ1ボクセルの時系列データに対して、線型予測による外挿処理、指数関数フィッティングによる正規化処理、窓関数フィルタリング処理、1次元逆FFT、位相補正を施して得られた1次元スペクトルを、リアル成分表示したものである。同図(a),(b)と比較して明らかなように、時系列データの外挿による拡張及び正規化により、トランケーションアーチファクトが低減され、かつスペクトル分解能の低下が生じていない。つまり、ChoとCrのピークが分離されていることが解る。
【００３６】
図８に、本発明のケミカルシフト情報処理法を4次元データに適用して得られるスペクトルを示す。すなわち、４次元MRSIで計測した4次元時系列データに３次元(X、Y、Z方向の3軸)の逆FFTを施した後、全voxelの時系列データに対し、図７(c)作成の際と同等の新規再構成処理を追加、つまり全ボクセルのＮ−アセチルアスパラギン酸（NAA）ピーク位置検索及びNAAピーク位置揃え処理を追加して行った。計測した8スライスエンコード分の4次元データから中央のスライスを抽出し、抽出した画像の脳内に位置する一ボクセルのスペクトルを重ね書きしたものが図８のスペクトルである。このスペクトルにおいても、ChoとCrピークがよく分離できていることがわかる。
【００３７】
上述したように、図１の実施形態によれば、代謝物質信号のスペクトル分解能を低下させることなく、残留水信号によるリンギングアーチファクトを低減させることができるため、近接した複数のピークの分離が容易となり、データの定量性を向上させることが可能となる。
【００３８】
なお、本発明の処理によりトランケーションアーチファクトの大部分を抑えることが可能となるが、位相補正後のリアル（real）成分に小さなトランケーションアーチファクトが残留する場合がある。これは、時系列データのアブソリュート（absolute）成分では始点と終点を滑らかに接続できていても、そのリアル成分及びイマジナリ−（imaginary）成分だけをみた場合、始点と終点にギャップが残留していることがあるためである。この小さなトランケーションアーチファクトを無くすためには、正規化後データのリアル成分及びイマジナリ−成分毎に、両端の数点のデータに対して線形（linear）フィルタを施すことが有効である。
【００３９】
ここで、図１のステップＳ４、Ｓ５の線型予測による外挿処理及び指数関数フィッティングによる正規化処理の具体例について説明する。まず、1次元の線型予測法について説明する。測定データをS(n), (n = 1, 2, 3, ... , M)と表す。ここでMはデータ点数である。これに対してN > Mなる整数Nをとってデータ外挿を行うことにより、s(n), (n = 1, 2, 3, ..., N)を作成する。この際、次の漸化式により外挿を行う。
【００４０】
【数１】

【００４１】
ここで、pは予測次数、xｋ (1≦k≦p)は予測係数と呼ばれるものである。この式は、あるデータが前p個のデータの線型結合で表されることを意味している。予測係数xｋは、シンギュラーヴァリューデコンポジション（Singular ValueDecomposition (SVD)）を使う方法により求める。SVDの概要を以下に説明する。
【００４２】
Bを(M - ｐ)行ｐ列の行列で、要素Bijが
【００４３】
【数２】
Bij＝ s (M - I- j + 1), (1≦j≦M - p, 1≦j≦p) (式2)
で与えられるものとし、bを(M - p)次ベクトルで、要素biが
【００４４】
【数３】
bi= s (M i +1) (式3)
で与えられるものとする。予測係数x(p次ベクトル)を求めるため次の方程式をたてる。
【００４５】
【数４】
Bx≒- b (式４)
ここで≒は最小二乗ノルムの意味で、‖Bx + b‖を最小にするxを求めることになる。ただしm次ベクトルvに対して、‖v‖とは
【００４６】
【数５】

【００４７】
で定義されるものである。
【００４８】
正規方程式にするために、Bのエルミート共益B＊を両辺に前から掛けて、
【００４９】
【数６】
B＊Bx＝−B＊ｂ (式6)
を得る。B＊BをA、 B＊bをaと書くと、
【００５０】
【数７】
Ax=−a (式7)
となる。ここでAはp行p列エルミート行列、aはp次ベクトルである。このAに対してSVDを用いれば、あるp行p列ユニタリ行列Vと、p行p列対角行列Lが存在して、
【００５１】
【数８】
A = V L V＊ (式8)
が成立する。よって、
【００５２】
【数９】
x = −V L~¹V＊a (式9)
が成立する。ただしL~¹はp行p列対角行列で、
【００５３】
【数１０】
Λ~¹ii＝１／Λ~¹ii，但し、Λii≠０又は
＝０，但し、Λii＝０、（１≦ｉ≦ｐ） (式10)
で定義されるものである。ここで，Ｌ、Vはエルミート行列Aの個有値、固有ベクトルに対応するものである。
【００５４】
次に、指数関数フィッティングによる正規化処理について説明する。外挿処理を施したデータのアブソリュート（absolute）成分の平方根をS(n), (n = 1, 2, 3, ... , N) )と表す。ここでNは外挿後のデータ点数である。S(n)の近似関数が次の式１１となる場合を考える。
【００５５】
【数１１】
F(n) = a×exp(−b×n) (式11)
式１１の両辺の自然対数をとり、式１２とおく。
【００５６】
【数１２】
log F(n) = log a − b×n (式12)
更に次の式１３に示す変換を行うことにより、
【００５７】
【数１３】
a ₀ = log a, a ₁ = b, f(n) = log F(n) (式13)
式１４の１次式に当てはめることができる。
【００５８】
【数１４】
f(n) = a ₀ − a ₁×n (式14)
式１４に対して最小二乗法を用いてフィッティングを行い、指数をとってF(n)を求める。そして、S(n)のリアル成分Re(n)と、イマジナリ−成分Im(n)に、F(n)の逆数を乗じることにより、始点と終点を滑らかに接続することが可能になる。また更に、逆FFT時の中心データとなるF(N / 2)の逆数を両成分に乗じておくことにより、代謝物質信号(T2の長い物質)のピーク面積の処理前後における変化を抑えることができる。これらの処理により最終的な両成分の式は、次の式１５、１６となる。
【００５９】
【数１５】

【００６０】
【数１６】

【００６１】
以上説明した実施形態では、MRSI測定手法としてエコープラナースペクトロスコピックイメージングを用いた。しかし、本発明はこれに限らず、一般的なMRSIの一手法であるジャーナルオブマグネティックレゾナンス(Journal of Magnetic Resonance)第61巻、第188-191頁 (1985年) 等に記載された4-Dimensional Chemical Shift Imaging (4D-CSI)法や、エコープラナースペクトロスコピックイメージングとは異なる高速MRSIの一手法であるマグネティックレゾナンスインメディスン(Magnetic Resonance in Medicine)第1巻、第370-386頁(1984年)等に記載されたマンスフィールド(Mansfield)による提案のEPSM法及びPREP法を用いた場合も、同様の効果が得られる。
【００６２】
また、上述した実施形態では、¹H(プロトン)を含む化合物を測定対象とする¹H−MRSIの場合について述べたが、水素原子以外の原子核を含む化合物からの磁気共鳴信号を観測する場合にも同様の効果が得られる。
【００６３】
【発明の効果】
以上述べたように、本発明によれば、ケミカルシフト情報のスペクトル分解能を向上させることができる。
【図面の簡単な説明】
【図１】本発明に係るケミカルシフト情報処理方法の一実施形態の処理手順を説明する図である。
【図２】本発明のケミカルシフト情報処理方法を適用してなる磁気共鳴撮影装置の全体構成図である。
【図３】本発明の説明に用いたエコープラナースペクトロスコピックイメージングの一例のパルスシーケンスを示す図である。
【図４】エコープラナースペクトロスコピックイメージングにおける信号処理方法を説明する図である。
【図５】本発明のケミカルシフト情報処理方法を説明する図であり、（ａ）は計測で得られる元時系列データの例、（ｂ）は元時系列データの時間軸の両側のデータを外挿で補完してなる拡張時系列データの例、（ｃ）は拡張時系列データを近似関数で正規化した結果の拡張時系列データの例、（ｄ）は正規化に用いた近似関数の例をそれぞれ示す。
【図６】窓関数フィルタリング処理に用いたフィルタリング関数の例を示す。
【図７】ケミカル情報処理により得られるスペクトルの例を示すものであり、（ａ）は従来法によりフィルタリング処理をしない場合の例、（ｂ）は従来法によりフィルタリング処理を行った場合の例、（ｃ）は本発明方法により得られるスペクトルの例である。
【図８】本発明方法により得られる4次元データのスペクトルの例である。
【符号の説明】
RF 高周波磁場
Ｇｚ Z軸方向の傾斜磁場
Ｇｘ X軸方向の傾斜磁場
Gｙ Y軸方向の傾斜磁場
RF１ 90°パルス
RF２ 180°パルス
Ｇｓ１スライス磁場
Ｇｓ２スライス磁場
Gｒリードアウト磁場
Gｅ、Ge1、Ge2 エンコード磁場
Sig エコートレイン信号（磁気共鳴エコー信号）
TR 繰返し時間
TE エコー時間
２静磁場発生マグネット
３傾斜磁場発生コイル
４プローブ
５計算機
６ディスプレイ
７、１２駆動用電源部
８送信器
９受信器
１０シーケンス制御装置
１１シムコイル
１３記憶媒体[0001]
BACKGROUND OF THE INVENTION
The present invention relates to a magnetic resonance imaging apparatus, and more particularly to a chemical shift information processing method for measuring a spatial distribution of information related to chemical shift.
[0002]
[Prior art]
A magnetic resonance imaging apparatus irradiates a subject placed in a static magnetic field with a high-frequency magnetic field of a specific frequency to cause a magnetic resonance phenomenon, thereby using a magnetic resonance signal generated from the subject to For example, an image is created based on the measurement of chemical and physical information of a substance constituting the specimen. In particular, Magnetic Resonance Spectroscopic Imaging (hereinafter abbreviated as MRSI) measures the difference in magnetic resonance frequency due to the difference in chemical bonds of various molecules as chemical shift information, and spatial distribution for each molecule. This method is known as a method for obtaining an image as a chemical shift image.
[0003]
Information processing methods for acquiring chemical shift images from magnetic resonance signals generally include a series of magnetic resonance signals (echo train signals) measured by adding spatial information encoded in frequency and phase, including the chemical shift axis. After arranging for each measurement axis, inverse Fourier transform is performed for each axis. That is, for example, when the space to be measured is a four-dimensional real space of chemical shift axis (δ = time axis), X-axis, Y-axis, and Z-axis, each axis of the four-dimensional Fourier space (kδ axis, kx axis, ky axis, and kz axis) correspond to measurement axes. In order to improve the image quality and spectral accuracy of the chemical shift image, pre-processing such as window function filtering and zero filling is performed before inverse Fourier transform, and post-processing such as phase correction is performed after inverse Fourier transform. That is done as appropriate.
[0004]
[Problems to be solved by the invention]
By the way, proton ( ¹ H) ¹ In H-MRSI, in order to detect a signal from a very small amount of a metabolite, it is necessary to suppress the generation of a water signal or a subcutaneous fat signal having a very large signal level relative to the signal level of the metabolite. For example, the water signal reaches several thousand times the metabolite signal, and the subcutaneous fat signal reaches several hundred times the metabolite signal. In particular, suppression of the water signal is important, and usually three Gaussian high-frequency magnetic field pulses with a limited excitation band (with the band limited to about 1 to 2 ppm and the action of selectively exciting only the water signal) A water signal suppression sequence comprising a spoiler gradient magnetic field (having the effect of suppressing the generation of signals by causing the phases of excited magnetization to fall apart) is used. With this processing, the water signal intensity can be suppressed to about one-hundredth of the signal intensity at the time of non-suppression, but the water signal that could not be suppressed remains as a residual water signal, which may cause truncation artifacts.
[0005]
Here, the reason why the generation of truncation artifacts due to the residual water signal becomes a problem will be described. All magnetic resonance signals generated by the excitation of nuclear magnetization decay with time. This process is called relaxation. Relaxation includes relaxation of longitudinal magnetization (T1 relaxation) and relaxation of transverse magnetization. Transverse magnetization relaxation is further divided into T2 * relaxation due to static magnetic field inhomogeneity and T2 relaxation due to the internal state of nuclear magnetization. When the static magnetic field inhomogeneity is large and 180 ° pulse irradiation is applied, T2 * relaxation is dominant, and an echo train signal with a bell-shaped envelope having a maximum echo time (TE) is generated. On the other hand, in the ideal state with no static magnetic field inhomogeneity, signal attenuation during one sampling period is governed only by T2 relaxation, and is always exponential (exp (−t / T2) regardless of the presence or absence of a 180 ° pulse. Where t is the sampling time and T2 is the transverse relaxation time, ie the envelope is an exponential echo train signal, so the signal is usually a combination of bell and exponential types However, in general, when MRSI is executed, in order to improve the spectral resolution and improve the quality of metabolite images, the uniformity of the static magnetic field is often improved by auto shimming etc. In this case, T2 relaxation is dominant. Thus, it approaches an exponential signal shape.
[0006]
When inverse FFT is performed on an exponential signal, the signal level at the start and end of the signal is different, so the Gibbs phenomenon due to discontinuity (gap) (periodic increase or decrease near the discontinuity), that is, Truncation artifacts will occur. In addition, the signal strength of this truncation artifact is proportional to the gap size, so truncation artifacts due to high signal strength and short T2 residual water signal (T2 ≈ 100 ms) The effect on (T2 ≈ 250 to 400 ms) is a very big problem.
[0007]
In order to reduce the truncation artifact, it is effective to perform the window function filtering process as described above and smoothly connect the start point and the end point of the signal. However, since this filtering process has a broadening effect that increases not only the residual water signal but also the peak width of the metabolite signal, there is a problem of reducing the spectral resolution. In other words, the window function filtering process has an effect of transforming the exponential signal shape into a bell shape, and this effect has the same effect on the spectrum as when the static magnetic field uniformity is lowered.
[0008]
The problem to be solved by the present invention is to improve the spectral resolution of chemical shift information.
[0009]
[Means for Solving the Problems]
The above problem is that an echo train signal to which encode information measured by magnetic resonance imaging is added is captured, the echo train signal is arranged along the chemical shift axis, and then the inverse Fourier transform in the axial direction related to the encode information is performed. Based on the time series data obtained by the inverse Fourier transform, the time series data is supplemented by extrapolation over a predetermined range in the time axis direction to generate extended time series data. A chemical shift information processing method that obtains chemical shift information by normalizing by dividing by the approximate function of the extended time series data and performing inverse Fourier transform in the chemical shift axis direction based on the normalized extended time series data This can be solved.
[0010]
That is, before obtaining chemical shift information by inverse Fourier transform in the chemical shift axial direction, time series data obtained by applying the inverse Fourier transform in the axial direction related to spatial information is extrapolated to the time zone that is not measured. Complement data. By normalizing the extended time-series data obtained in this way, the signal levels at the start and end points of the extended time-series data are made equal, that is, the attenuation due to T2 relaxation is compensated, and the discontinuity (gap) is reduced. It can be greatly eased or reduced. As a result, truncation artifacts generated by inverse Fourier transform in the chemical shift axis direction can be reduced, and the spectral resolution of chemical shift information can be improved.
[0011]
In addition to the above processing, the normalized extended time series data is filtered by a window function filter including the time range of the time series data, and the chemical shift axial direction is based on the filtered extended time series data. It is preferable to obtain the chemical shift information by performing the inverse Fourier transform. According to this, the spectral resolution of chemical shift information can be further improved.
[0012]
A magnetic resonance imaging apparatus according to the present invention includes a static magnetic field generation means, a gradient magnetic field generation means for generating gradient magnetic fields in three different directions, a high-frequency magnetic field generation means, a signal detection means for detecting a magnetic resonance signal, Comprising arithmetic means for performing processing arithmetic based on the magnetic resonance signal obtained by the signal detection means, storage means, and sequence control means for controlling the operation of each means, the arithmetic means includes the above-mentioned This can be realized by incorporating a chemical shift information processing function.
[0013]
In addition, the chemical shift information processing according to the present invention captures a series of magnetic resonance signals (echo train signals) measured by a magnetic resonance imaging apparatus into a computer provided separately via a signal transmission path or a recording medium, Can be processed by computer.
[0014]
A linear prediction method can be used when extrapolating unmeasured time series data. A function that approximates time-series data is an exponential function, and the least square method can be used when obtaining this approximate exponential function. The control means preferably performs control to detect the magnetic resonance signal while applying the gradient polarity of the gradient magnetic field in at least one direction out of the three different gradient magnetic fields.
[0015]
[Embodiment]
Hereinafter, embodiments of the present invention will be described with reference to the drawings. FIG. 1 shows the procedure of the chemical shift information processing method of the present invention, and FIG. 2 shows the overall configuration of a magnetic resonance imaging apparatus of one embodiment to which the chemical shift information processing method of the present invention is applied. FIG. 3 shows an example of a measurement sequence to which the present invention can be applied.
[0016]
In the magnetic resonance imaging apparatus of FIG. 2, the subject 1 is placed in a space to which a static magnetic field generated by the static magnetic field generating magnet 2 and a gradient magnetic field in three different directions generated by the gradient magnetic field generating coil 3 are applied. The shim coil 11 includes a plurality of coils, and the uniformity of the static magnetic field can be adjusted by changing the current flowing through each coil. However, the shim coil 11 is not necessarily provided. The probe 4 generates a high-frequency magnetic field that irradiates the subject 1 to cause a magnetic resonance phenomenon in the subject 1 and detects the magnetic resonance signal generated from the subject 1 by the probe 4. The computer 5 processes the magnetic resonance signal to generate necessary information, such as a chemical shift spectrum and a chemical shift image, and displays them on the display 6. The sequence power supply unit 12 for the shim coil 11, the drive power unit 7 for the gradient magnetic field generating coil 3, the transmitter 8 and the receiver 9 are controlled by the sequence control device 10. In addition, measurement signals, measurement conditions, and the like are stored in the storage medium 13 as necessary.
[0017]
The pulse sequence shown in FIG. 3 is echo planar spectroscopic imaging, which is one method of high-speed MRSI that can obtain a magnetic resonance spectrum image at high speed. This is a proposal proposed by Matsui et al. It is an example. In the figure, in order from the top, the excitation high-frequency magnetic field pulse RF, the gradient magnetic field Gz in the Z-axis direction, the gradient magnetic field Gx in the X-axis direction, the gradient magnetic field Gy in the Y-axis direction, and the magnetic resonance signal Sig are shown. In this pulse sequence, an excitation high-frequency magnetic field pulse RF1 (90 ° pulse) and a slice (volume) selection gradient magnetic field Gs1 in the Z-axis direction are applied simultaneously. As a result, the nuclear spins included in the region to be measured are selectively excited. Next, by applying phase encoding gradient magnetic fields Ge1 and Ge2 in the Y-axis and Z-axis directions to the magnetization generated by excitation, spatial information in the Y-axis and Z-axis directions is given to the phases of the respective magnetizations. Further, TE / 2 after applying the excitation pulse RF1 (TE: echo time), the inversion high-frequency magnetic field pulse RF2 (180 ° pulse) and the slice (volume) selection gradient magnetic field Gs2 in the Z-axis direction are simultaneously applied. As a result, the nuclear spins whose phases are shifted from each other due to the static magnetic field inhomogeneity are converged again, and a spin echo signal is generated after TE / 2 from the application of RF2. When measuring this spin echo signal, the gradient of the readout gradient magnetic field Gr in the X-axis direction is periodically reversed and applied, thereby providing an echo train that includes the chemical shift and spatial information in the X-axis direction. A signal can be generated. A trapezoidal wave or a sine wave (sin wave) is used as the waveform of Gr. The illustrated example is an example in which a trapezoidal wave is used for Gr, and the polarity switching unit is not smooth, but may be switched so as to shift continuously.
[0018]
This series of pulse sequences is changed stepwise as indicated by broken lines in the drawing, and the Ge1 and Ge2 application intensities are repeated at intervals of TR by the number of times corresponding to the number of pixels in the Y-axis and Z-axis directions. For example, in order to obtain an image of the number of spatial matrices mx × my × mz, it is necessary to repeat excitation and measurement mx × my × mz times in the most widely used conventional method (4D-CSI). However, the high-speed MRSI in FIG. 3 only requires my × mz excitation / measurement, so the measurement time can be reduced to 1 / mx.
[0019]
Next, a data processing method in this high-speed MRSI will be described. The echo train signal obtained by the above measurement is frequency-modulated by both the external magnetic interaction due to the readout gradient magnetic field Gr and the internal magnetic interaction corresponding to the chemical shift. That is, when attention is paid to each echo signal one by one, the spatial information in the readout axis direction is the same as each echo signal obtained by ultra-high-speed imaging method EPI (Echo Planar Imaging) that is widely used in MRI at present. It can be seen that is frequency encoded. Further, when attention is paid to the center points of the echo signals, the influence of Gr is canceled because the gradient polarity of the gradient magnetic field of Gr is periodically reversed and applied at these points. For this reason, the envelope of the echo train signal connecting the center points of the echo signals is modulated only by the chemical shift. Therefore, by performing the data processing described below, it is possible to extract chemical shift information and spatial information in the X-axis direction from the echo train signal.
[0020]
That is, as shown in FIG. 4, one phase-modulated echo train signal E1, E2, E3, E4, E5, E6,... Is divided into echoes and grouped into odd-numbered echoes and even-numbered echoes. After that, two-dimensional rearrangement is performed along the time axis kδ axis. The reason for grouping into odd-numbered and even-numbered groups is to reduce artifacts by performing processing for each echo having the same Gr gradient polarity. In the rearranged two-dimensional data array, the time change of the signal along the newly defined time axis kδ axis is determined purely by the internal magnetic interaction corresponding to the chemical shift. This is because the external magnetic interaction due to Gr disappears along the kδ axis due to the periodic inversion of Gr. On the other hand, the shape of each echo signal along the kx axis is governed by the external magnetic interaction due to Gr. The spatial information in the Y-axis and Z-axis directions is Ge ₁ And Ge ₂ Thus, phase encoding is performed along the ky and kz axes perpendicular to the kx axis. By acquiring a series of echo signals in this way, a four-dimensional Fourier space (kδ axis, kx axis, connected by a Fourier transform with a four-dimensional real space (chemical shift δ axis, X axis, Y axis, and Z axis), ky axis and kz axis) can be measured.
[0021]
Here, the number of matrices on each axis is determined by the number of sampling points per echo in the readout axis (X) direction, as in EPI, and the number of spectral points in the chemical shift axis direction is the echo contained in one echo train. The number of pixels in the phase encoding axis direction is determined by Ge ₁ And Ge ₂ Determined by the number of steps. A chemical shift image can be obtained by performing a four-dimensional inverse Fourier transform on these four-dimensional data arrays. Two chemical shift images created from the odd-numbered and even-numbered echo groups are added after phase correction. When sine wave drive is used for Gr, the trajectory for scanning the Fourier space has a sine wave shape, and it is necessary to perform a four-dimensional inverse Fourier transform including sine correction.
[0022]
Next, the chemical shift information processing method according to the features of the present invention will be described in detail. The computer 5 takes in an echo train signal Sig in which spatial information measured by magnetic resonance imaging is encoded, and arranges the echo train signal along the chemical shift axis δ, and then the axial direction (X, Y, Z) and the chemical shift axis direction are subjected to inverse Fourier transform to obtain chemical shift information. At this time, when the window function filtering process is directly performed on the time series data (original data) obtained by the inverse Fourier transform, as described above, the broadening effect given to the metabolite signal is too large.
[0023]
Therefore, as shown in FIG. 1, the present invention performs inverse Fourier transform in the axial direction (X, Y, Z) related to spatial information encoded in frequency and phase before performing inverse Fourier transform in the chemical shift axial direction. The time series data obtained in this way is supplemented by extrapolation over a predetermined range in the time axis direction to generate extended time series data, and the extended time series data obtained thereby is divided by its approximate function and normalized, Chemical shift information is obtained by performing inverse Fourier transform in the chemical shift axis direction on the normalized extended time series data.
[0024]
(Step S1)
The pulse sequence of FIG. 3 is performed to measure the four-dimensional MRSI, and an echo train signal Sig (magnetic resonance signal) to which spatial information frequency-encoded and phase-encoded is added is collected.
[0025]
(Step S2)
The echo train signals are arranged along the chemical shift axis δ and rearranged into a four-dimensional data array.
[0026]
(Step S3)
Time series data is obtained by performing an inverse Fourier transform in the axial direction (X, Y, Z) related to each spatial information.
[0027]
(Step S4)
Extrapolation processing using linear prediction is performed on the time-series data of each voxel obtained in step S3. For example, with respect to the time-series data (original data) of the one-dimensional data (1 voxel) shown in FIG. 5A, data at a time earlier than the original data acquisition time and data at a later time than the original data acquisition time are By extrapolating from the data and predicting, the expanded time series data shown in FIG. 5B is generated.
[0028]
(Step S5) An exponential function fitting process is performed on the expanded time series data of each voxel obtained in Step S4. That is, an approximate function that approximates the attenuation of the extended time series data is obtained. This approximate function (fitting function ) Is, for example, an exponential function shown in FIG.
[0029]
(Step S6)
Normalization processing is performed by multiplying the extended time series data by the reciprocal of an exponential function that is an approximate function. That is, the extended time series data is divided by the approximate function to perform processing for compensating (normalizing) attenuation of T2 relaxation. Thereby, for example, the normalized extended time series data shown in FIG. 5C is obtained.
[0030]
(Step S7)
A window function filtering process is performed on the extended time-series data of each normalized voxel. The filtering function used in the window function filtering process has a pattern that selectively attenuates the signal intensity in the vicinity of the start point and end point of the sampling time axis, as shown in FIG.
[0031]
(Step S8)
An inverse Fourier transform (FFT) process in the δ-axis direction is performed on the expanded time series data of each filtered voxel. Thereby, spectrum data of chemical shift information is obtained.
[0032]
Spectral data of chemical shift information obtained by performing the chemical shift information processing of FIG. 1 will be described with reference to FIG. 7 while comparing with spectral data obtained by a conventional method. FIG. 7 shows the result of application to one-dimensional data (1 voxel spectrum) as in FIG.
[0033]
FIG. 7 (a) shows a conventional example, which is obtained by performing zero filling, one-dimensional inverse FFT, and phase correction on time series data of one voxel obtained by measuring a part of the human head. The one-dimensional spectrum is displayed as a real component and is not subjected to window function filtering. The zero filling is a process for extending the sampling time band of FIG. 5 (a) to 64 points of time series data centering on the echo time TE, so that 64 points. In other words, zeros are added to both sides of the time axis of the time series data. According to the spectrum obtained by this conventional example, as shown in FIG. 7A, it can be seen that the truncation artifact due to the residual water signal extends to the vicinity of the signals of choline (Cho) and creatine (Cr).
[0034]
FIG. 7 (b) also shows a conventional example. For the same 1 voxel time-series data as in FIG. 7 (a), window function filtering processing, zero filling (64 points → 256 points), one-dimensional inverse FFT, A one-dimensional spectrum obtained by performing phase correction is displayed as a real component. The window function of the filtering process is the function shown in FIG. 6, and is a function obtained by multiplying the cosine filter function and the linear filter function in the time axis direction. Although this truncation artifact is reduced by this window function filtering processing, it can be seen that the peak width is widened and the peaks of Cho and Cr are superimposed.
[0035]
FIG. 7 (c) is a spectrum obtained by applying the processing of the present invention of FIG. 1, and extrapolation processing and exponential function fitting by linear prediction are performed on the same time series data of 1 voxel as FIG. The one-dimensional spectrum obtained by performing normalization processing, window function filtering processing, one-dimensional inverse FFT, and phase correction is displayed as a real component. As is clear from comparison with FIGS. 9A and 9B, the truncation artifacts are reduced and the spectral resolution is not reduced by the extension and normalization by extrapolation of the time series data. That is, it can be seen that the peaks of Cho and Cr are separated.
[0036]
FIG. 8 shows a spectrum obtained by applying the chemical shift information processing method of the present invention to four-dimensional data. In other words, three-dimensional (three axes in the X, Y, and Z directions) inverse FFT is applied to the four-dimensional time series data measured by the four-dimensional MRSI, and then the time series data for all voxels is created as shown in FIG. A new reconstruction process equivalent to that in the above was added, that is, N-acetylaspartic acid (NAA) peak position search and NAA peak alignment process of all voxels were added. The spectrum shown in FIG. 8 is obtained by extracting the center slice from the measured four-slice encoded 4-slice data and overwriting the spectrum of one voxel located in the brain of the extracted image. Also in this spectrum, it can be seen that the Cho and Cr peaks are well separated.
[0037]
As described above, according to the embodiment of FIG. 1, ringing artifacts due to the residual water signal can be reduced without reducing the spectral resolution of the metabolite signal, so that it is easy to separate a plurality of adjacent peaks. It becomes possible to improve the quantitativeness of data.
[0038]
Although the processing of the present invention can suppress most of the truncation artifacts, small truncation artifacts may remain in the real component after phase correction. This is because even if the absolute component of the time series data can smoothly connect the start point and the end point, when only the real component and the imaginary component are seen, a gap remains at the start point and the end point. Because there are things. In order to eliminate this small truncation artifact, it is effective to apply a linear filter to several data at both ends for each real component and imaginary component of the normalized data.
[0039]
Here, specific examples of extrapolation processing by linear prediction and normalization processing by exponential function fitting in steps S4 and S5 in FIG. 1 will be described. First, a one-dimensional linear prediction method will be described. The measured data is represented as S (n), (n = 1, 2, 3, ..., M). Here, M is the number of data points. On the other hand, s (n), (n = 1, 2, 3,..., N) is created by extrapolating data by taking an integer N where N> M. At this time, extrapolation is performed by the following recurrence formula.
[0040]
[Expression 1]

[0041]
Here, p is a prediction order, and xk (1 ≦ k ≦ p) is a prediction coefficient. this The expression means that certain data is represented by a linear combination of the previous p pieces of data. The prediction coefficient xk is obtained by a method using a Singular Value Decomposition (SVD). The outline of SVD is explained below.
[0042]
B is a matrix of (M-p) rows and p columns, and element Bij is
[0043]
[Expression 2]
Bij = s (M-I- j + 1), (1 ≤ j ≤ M-p, 1 ≤ j ≤ p) (Formula 2)
Where b is an (M-p) degree vector and the element bi is
[0044]
[Equation 3]
bi = s (M i +1) (Equation 3)
It shall be given by In order to obtain the prediction coefficient x (p-th order vector), the following equation is established.
[0045]
[Expression 4]
Bx ≒ -b (Formula 4)
Here, ≒ means the least square norm, and x that minimizes ‖Bx + b‖ is obtained. However, for m-order vector v, ‖v‖ is
[0046]
[Equation 5]

[0047]
Is defined by
[0048]
To make a normal equation, multiply both sides by Hermitian common benefit B * from the front,
[0049]
[Formula 6]
B * Bx = -B * b (Formula 6)
Get. If you write B * B as A and B * b as a,
[0050]
[Expression 7]
Ax = −a (Formula 7)
It becomes. Here, A is a p-row p-column Hermitian matrix, and a is a p-order vector. If SVD is used for this A, there exists a p-row p-column unitary matrix V and a p-row p-column diagonal matrix L,
[0051]
[Equation 8]
A = VLV * (Formula 8)
Is established. Therefore,
[0052]
[Equation 9]
x = −VL ~ ¹ V * a (Formula 9)
Is established. L ~ ¹ Is a p-by-p diagonal matrix,
[0053]
[Expression 10]
Λ ~ ¹ ii = 1 / Λ ~ ¹ ii, where Λii ≠ 0 or
= 0, Λii = 0, (1 ≦ i ≦ p) (Equation 10)
Is defined by Here, L and V correspond to the unique values and eigenvectors of the Hermitian matrix A.
[0054]
Next, normalization processing by exponential function fitting will be described. The square root of the absolute component of the extrapolated data is represented as S (n), (n = 1, 2, 3, ..., N)). Here, N is the number of data points after extrapolation. Consider the case where the approximate function of S (n) is given by the following equation (11).
[0055]
## EQU11 ##
F (n) = a × exp (−b × n) (Equation 11)
Taking the natural logarithm of both sides of Equation 11, Equation 12 is given.
[0056]
[Expression 12]
log F (n) = log a − b × n (Equation 12)
Furthermore, by performing the conversion shown in the following equation 13,
[0057]
[Formula 13]
a ₀ = log a, a ₁ = b, f (n) = log F (n) (Equation 13)
It can be applied to the linear expression of Equation 14.
[0058]
[Expression 14]
f (n) = a ₀ − A ₁ × n (Formula 14)
Fitting is performed on the equation 14 using the least square method, and an exponent is taken to obtain F (n). Then, by multiplying the real component Re (n) of S (n) and the imaginary component Im (n) by the reciprocal of F (n), the start point and the end point can be smoothly connected. Furthermore, by multiplying both components by the reciprocal of F (N / 2), which is the central data during inverse FFT, it is possible to suppress changes in the peak area of the metabolite signal (a substance with a long T2) before and after processing. it can. By these processes, the final formulas of both components become the following formulas 15 and 16.
[0059]
[Expression 15]

[0060]
[Expression 16]

[0061]
In the embodiment described above, echo planar spectroscopic imaging is used as the MRSI measurement method. However, the present invention is not limited to this, and 4-Dimensional described in Journal of Magnetic Resonance, Vol. 61, pp. 188-191 (1985), which is a general MRSI technique. Magnetic Resonance in Medicine, Volume 1, pages 370-386 (1984) Similar effects can be obtained by using the EPSM method and the PREP method proposed by Mansfield described in the above.
[0062]
In the embodiment described above, ¹ Measuring compounds containing H (proton) ¹ Although the case of H-MRSI has been described, the same effect can be obtained when observing magnetic resonance signals from compounds containing nuclei other than hydrogen atoms.
[0063]
【The invention's effect】
As described above, according to the present invention, the spectral resolution of chemical shift information can be improved.
[Brief description of the drawings]
FIG. 1 is a diagram illustrating a processing procedure of an embodiment of a chemical shift information processing method according to the present invention.
FIG. 2 is an overall configuration diagram of a magnetic resonance imaging apparatus to which the chemical shift information processing method of the present invention is applied.
FIG. 3 is a diagram showing a pulse sequence of an example of echo planar spectroscopic imaging used for explaining the present invention.
[Figure 4] Echope Lana It is a figure explaining the signal processing method in spectroscopic imaging.
5A and 5B are diagrams for explaining the chemical shift information processing method of the present invention, in which FIG. 5A shows an example of original time series data obtained by measurement, and FIG. 5B shows data on both sides of the time axis of the original time series data. Example of extended time series data complemented by extrapolation, (c) is an example of extended time series data obtained by normalizing extended time series data with an approximation function, and (d) is an approximation function used for normalization. Each example is shown.
FIG. 6 shows an example of a filtering function used in the window function filtering process.
FIGS. 7A and 7B show examples of spectra obtained by chemical information processing, where FIG. 7A shows an example when filtering is not performed by a conventional method, and FIG. 7B shows an example when filtering is performed by a conventional method; (C) is an example of a spectrum obtained by the method of the present invention.
FIG. 8 is an example of a spectrum of four-dimensional data obtained by the method of the present invention.
[Explanation of symbols]
RF high frequency magnetic field
Gz Z-axis gradient magnetic field
Gx X-axis gradient magnetic field
Gy Y-axis gradient magnetic field
RF1 90 ° pulse
RF2 180 ° pulse
Gs1 slice magnetic field
Gs2 slice magnetic field
Gr lead-out magnetic field
Ge, Ge1, Ge2 encoding magnetic field
Sig echo train signal (magnetic resonance echo signal)
TR repetition time
TE echo time
2 Static magnetic field generation magnet
3 Gradient magnetic field generating coil
4 Probe
5 Calculator
6 Display
7, 12 Power supply for driving
8 Transmitter
9 Receiver
10 Sequence control device
11 Shim coil
13 Storage media

Claims

Static magnetic field generating means, gradient magnetic field generating means for generating gradient magnetic fields in three different directions, high frequency magnetic field generating means, signal detecting means for detecting magnetic resonance signals, and magnetic resonance signals obtained by the signal detecting means In a magnetic resonance imaging apparatus including a calculation unit that performs processing calculation based on the above, a storage unit, and a sequence control unit that controls the operation of each unit,
The calculation means takes in an echo train signal in which position information measured by magnetic resonance imaging is encoded, arranges the echo train signal along a chemical shift axis, and then reverses the axial direction related to the position information. Based on the time series data obtained by performing Fourier transform and the inverse Fourier transform, the time series data is supplemented by extrapolation over a predetermined range in the time axis direction to generate expanded time series data, Dividing the series data by the approximate function of the extended time series data to normalize, and applying the inverse Fourier transform in the chemical shift axis direction based on the normalized extended time series data to obtain chemical shift information Magnetic resonance imaging device.

Different from static magnetic field generation means Three Gradient magnetic field generating means for generating a gradient magnetic field in the direction, high frequency magnetic field generating means, signal detecting means for detecting a magnetic resonance signal, and arithmetic means for performing a processing operation based on the magnetic resonance signal obtained by the signal detecting means And a magnetic resonance imaging apparatus including storage means and sequence control means for controlling the operation of each means.
The computing means generates extended time series data by supplementing the measured time series data by extrapolation over a predetermined range in the time axis direction, and generates the extended time series data. Magnetic resonance imaging, wherein normalization is performed by dividing by an approximate function of time series data, and chemical shift information is obtained by performing inverse Fourier transform in the chemical shift axis direction based on the normalized extended time series data apparatus.

The magnetic resonance imaging apparatus according to claim 1 or 2,
The calculation means filters the normalized extended time series data with a window function filter including a time range of the time series data, and the filtered extended time series data is subjected to inverse Fourier in the chemical shift axis direction. A magnetic resonance imaging apparatus characterized by performing conversion.