(57)【要約】【課題】  ベースライン補正された周波数領域のＭＲＳデータ解析のめのフィッティング処理時間を短縮しかつフィッティング精度を向上したＭＲＳ計測および計測されたＮＭＲのスペクトルデータを処理，解析する機能を備えたＭＲＩ装置を提供する。【解決手段】ベースライン補正されたＭＲＳのデータに対して，判別分析法（多変量解析）を適用し，周波数領域のＭＲＳデータをノイズ領域と，信号領域に分離する。この結果，周波数領域は，幾つかの信号領域と，ノイズ領域に分割される。フィッティングを行なう領域を，信号領域のみとし，さらに，モデル関数の和に参加する関数は，共鳴周波数がフィッティングを行なう信号領域に存在するガウス関数または，ローレンツ関数だけに制限する。この制限により，フィッティングするパラメータの数を少なくし，フィッティング処理の時間を短縮する。【選択図】  図２
2008048878.tif 000002
【特許請求の範囲】
【請求項１】
被検体の特定領域についての磁気共鳴スペクトロスコピ（以下ＭＲＳと称す）計測から得られる核磁気共鳴（以下ＮＭＲと称す）信号データを処理，解析するＭＲＳデータ処理部を備えた磁気共鳴イメージング装置において，上記ＭＲＳデータ処理部は，ＭＲＳ計測から得られた時間領域のＮＭＲ信号データを周波数領域のＮＭＲ信号にフーリエ変換した後、特定領域に含まれる水および巨大分子の寄与によるベースライン分を差し引いて得られた局所データであるベースライン補正されたＮＭＲ信号データに判別分析法を適用してノイズ領域と信号領域に区分する機能および区分されたそれぞれの信号領域についてのみその領域に共鳴周波数を持つ特定領域中の検出可能な物質の共鳴周波数をピークに持つガウス関数またはローレンツ関数の和から構成されるモデル関数を生成し，その領域の局所データとのフィッティングを行ない，モデル関数のパラメータを決定し，決定されたパラメータに基づいてそれぞれの信号領域がカバーしている特定領域中の検出可能な物質の存在比を演算する機能を有することを特徴とする磁気共鳴イメージング装置。
【請求項２】
上記ＭＲＳデータ処理部はさらにそれぞれの信号領域のデータについてのフィッティングに先立って，それぞれの信号領域の局所データについてスムーシングを行なう機能を有していることを特徴とする請求項１記載の磁気共鳴イメージング装置。
【請求項３】
上記スムーシングは，信号領域のデータを，グラフにし，このグラフを画像と見做して，画像の拡大縮小を行い，その結果滑らかになった画像を，グラフに戻すことで，周波数領域の画像をスムーシングすることを特徴とする請求項２記載の磁気共鳴イメージング装置。
【請求項４】
それぞれの信号領域の局所データとそれに対応するモデル関数についての上記フィッティングは最小自乗法あるいは最大の相互情報量を利用して行なわれることを特徴とする請求項１から３のいづれかに記載の磁気共鳴イメージング装置。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は，磁気共鳴スペクトロスコピ（以下ＭＲＳと称す）計測および計測された核磁気共鳴（以下ＮＭＲと称す）信号のスペクトルデータを処理，解析する機能を備えた磁気共鳴イメージング（以下ＭＲＩと称す）装置の改良に関する。
【背景技術】
【０００２】
ＭＲＳ計測では，ＮＭＲ信号を収集したい，例えば，生体内の特定の領域のみをＲＦパルスで励起しＮＭＲ信号を収集する。この収集したＮＭＲ信号は，各時刻における種々のＮＭＲ信号の信号値をまとめたデータから構成されている。
【０００３】
このデータは，時間領域のデータと呼ばれる。時間領域のデータは，一次元フーリエ変換を行なうことで，周波数におけるＮＭＲ信号のデータとなる。このデータは周波数領域のデータと呼ばれる。
【０００４】
通常のＭＲＳ計測では，水からのＮＭＲ信号は抑制されて計測されているが，完全には水からのＮＭＲ信号を抑制することはできない，また，生体内に存在する，例えば，たんぱく質，脂質，アミノ酸等の巨大分子からのＮＭＲ信号は，広い周波数領域に渡って分布する。
【０００５】
従って，水の信号成分を除去し，更に，生体内巨大分子からの信号の寄与からなるベースラインを決定し，このデータを原データから引くことで，ベースライン補正を行なう。
【０００６】
このベースライン補正されたデータは，ＲＦパルスで励起した領域に存在する検出可能な物質，例えば，ＮＡＡ，コリン等の代謝物質の共鳴周波数に，その存在量に比例した振幅を持つデータとなる。従って，このベースライン補正されたデータを，検出可能な物質の共鳴周波数をピーク位置にもつ，ガウス関数またはローレンツ関数の和から構成されるモデル関数で，フィティングを行い，各モデル関数のパラメータを算出すると，各物質の存在比が算出される。
【０００７】
頭部に話を限ると，正常な神経細胞は，ＮＡＡを多量に含むが，異常な状態の神経細胞は，ＮＡＡの存在量が低下する。また，脳腫瘍の場合には，細胞膜の前駆体であるコリンが大量に含まれることになる。このことを利用し，ＮＡＡとコリンの存在比を調べると，脳腫瘍の悪性度を判定することができる。
【０００８】
検出可能な物質の数に対応したガウス関数または，ローレンツ関数を準備し，この関数の和からなるモデル関数を用いてフィッティングを行なう場合，決定しなければならないパラメータの数が多くなってしまい，フィッティングの処理に時間がかかる。
【０００９】
このようなフィッティング処理にかかる時間を短縮するためには，フィティングに用いるデータの領域と，また，モデル関数に加わるガウス関数または，ローレンツ関数の数を少なくする必要がある。そのため，例えば，特開２０００−１４７０８１（特許文献１）は，検出可能な物質の周波数情報を用いて予め計算機システムにカーブフィッティング範囲を記憶させておき，カーブフィッティングの計算の範囲を限定することによって，カーブフィッティング処理時間の短縮を計った磁気共鳴装置を開示している。
【発明の開示】
【発明が解決しようとする課題】
【００１０】
しかし，この磁気共鳴装置で設定されるカーブフィッティングの範囲は，実際のＭＲＳ計測で得られるＮＭＲ信号データとは，独立して事前に設定されたものなので，事前に設定されたある範囲について計測されたＮＭＲ信号データ中に大量なノイズが含まれていた場合でもカーブフィッティングの処理は行なわれるため，それは処理の時間を引き延ばすと共に計算の精度，即ちＭＲＳ解析の精度を低減する恐れが有った。
【００１１】
本発明の目的は，ベースライン補正された周波数領域のＭＲＳデータ解析のためのフィッティング処理時間を短縮しかつフィッティング精度を向上したＭＲＳ計測および計測されたＮＭＲのスペクトルデータを処理，解析する機能を備えたＭＲＩ装置を提供することである。
【課題を解決するための手段】
【００１２】
本発明では，ベースライン補正されたＭＲＳのデータに対して，判別分析法（多変量解析）を適用し，周波数領域をノイズ領域と，信号領域に分離する。この結果，周波数領域は，幾つかの信号領域と，ノイズ領域に分割される。フィッティングを行なう領域を，信号領域のみとし，さらに，モデル関数の和に参加する関数は，共鳴周波数がフィッティングを行なう信号領域に存在するガウス関数または，ローレンツ関数だけに制限する。この制限により，フィッティング時の未知数を減らし，また，ノイズに対応する信号を計算から除外するので，高速に，高精度のフィッティング処理が実行できる。
【発明の効果】
【００１３】
本発明によれば，周波数領域のＭＲＳのデータについて，ノイズ領域を除いた領域で，少ないパラメータに対するパラメータフィッティングを行なうので，より速くまたより高い精度でＭＲＳのデータの解析が行なえる。
【発明を実施するための最良の形態】
【００１４】
以下図面を参照して本発明の一実施例であるＭＲＩ装置を説明する。
【００１５】
図１は本発明の一実施例であるＭＲＩ装置をブロック図で示したもので，被検体等を受けいれることが出来る均一静磁場の検査空間Ｐを形成するように空間Ｐを挟んで上下に対向して配置された一対の静磁場コイル１，静磁場コイル１の内側に同様に上下に対向して配置された一対の傾斜磁場コイル２，傾斜磁場コイル２の内側に配置され検査空間Ｐ内に置かれる被検体に向けて励起用の高周波磁場を照射しそれによって励起された被検体からのＮＭＲ信号を受信するＲＦコイル３，傾斜磁場コイル２に傾斜磁場形成のためのパルス電流を供給する傾斜磁場コイル電源４，ＲＦコイル３に励起用高周波磁場発生のための高周波電流を供給するための送信部５，ＲＦコイル３が受信した微弱なＮＭＲ信号を増幅，検波する受信部６，受信部６で増幅，検波されたＮＭＲ信号を収集するデータ収集部７，ＭＲＳ計測のための，例えば，ＳＴＥＡＭ法あるいはＰＲＥＳＳ法等によるパルスシーケンスあるいはＭＲＩのためのパルスシーケンスに従って高周波磁場印加，傾斜磁場印加およびデータ処理が行なわれるように，傾斜磁場コイル電源４，送信部５，受信部６およびデータ収集部７を制御するシーケンサ８，データ収集部７が収集したＭＲＳデータを処理，解析して生体内の特定領域における目的とする検出対象の代謝物質の存在比を演算するＭＲＳデータ処理部９，ＭＲＳデータ処理部９で演算された結果等を表示する表示部１０，シーケンサ８，ＭＲＳデータ処理部９，表示部１０等の動作を統括するＣＰＵ１１，ＣＰＵ１１に接続され，検者によって検査対象の特徴領域の位置の設定，検出対象の代謝物質の特定，その他ＭＲＳ計測に必要な種々のパラメータが入力される入力部１２，およびシーケンサ８，ＭＲＳデータ処理部９，表示部１０およびＣＰＵ１１による種々の処理動作に必要な種々のプログラムを格納したＲＯＭおよびＲＡＭ等を有する記憶部１３から構成されている。
【００１６】
ＭＲＳデータ処理部９で実行されるＭＲＳデータの処理，解析について図２のフロー図に沿って説明する。
【００１７】
データ収集部７が収集したＭＲＳデータは時間領域のデータなので，まずこれをフーリエ変換して周波数領域のデータとする（２１）。この周波数領域のデータを，周波数に対してプロットしたものが図３に示されている。ここでＭＲＳのデータは点線で，このデータをｌｏｗｅｓｓフィルタでフィッティングしたデータが実線で示されている。
水からの信号と，巨大分子からの信号は，広い周波数領域に渡って分布し，それ以外の物質からの信号は，限られた周波数帯域に分布している。このため，図３に示したグラフの広域的な構造は，水信号と巨大分子からの信号の寄与によって決定され，これらの信号の寄与は，ベースラインと呼ばれる。ベースラインは，特定の物質からの信号のピークが存在する周波数帯域にも存在するため，これらの局所的な構造を隠してしまう。局所的な構造を解析するためには，ベースラインを特定し，この成分を原データから除去し，局所的に存在する信号のみのデータを作成する必要が有る。この処理は，ベースライン補正と呼ばれる。ベースラインの決定には，主に，ｌｏｗｅｓｓフィルタによる広域的なデータフィッティングか，ｗａｖｅｌｅｔ変換を原データに適用し，その後，大きいスケールに関係する成分のみを用いて，逆ｗａｖｅｌｅｔ変換を行なうことで，局所的な構造を除去するｗａｖｅｌｅｔフィルタが用いられる（２２）。これらの処理を用いて作成されたデータを，ベースラインデータと呼ぶことにする。原データから，ベースラインデータの値を引き算して，局所的構造のみのデータが得られる（２３）。このデータを局所データと呼ぶことにする。
【００１８】
先にも言及したように局所データの全領域に対し、検出可能な物質に対応したピークのフィッティングを直接行うと、フィッティングパラメータの数が膨大に成ってしまうので、フィッティングを行う領域を限定する。この目的のために、本発明では、判別分析法を用いて、閾値を決定する（２４）。具体的には次の通り。局所データを信号値の大きさに従って2分割する。基準値を決め、信号値が、基準値より大きい集団と、信号値が、基準値より小さい集団に分割する。分割する際、2つの集団それぞれの内部では、分散がなるべく小さくなるように、そして、2つの集団の間の分散はなるべく大きくなるように、基準値を決めることを考える。ここで、基準値より信号値が大きい集団の、データ数をｎ1、平均値をｍ1、分散をσ1と置き、基準値より信号値が小さい集団の、データ数をｎ2、平均値をｍ2、分散をσ2と置く。また、分割する前の母集団の、平均値をｍＴ、分散をσＴと置く。これらの値を用いて、集団内部の分散を表す指標、クラス内分散σＷと、2つの集団の間の分散を表す指標、クラス間分散σＢを次の様に定義する。
2008048878.tif 000003
2008048878.tif 000004
この指標を用いて、クラス内分散を小さくし、クラス間分散を大きくする基準値を探す。これは、クラス間分散をクラス内分散で割った値を最大とする基準値である。即ち、基準値を変えながら、次の量を計算し、この量を最大とする基準値が、求める基準値である。
2008048878.tif 000005
この判別分析法を、データの存在する点に直接適用しても良いが、周波数領域を、予め決められた数に分割し、その分割した領域内のデータについて平均値を計算し、平均値からなる集団を作成し、この集団について、判別分析法を適用して、基準値を決定しても良い。この基準値を閾値に設定し、信号値が、閾値より大きい領域と、信号値が閾値より小さい領域とに、周波数の領域を分ける。信号値が閾値以上の領域を信号領域、信号値が閾値以下の領域をノイズ領域とする。この処理を行うことで、周波数の領域は、図４のように、信号領域と、ノイズ領域に分割される（２５）。信号領域に番号を振っておく。
【００１９】
フィッティングを行う信号領域のデータに対して、フィッティングに先立ってスムーシングを行うと、最小値が、速く見つかる。このための処理として、lowessフィルタを用いても良いし、waveletフィルタを用いても良い。
また、以下に説明するスムージングの方法を用いても良い。即ち、局所データを周波数に対してプロットし、このグラフを画像として考え、図５に示すように、データを表すグラフより上の領域の画素値を0に設定し、グラフより下の領域の画素値に0でない、一定値を設定する。この関係は、逆でもよい。また、二つの領域の画素値が異なっていれば良いので、必ずしも0を設定する必要はない。この処理によって、得られた画像に対して、図６に示す拡大処理を行う。画素値が0でない画素に隣接する点の画素値を0でない値に設定する。次に、拡大した画像に対して、図７に示す縮小処理を行う。つまり、画素値が0の点に隣接する4点で、画素値が0でない点の画素値を0にする。拡大する画像は、1点ずつでも良いし、2点以上の点ずつでもよい。縮小する場合も同様である。ただし、拡大する点の数と、縮小する点の数は必ず同じでなければならない。この処理を行うことによって、この画像の輪郭は滑らかにスムーシングされる。即ち，図７の右側の図で，格子で示した部分だけ元に戻っていない。この結果，この部分だけ，構造がスムーシングされている。この画像に対し、図８に示すように、横軸に沿って、画素値が0でない点から、画素値が0に変わる点を探し、この点の縦の位置を、その横座標の信号値に変換する。これにより、グラフのスムーシングが実現できる。
【００２０】
フィッティングに用いるモデル関数は、検出可能な物質の共鳴周波数を、ピークに持つ、ガウス関数、または、ローレンツ関数の和で構成される。検出可能な物質の共鳴周波数を、Ｃｉとし、また、その物質に対応する振幅と、半値幅を、それぞれ、ＡｉとＢｉと表す。ここで、添え字iは、それぞれの物質に対応した番号である。モデル関数は、ガウス関数を用いた場合、(1)式で、ローレンツ関数を用いた場合には、(2)式で表される。
2008048878.tif 000006
2008048878.tif 000007
検出可能な物質の共鳴周波数は、予めテーブルに纏められて記憶部１３に格納されているものとする。また、フィッティングの対象としたい物質を検者が入力部１２から予め入力しておきその物質の共鳴周波数のみを、モデル関数に用いてフィッテイングを行っても良い。フィッティングは、信号領域でのみ行う。また、フィッティングは順番に、1つの領域ごとに行う。まず対象となる、信号領域内に、検出可能な物質の共鳴周波数が存在するかを調べる（２７）。共鳴周波数が存在しない場合には、処理を中止し、次の領域のフィッティングを行う。共鳴周波数が存在する場合には、存在する共鳴周波数を用いて、モデル関数を作成する。例えば、3つの共鳴周波数が存在し、ガウス関数を用いたフィッティングの場合、共鳴周波数を、Ｃ１、Ｃ２、Ｃ３とすると、（３）式に示すモデル関数が作成される。
2008048878.tif 000008
この場合、フィッティングパラメータは、振幅Ａ１、Ａ２、Ａ３そして、半値幅Ｂ１、Ｂ２、Ｂ３である。これらのパラメータを、対象となっている、信号領域において、局所データとの差の自乗の和が最小となる条件で決定する（２８）。最小値を求めるアルゴリズムは、Levenberg-Marquardt法(Ｗｉｌｌｉａｍ Ｈ. Press 他、Numerical Recipes in C  pp681-688)を用いて行う。
【００２１】
このようにして決定された，それぞれの該当物質の振幅値Ａ１、Ａ２、Ａ３および半値幅Ｂ１、Ｂ２、Ｂ３からなるパラメータに基づいて，それぞれ該当する物質の存在比を求めその結果を表示部１０に表示する（２９）。
【００２２】
以上説明した実施例では、データのフィッティングにおいて、各点での局所データと、モデル関数の値の差の自乗が最小になるように、各パラメータを決定したが、最小自乗法を使ったフィッティングに替えて，相互情報量を、局所データのグラフと、モデル関数のグラフで計算し、この相互情報量が最大となるようにパラメータフィッティングを行っても良い。具体的には、対象となる信号領域内に均等に分布する点を複数個(ＮＡ)選ぶ、この点をｘｉとする。iは、0から、ある決められた数ＮＡまでの値を取る。これらの点に対し、信号値をサンプリングする。得られたデータを、局所データに対しては、ｕｉ、モデル関数に関しては、ｖｉとする。また、それらの標準偏差をそれぞれσｕとσｖとする。サンプル点の数は、局所データと、モデル関数で異なってもよい。モデル関数の計測点を、ＮＢとする。これらの信号値
ｚ についてPazen window を用いて、確率を近似する。
局所データでの確率
2008048878.tif 000009
モデル関数での確率
2008048878.tif 000010
局所データとモデル関数の重なり確率
2008048878.tif 000011
この確率を用いて、相互情報量Iが次の様に計算できる。
2008048878.tif 000012
この相互情報量については、例えば、P.Viola
and W.M.Wells III “Alignment by Maximization of Mutual
Information” IJCV 24 (2) 137-154, 1997 に記載されている。また、この文献に記載されているように解析的に相互情報量の微分が計算できるので、最急勾配法を用いて最大値を検索し、その時のモデル関数のパラメータを求めてそれぞれの物質の存在比を求めることが出来る。
【図面の簡単な説明】
【００２３】
【図１】本発明ＭＲＩ装置の一実施例のブロック図。
【図２】図１のＭＲＳデータ処理部で実行される機能を説明するフロー図。
【図３】図１のＭＲＳデータ処理部で処理される周波数領域でのＭＲＳのデータの一例を示す図。
【図４】図１のＭＲＳデータ処理部で判別分析法を適用して周波数領域でのＭＲＳデータを分割する様子を説明する図。
【図５】図１のＭＲＳデータ処理部で、フィッティング処理に先立って分割された信号領域のＭＲＳデータスムージングのため画像化する手順を説明する図。
【図６】図５で画像化された画像を拡大する手順を説明する図。
【図７】図６で拡大された画像を縮小する手順を示す図。
【図８】図７で縮小された画像データからスムージングされた周波数領域のデータへの変換の手順を説明する図。
【符号の説明】
【００２４】
１  …  静磁場コイル
２  …  傾斜磁場コイル
３  …  ＲＦコイル
４  …  傾斜磁場コイル電源
５  …  送信部
６  …  受信部
７  …  データ収集部
８  …  シーケンサ
９  …  ＭＲＳデータ処理部
10  …  表示部
11  …  ＣＰＵ
12  …  入力部
13  …  記憶部
【図１】
2008048878.tif 000013
【図２】
2008048878.tif 000014
【図３】
2008048878.tif 000015
【図４】
2008048878.tif 000016
【図５】
2008048878.tif 000017
【図６】
2008048878.tif 000018
【図７】
2008048878.tif 000019
【図８】
2008048878.tif 000020