JP4019146B2

JP4019146B2 - 拡散テンソル場の構造強調濃淡画像表示方法及び装置

Info

Publication number: JP4019146B2
Application number: JP2003294707A
Authority: JP
Inventors: 茂村木
Original assignee: National Institute of Advanced Industrial Science and Technology AIST
Current assignee: National Institute of Advanced Industrial Science and Technology AIST
Priority date: 2003-08-18
Filing date: 2003-08-18
Publication date: 2007-12-12
Anticipated expiration: 2023-08-18
Also published as: JP2005058570A

Description

本発明は、拡散テンソル画像に基づいてその構造を強調した構造強調濃淡画像を形成して、それを画像表示手段に表示する拡散テンソル場の構造強調濃淡画像表示方法及び装置に関するものである。

直接的には視ることができない体内や脳内部を可視化する画像診断は、医療現場において大きな役割を担っている。そのなかでもＭＲＩ（ＭａｇｎｅｔｉｃＲｅｓｏｎａｎｃｅＩｍａｇｉｎｇ：核磁気共鳴撮影法）は、骨に遮られることなく脊髄や軟骨までも明確に可視化することができ、さまざまな病気の早期発見に貢献している。またＭＲＩは撮影条件により強調される組織が異なり、多様な強調画像を作ることができる。これにより、特定の組織の特性量などを計測することができるようになった。

ＭＲＩにより計測されたテンソルデータを効果的に可視化することは、生理学や神経分析などさまざまな研究分野の発展促進につながると考えられる。しかし、これまでにスカラ場やベクトル場の効果的な可視化法は数多く提案されているものの、テンソル場のもつ多様な情報を、目的に応じて柔軟に可視化できる手法はほとんど開発されていない。

例えば、非特許文献１に示されるように２次元流れ場の効果的な可視化法として知られるＬＩＣ（ＬｉｎｅＩｎｔｅｇｒａｌＣｏｎｖｏｌｕｔｉｏｎ）のボリューム版を用いて、拡散テンソル場を可視化する方法について、非特許文献２に開示されている。また非特許文献３には、ノンフォトリアリスティックレンダリングで用いられる非一様拡散シミュレーションを応用した、画像特徴の強調法、および、平滑化法に関する技術が開示されている。
Ｂ．ＣａｂｒａｌＬ．Ｌｅｅｄｏｍの論文"Ｉｍａｇｉｎｇｖｅｃｔｏｒｆｉｅｌｄｓｕｓｉｎｇｌｉｎｅｃｏｎｖｏｌｕｔｉｏｎ，" ＩｎＣｏｍｐｕｔｅｒＧｒａｐｈｉｃｓＰｒｏｃ．，ＡｎｎｕａｌＣｏｎｆｅｒｅｎｃｅＳｅｒｉｅｓ，ＡＣＭＳＩＧＧＲＡＰＨ，Ａｕｇ．１９９３，第２６３頁−第２７０頁．安藤祥子，藤代一成及び村木茂の論文「ＶＬＩＣによる拡散テンソル場可視化法の検討」，情報処理学会第６５回全国大会，２００３年３月，第１１９頁−第１２０頁ＪｏａｃｈｉｍＷｅｉｃｋｅｒｔの論文"Ｃｏｈｅｒｅｎｃｅ−ｅｎｈａｎｃｉｎｇｄｉｆｆｕｓｉｏｎｏｆｃｏｌｏｒｉｍａｇｅｓ，" ＩｍａｇｅａｎｄＶｉｓｉｏｎＣｏｍｐｕｔｉｎｇ，ｖｏｌ．１７，１９９９，第２０１頁−第２１２頁

ここでテンソルとは、ベクトルや行列を一般化した数学的な概念で、０階のテンソルはスカラ、１階のテンソルはベクトル、２階のテンソルは行列である。各点にテンソルが与えられた空間をテンソル場と呼び、例えば地中の応力分布は２階３次元テンソル場と考えられる。一般相対性理論によれは、我々が住む世界は、各点に計量テンソル（無限小線素の長さの単位）が与えられた、２階４次元テンソル場になる。拡散テンソルは、ブラウン運動する微小粒子（分子など）の単位時間後の分布を示す物理量で下記［数１］に示すような３次元の対角行列で示される。

拡散テンソルＤの固有値λ_ｉ（ｉ＝０，１，２）及び固有ベクトルｖ_ｉ（ｉ＝０，１，２）は、下記［数２］の方程式を解くことにより求められる。

しかしながら、拡散テンソルＤによっては固有値が３個未満である場合もある。固有値λ_ｉの中で絶対値が最大のものを最大固有値、それに対応する固有ベクトルｖ_ｉを最大固有ベクトルとすると、最大固有ベクトルは最も拡散の大きい方向を示し、最大固有値はその拡散係数を示すことになる。ただし、上記［数２］の式が、両辺のｖ_ｉの符号を逆転しても成立することから、固有ベクトルの方向の正逆は不定であることに注意が必要である。前述のように、スカラ場やベクトル場の可視化法は数多く提案されているが、各点の自由度の多いテンソル画像の場合、そのすべてを２次元のディスプレイ上に可視化するのは困難である。

核磁気共鳴イメージング（ＭＲＩ）の一つである拡散強調像（diffusion-weighted imaging:ＤＷＩ）は、プロトンの拡散を強調する傾斜磁場（motion probing gradient:ＭＰＧ）を加えて撮影した画像である。ＭＰＧの方向や大きさを変えて撮影した複数のＤＷＩから、拡散テンソル画像を作ることができる。つまり、ＭＰＧの印加方向をＧ_ｉ＝（Ｇ_ｉｘ，Ｇ_ｉｙ，Ｇ_ｉｚ）^Ｔとすれば、拡散強調画像の信号強度は、下記［数３］のように表される。

ここでＳ_０はＭＰＧを印加しない画像（Ｔ２強調像）の信号強度、上付きのＴは転置を意味するものとする。対角行列である拡散テンソルＤの未知数は６個なので、６方向以上のＭＰＧを用いて上記［数３］の式を連立させれば、拡散テンソルを求めることができる。

脳について計測された拡散テンソル画像では、細胞膜により水の拡散が制限されているために、最大固有ベクトルが神経繊維の走行方向とほぼ一致する。したがって、脳の拡散テンソル画像の可視化は、脳の神経走行の可視化と密接に関連し、いくつかの手法が提案されている。たとえば、図１に示すように、拡散テンソルの３組の固有値と固有ベクトルを、主軸の長さと方向に対応させた楕円体をテンソル画像の各点に置く楕円体表記法がある。この方法では、テンソルから３個の固有値とそれに対応する固有ベクトルを計算する。拡散テンソルにおける固有値とその対応固有ベクトルは、それぞれ拡散の大きさとその方向を示しているので、楕円体表記法では、固有ベクトルを楕円体の各主軸とし、固有値をその主軸方向の大きさとすることにより、単一のテンソル量を描くことになる。しかし、単一のテンソルを楕円体などのグリフで表現した可視化法をテンソル場に適用すると、隣接するテンソル間が不連続となり、データの全体像の把握が困難となってしまう。また、ディスプレイの制約により一度に表示できるテンソル数が制限されてしまうため、サイズの大きなボリュームデータなどを可視化することは難しい。また、隣接する楕円体が重なりあうため、密なテンソル場を可視化することが困難である。なお楕円体表記法については、株式会社秀潤社発行の「これでわかる拡散ＭＲＩ」に詳しく説明されている。図２には、楕円体表記法により表示した脳の断層画像とその要部の拡大図を横に並べて示してある。また図３も楕円体表記法による表示例を示している。その他の表記法としては、図４に示すような、最大固有値の固有ベクトル方向を追跡した曲線を描画するトラクトグラフィー（tractography）などがある。しかし、トラクトグラフィーでは、最大固有方向のみを追跡するため、テンソルの持つ情報の一部しか利用できず、先に述べた固有ベクトルの正逆方向の不定性により、原理的に分岐や交差する神経の走行を正しく可視化することができない。また、隣接する線素が重なりあうため、密なテンソル場を可視化することが困難である。

他の表記としては、固有値を用いて、不等方度の一つであるＦＡ（Fractional Anisotropy)を各点について下記の［数４］で算出し、濃淡値や色で表示する方法などがある。図５にはこの方法で画像を表記した一例を示してある。

ここで〈Ｄ〉は［数５］に示すみかけの拡散係数（Apparent Diffusion Coefficient:ＡＤＣ）である。

この方法では、図６のＴ２強調像にくらべて、異方性の高い神経線維が強調されていることがわかるが、方向に関する情報が失われる問題がある。

本発明の目的は、拡散テンソル画像が持つ、密で連続なテンソル場の構造を強調して表示することができる拡散テンソル場の構造強調濃淡画像表示方法および装置を提供することにある。

本発明のより具体的な目的は、３次元の各点に拡散テンソル（２階３次元テンソル）が与えられた一次の拡散テンソル画像に対して、固有値分解を利用して生成した二次の拡散テンソル画像を用いて、二次の拡散テンソル画像と同じ解像度の濃淡画像を拡散させ、それにより発生するパターンを用いて、一次の拡散テンソル画像が持つ、密で連続なテンソル場の構造を強調して表示することができる拡散テンソル場の構造強調濃淡画像表示方法および装置を提供することにある。

本発明は、画像を表示する各点に拡散テンソルが与えられて構成される一次の拡散テンソル画像に基づいてその構造を強調した構造強調濃淡画像を形成してそれを画像表示手段に表示する拡散テンソル場の構造強調濃淡画像表示方法を対象とする。ここで一次の拡散テンソル画像とは、ＭＲＩ（ＭａｇｎｅｔｉｃＲｅｓｏｎａｎｃｅＩｍａｇｉｎｇ：核磁気共鳴撮影法）で計測された複数の拡散強調像（diffusion-weighted imaging:ＤＷＩ）等から形成される拡散テンソルデータを各点に与えた画像データを意味する。

本発明の方法では、第１のステップとして、まず一次の拡散テンソル画像を構成する各点の一次の拡散テンソルについて固有値分解を行い、固有値分解により得た各点の二次の拡散テンソルに基づいて二次の拡散テンソル画像を生成して保存する。ここで二次の拡散テンソル画像とは、各点の二次の拡散テンソルデータによって形成される画像データを意味する。具体的に第１のステップでは、一次の拡散テンソル画像を構成する各点の一次の拡散テンソルについてそれぞれ固有値分解を行い、各点についての固有値と固有ベクトルの画像データを保存する。次に第１の制御値と固有値とから各点の固有値の重みを決定する。この決定した重みと、固有値と固有ベクトルとから各点の二次の拡散テンソルを計算する。具体的には、各点の固有値の重みは、二次の拡散テンソルの一部を構成する３つの固有値の関数と定義することができる。次に各点の二次の拡散テンソルに基づく二次の拡散テンソル画像を保存する。なお前述の第１の制御値とは、どのような二次の拡散テンソル画像を生成するかをユーザーが選択するために定められる値であり、この値を変えることにより、二次の拡散テンソルの異方性や、二次の拡散テンソル画像の大きさ（解像度）をユーザーの好みに合わせて調整することができる。なお各点の固有値の重みを、その絶対値が最大の固有値に対して１となり、その絶対値が最大の固有値以外の値に対しては０になるようにすると、最も拡散係数の大きい方向にのみ拡散が生じ、神経線維などの線状の構造を強調することができる。

次に第２のステップでは、二次の拡散テンソル画像と同じ解像度の濃淡画像を生成する。具体的には、二次の拡散テンソル画像の各点について、第２の制御値と、固有値とから濃淡値を決定する。そして各点について決定した濃淡値を有する二次の拡散テンソル画像と同じ解像度になる濃淡画像を生成する。濃淡画像とは、拡散される粒子の密度を値とする画像データである。またここで使用する第２の制御値とは、各点にどのような濃淡値を生成するかをユーザーが選択するための値であり、この値を変えることにより、固有値の小さい点の濃淡値をゼロにして拡散によるパターンが生じないようにしたり、パターンの密度を変えたりできる。具体的に、第２のステップでは、二次の拡散テンソル画像の各点について、濃淡値を調整するための制御値（第２の制御値）と固有値分解により得た固有値とから濃淡値を決定して、二次の拡散テンソル画像と同じ解像度になる濃淡画像を生成する。なお濃淡画像の点各の濃淡値は、構造強調を行う濃淡画像の点においては乱数であり、それ以外の濃淡画像の点では一定値であることが好ましい。このようにすると、構造強調を行う領域のみにコントラストの強いパターンが生じ、その構造が強調される。なお二次の拡散テンソル画像の解像度は、一次の拡散テンソル画像の解像度と同じでもよいが、それよりも大きくしてもよい。二次の拡散テンソル画像の解像度を一次の拡散テンソル画像の解像度よりも大きくすると、より細かいパターンが発生し、一次の拡散テンソル画像に含まれる細部の構造が鮮明に表示できる。

第３のステップでは、濃淡画像の各点において、濃淡値の時間微分を濃淡画像と各点の二次の拡散テンソルを用いて算出する。ここで濃淡値の時間微分とは、拡散方程式によって計算される値を意味する。そして第４のステップでは、この時間微分を設定された時間幅で積分して得た積分値を濃淡画像の対応する各点の濃淡値に加算して更新された濃淡画像を生成する。この時間幅は、流入量が大きくなり過ぎず、かつ、拡散シミュレーションが十分高速に行えるように、使用するデータに応じて適宜に定める。

そして第５のステップでは、更新された濃淡画像に対して、第３及び第４のステップを繰り返して実行して構造強調濃淡画像を作成する。第３及び第４のステップを繰り返す反復回数は、任意に定めることになる。反復回数が多いほど、各点により広い範囲の拡散テンソル情報が反映されることになるが、ある程度の反復回数を経ると濃淡画像がぼけて情報が失われるため、反復回数は適宜に定めることになる。

また画像を表示する各点に拡散テンソルが与えられて構成される一次の拡散テンソル画像に基づいてその構造を強調した構造強調濃淡画像を形成してそれを画像表示手段に表示することをコンピュータを利用して実現するために用いる本発明の構造強調濃淡画像表示用プログラムは、一次の拡散テンソル画像を構成する各点の一次の拡散テンソルについて固有値分解を行って得た各点の二次の拡散テンソルに基づいて二次の拡散テンソル画像を生成して保存する第１の機能と、二次の拡散テンソル画像と同じ解像度の濃淡画像を生成する第２の機能と、濃淡画像の各点において、濃淡値の時間微分を濃淡画像と各点の二次の拡散テンソルを用いて算出する第３の機能と、時間微分を設定された時間幅で積分して得た積分値を濃淡画像の対応する各点の濃淡値に加算して更新された濃淡画像を生成する第４の機能と、更新された濃淡画像に対して、第３及び第４の機能を繰り返して実行することにより構造強調濃淡画像を作成する第５の機能とをコンピュータに実現させる構成を備えている。

本発明の方法を実施する構造強調濃淡画像表示装置は、一次の拡散テンソル画像を構成する各点の一次の拡散テンソルについて固有値分解を行う固有値分解手段と、固有値分解により得た各点についての固有値と固有ベクトルの画像を保存する固有値分解保存手段と、制御値と固有値とから各点の固有値の重みを決定する重み決定手段と、重みと、固有値と固有ベクトルとから二次の拡散テンソルを計算する二次の拡散テンソル計算手段と、二次の拡散テンソルに基づく二次の拡散テンソル画像を保存する二次の拡散テンソル画像保存手段と、二次の拡散テンソル画像と同じ解像度の濃淡画像を生成する濃淡画像生成手段と、濃淡画像の各点において、濃淡値の時間微分を濃淡画像と各点の二次の拡散テンソルを用いて算出する時間微分量算出手段と、時間微分を設定された時間幅で積分して得た積分値を濃淡画像に加算して更新された濃淡画像を生成する更新濃淡画像生成手段と、更新された濃淡画像の各点について、繰り返し時間微分算出手段及び更新濃淡画像生成手段による処理を行わせて構造強調濃淡画像を作成する構造強調濃淡画像作成手段とから構成することができる。

本発明によれば、一次の拡散テンソル画像では、はっきりと見ることができなかった一次の拡散テンソル画像に含まれる構造の特徴を強調して表示することができるようになる利点が得られる。その結果、本発明によれば、従来の表示法によれば不明瞭にならざるを得なかった、例えば、脳の神経路部分も明確に表示できるようになる。

図７は、本発明の方法を実施するために用いる構造強調濃淡画像表示装置の実施の形態の一例の構成を示すブロックであり、図８は図７の実施の形態を、コンピュータを用いて実現する場合に用いるソフトウエアのアルゴリズムを示すフローチャートである。

図７において、符号１で示したブロックは、ＭＲＩのように拡散テンソルデータ（一次の拡散テンソル画像）を出力する一次の拡散テンソル画像出力装置である。本実施の形態の構造強調濃淡画像表示装置３では、この装置１から出力される一次の拡散テンソルデータまたは拡散テンソル画像を入力として、この一次の拡散テンソル画像に含まれる構造を強調した濃淡画像として画像表示手段２に表示する。

この構造強調濃淡画像表示装置３は、一次の拡散テンソル画像を構成する各点の一次の拡散テンソルについて固有値分解を行う固有値分解手段４を有している。固有値分解手段４による固有値分解により得た一次の拡散テンソル画像の各点についての固有値と固有ベクトルの画像（画像データ）は固有値分解保存手段５に保存される。固有値分解保存手段５に保存された固有値は、重み決定手段６に入力される。重み決定手段６では、固有値と制御値（第１の制御値）とから各点の固有値の重みを決定する。重み決定手段６で用いる制御値とは、どのような二次の拡散テンソル画像を生成するかをユーザーが選択するための値であり、この値を変えることにより、二次の拡散テンソルの異方性や、二次の拡散テンソル画像の大きさ（解像度）をユーザーの好みに合わせて調整することができる。

二次の拡散テンソル計算手段７は、重み決定手段６により決定した重みと、固有値分解保存手段５に保存された固有値及び固有ベクトルとから、画像の各点についての二次の拡散テンソルを計算する。計算された二次の拡散テンソルから構成される二次の拡散テンソル画像は、二次の拡散テンソル画像保存手段８に保存される。

次に、濃淡画像生成手段９は、二次の拡散テンソル画像保存手段８に保存された二次の拡散テンソル画像と同じ解像度の濃淡画像を生成する。具体的には、濃淡画像生成手段９は、二次の拡散テンソル画像の各点について、制御値（第２の制御値）と、固有値とから二次の拡散テンソル画像と同じ解像度になる濃淡値を決定する。ここで使用する制御値とは、各点にどのような濃淡値を生成するかをユーザーが選択するための値であり、この値を変えることにより、固有値の小さい点の濃淡値をゼロにして拡散によるパターンが生じないようにしたり、パターンの密度を変えたりできる。そして各点について決定した濃淡値を有する濃淡値データから濃淡画像を生成する。

次に、時間微分算出手段１０は、濃淡画像生成手段９により生成した濃淡画像の各点において、濃淡値の時間微分を、濃淡画像と各点の二次の拡散テンソルとを用いて算出する。濃淡画像と各点の二次の拡散テンソルとを用いて濃淡値の時間微分（濃淡値の流入量）を算出する際に用いる式は後述する［数６］に示す式である。

更新濃淡画像生成手段１１は、時間微分算出手段１０により算出した時間微分を設定された時間幅Δｔで積分して得た積分値を濃淡画像に加算して更新された濃淡画像を生成する。この時間幅Δｔは、濃淡値の時間微分または濃淡値の流入量が大きくなり過ぎず、かつ、拡散シミュレーションが十分高速に行えるように、使用するデータに応じて適宜定めればよい。このように近傍の点から各点に拡散してくる濃淡値の流入量を考慮すれば、濃淡画像をより構造を強調したものにすることができる。この例では、濃淡画像生成手段９により生成された濃淡画像及び更新濃淡画像生成手段１１によって生成された更新された濃淡画像は、いずれも画像表示手段２に表示される。したがって処理過程を見ながら最適な反復回数を定めることができる。

最後に、構造強調濃淡画像作成手段１２は、更新濃淡画像生成手段１１により更新された濃淡画像の各点について、繰り返し時間微分算出手段１０及び更新濃淡画像生成手段１１による処理を行わせて構造強調濃淡画像を作成する。この処理の反復回数は、構造強調濃淡画像作成手段１２によって適宜に定められる。反復回数が多いほど、各点により広い範囲の拡散テンソル情報が反映されることになるが、ある程度の反復回数を経ると濃淡画像がぼけて情報が失われる。例えば、256×256画素程度ＭＲＩの脳の神経路部分を強調して表示する場合であれば、10回〜20回程度の反復回数を選択すればよい。

図７に示した実施の形態の装置をコンピュータを用いて実現する場合に用いるソフトウエアのアルゴリズムを示す図８のフローチャートを参照すると、固有値・固有ベクトル画像生成処理ステップＳＴ１によって、固有値分解手段４と固有値分解保存手段５とが実現されている。また拡散テンソル画像生成処理ステップＳＴ２により、重み決定手段６と、二次の拡散テンソル計算手段７と二次の拡散テンソル画像保存手段８とが実現されている。更に濃淡画像生成処理ステップＳＴ３により、濃淡画像生成手段９が実現されている。また主として拡散シミュレーション処理ステップＳＴ６により時間微分算出手段１０と、更新濃淡画像生成手段１１とが実現され、反復回数入力ステップＳＴ４と、濃淡画像表ステップＳＴ５及びＳＴ８と、反復終了確認ステップＳＴ９とにより、構造強調濃淡画像作成手段１２が実現されている。なお図８のステップＳＴ９に示されるように、重みや反復回数を変えて一連の処理を再施行するために、ステップＳＴ２からステップＳＴ８までの処理は繰り返される。

なお図９は、固有値・固有ベクトル画像生成処理ステップＳＴ１の詳細を示すフローチャートであり、図１０は拡散テンソル画像生成処理ステップＳＴ２の詳細を示すフローチャートであり、図１１は濃淡画像生成処理ステップＳＴ３の詳細を示すローチャートであり、図１２は拡散シミュレーション処理ステップＳＴ６の詳細を示すフローチャートである。なおこれらのフローチャートについては、より具体的に本発明の方法を実施する場合について説明した後に説明する。

本発明では、表示の対象が拡散テンソルであることを利用して、拡散テンソル情報に基づいた粒子の拡散シミュレーション（更新された濃淡画像の生成）を行うことにより（図８のステップＳＴ６により）、密で連続な拡散テンソル場の構造を強調して表示する。

拡散シミュレーションは、下記の［数６］に示された拡散方程式に従って行う。

上記式において、Ｕ_Ｔ（ｘ，ｙ，ｚ）は、時刻Ｔでの座標（ｘ，ｙ，ｚ）の粒子密度で、本実施の形態においては、シミュレーションに使用する拡散テンソル画像と同じ解像度に離散化された濃淡画像（濃淡画像生成手段９の出力）である。

初期の濃淡画像Ｕ_０（ｘ，ｙ，ｚ）として、例えば神経領域をランダムドットとした画像（図１３に示すＴ＝０の画像）を用いる場合を例にして説明する。なお図１３のような濃淡画像は、図５に示したようなＦＡ（Fractional Anisotropy)の低い画素や、図６に示すようなＴ２強調像の輝度の低い画素を背景画素とした上で、それ以外の画素に乱数を発生させることで、容易に作ることができる。この場合には、時間幅Δｔを使って、下記の［数７］のような積分を行うことにより、濃淡画像を拡散させて得られる更新された濃淡画像を得ることができる。

上記式において、更新を適当な回数反復することで、濃淡画像を拡散テンソル情報に基づいて拡散させたパターン即ち更新された濃淡画像を生じさせることができ、拡散テンソル場の構造を容易に把握できるようになる。

ここで、［数６］に示した拡散の方程式のＤ’として、核磁気共鳴イメージング（ＭＲＩ）の一つである拡散強調像ＤＷＩを複数用いて算出した拡散テンソルＤをそのまま一次の拡散テンソルとして用いると不都合を生じる。これは、計算により求められた拡散テンソルＤでは、負の固有値を含み拡散テンソルにならない場合があり、また、最大固有ベクトル方向以外の拡散成分が、濃淡画像を過度に拡散するからである。そこで本実施の形態では、各点の拡散テンソルの固有値分解を行い、以下のような方法で生成した第二の拡散テンソルを使用する。

固有値分解とは、一次の拡散テンソルをその固有値、固有ベクトルを用いて、下記の［数８］に示すように分解する方法である。

上記式において、ＡはＡ＝（ｖ_０ｖ_１ｖ_２）のような行列である。本実施の形態では、二次の拡散テンソルとして、重みｗ＝（ｗ_０ｗ_１ｗ_２）を用いて修正した下記の［数９］に示すものを用いる。

ここで重みｗは、制御値と固有値とから決定する。具体的には、各重みはｗ_i＝ｗ_ｘ（λ_０，λ_１，λ_２）のような固有値の関数として求める。ここで、i=0,1,2、（λ_０，λ_１，λ_２）が固有値である。また、どのような関数ｗ_ｘを選択するか、あるいは、二次の拡散テンソル画像をどのような解像度で生成するかは、制御値（第１の制御値）として入力する。

ここで重みが固有値に依存しない簡単な例としてｗ＝（１，１，１）の場合は、１次の拡散テンソル画像情報をそのまま使った拡散がシミュレーションされるが、最大固有ベクトル方向以外の拡散成分のために、図１４に示すようにぼけた画像が生成される。図１４の画像は、Ｔ＝１０（反復回数１０回）、ｗ＝（１，１，１）の画像である。ｗ＝（１，０，０）の場合は、最大固有値方向のみに拡散が生じるため、図１５に示す画像のように神経走行が強調される。図１５の画像は、Ｔ＝１０（反復回数１０回）、ｗ＝（１，０，０）の画像である。なお濃淡画像の表示方法としては、図１３〜図１６に示した断面表示以外に、ボリュームレンダリングによる３次元表示や、Ｔ２強調像の上に濃淡画像を重ねて表示する方法などが考えられる。

次に図９〜図１２に示したフローチャートを説明しながら、更に次のステップの説明を進める。図９の固有値・固有ベクトル画像生成処理ステップＳＴ１では、一次の拡散テンソル画像の生成処理ステップＳＴ１１で、拡散テンソル画像出力装置に依存する出力形態から実際に拡散テンソル画像を生成する。そしてステップＳＴ１２で、固有値・固有ベクトル画像を初期化する。ここで固有値・固有ベクトル画像を初期化するとは、画像メモリを確保し初期値を代入することである。その後、ステップＳＴ１３により、一次の拡散テンソル画像の処理すべき点（処理点）を選択する。なお画像の全ての点について処理が行われる。そしてステップＳＴ１４で、その点における一次の拡散テンソルを固有値分解する。その後ステップＳＴ１５で、固有値分解により得た固有値と固有値ベクトルとを保存する。

図１０の拡散テンソル画像生成処理ステップＳＴ２は、前述の［数９］に示した式に基づいて、第二の拡散テンソル画像を生成する処理である。まずステップＳＴ２１においては、前述の第１の制御値が入力される。ステップＳＴ２２では、二次の拡散テンソル画像の初期化を行う。この初期化とは、二次の拡散テンソル画像に必要なメモリを確保し初期値を代入することである。そしてステップＳＴ２３で処理点が選択される。次にステップＳＴ２４で、制御値と固有値とから各点の固有値の重みを決定する。次にステップＳＴ２５で、この決定した重みと、固有値と固有ベクトルとから各点の二次の拡散テンソルを計算する。そしてステップＳＴ２６で、各点の二次の拡散テンソルに基づく二次の拡散テンソル画像を保存する。ステップＳＴ２７により、画像の全ての点についてこの処理が繰り返される。ステップＳＴ２１の制御値入力では、どのような重みを選択するかといった情報が第１の制御値として与えられる。例えば、図１５では重みｗ＝（１，０，０）が得られる関数ｗ_ｘが制御値入力にあたる。この重みであれば、最大固有値方向のみの拡散を行わせる。しかし神経線維が分岐、交差する場所では、最大固有値以外も有意な大きさを持ち、図１の楕円体が扁平になることが知られており、ｗ＝（１，０，０）では不連続を生じる可能性がある。そのため、偏平度Ｃ_ｐを下記の［数１０］の式で算出する。

そして重みをｗ＝（１，Ｃ_ｐ，０）などとすることにより、神経線維が分岐、交差する場所での不連続を緩和することも可能である。図１６の画像は、Ｔ＝１０（反復回数１０回）、ｗ＝（１，Ｃ_ｐ，０）の画像である。なお図１０のフローチャートにおけるステップＳＴ２１において与える制御値入力として、二次の拡散テンソル画像の解像度を与え、補間処理により一次の拡散テンソル画像よりも高い解像度の、二次の拡散テンソル画像を生成することで、より解像度の高い濃淡画像を生成することも可能である。また、ＦＡなどの不当方度が適当な閾値より小さい場合にｗ＝（０，０，０）とすることで、そのような領域で拡散が生じないようにすることも可能である。

図１１の濃淡画像生成処理ステップＳＴ３は、図１３に示すような初期濃淡画像を生成するための処理である。ステップＳＴ３１で制御値（第２の制御値）を入力する。ここで使用する制御値とは、各点にどのような濃淡値を生成するかをユーザーが選択するための値であり、この値を変えることにより、固有値の小さい点の濃淡値をゼロにして拡散によるパターンが生じないようにしたり、パターンの密度を変えたりできる。そしてステップＳＴ３２では、二次の拡散テンソル画像の初期化が行われる。ここで行う初期化では、二次の拡散テンソルに必要なメモリの確保と初期値の代入を行う。次にステップＳＴ３３で、処理点の選択が行われ、ステップＳＴ３４で、二次の拡散テンソル画像の各点について、制御値と、固有値とから濃淡値を決定する。そしてステップＳＴ３５で、各点について決定した濃淡値を濃淡画像として確保したメモリに代入する。またステップＳＴ３６で、画像の全ての点で濃淡画像の処理が行われる。図１３の画像では、不等方度の一つであるＦＡ（Fractional Anisotropy)が大きく、Ｔ２強調像（プロトンの拡散を強調する傾斜磁場ＭＰＧを印加しない画像）の輝度が高い部分をランダムドットとし、それ以外を背景として黒にしている。しかしステップＳＴ３１において、制御値入力を変えることにより、ＦＡやＴ２強調像の閾値や背景の色を変えることも可能である。また、濃淡画像中の特定部分に、特徴のある輝度や色を与えることで、その位置からの拡散を強調して表示することも可能である。

図１２の拡散シミュレーション処理６では、ステップＳＴ６１で時間幅Δｔを設定する。この時間幅Δｔを変更することにより、シミュレーションのスピードを変えることが可能になる。本実施の形態では、この時間幅Δｔを0.001〜1の範囲の値を選択している。次にステップＳＴ６２で、処理点を選択し、ステップＳＴ６３で濃淡画像の処理点において、濃淡値の時間微分値を濃淡画像とその点の二次の拡散テンソルを用いて算出する。具体的には、前述の［数６］の拡散の方程式を用いる。次に、ステップＳＴ６４でこの時間微分を設定された時間幅で積分した上でその積分値を濃淡画像の対応する点の濃淡値に加算して更新された濃淡画像を生成する。このステップでは、前述の［数７］の式を用いて更新処理を行う。これらのステップはステップＳＴ６５により画像の全ての点について実行される。そして図８に示したステップＳＴ４で定めた反復回数分だけ上記のシミュレーション処理が繰り返される。反復による更新処理が終了したら、図８に示したステップＳＴ８へと進んで構造を強調した構造強調濃淡画像の表示が実行される。

本発明は、直接的には視ることができない体内や脳内部等を、従来よりも構造を強調して可視化することができるので、産業上の利用可能性が非常に高い。

楕円体表記法を説明するために用いる図である。楕円体表記法により表示した脳の断層画像とその要部の拡大図を横に並べて示した図である。楕円体表記法による表示例を示す図である。最大固有値の固有ベクトル方向を追跡した曲線を描画するトラクトグラフィーの表示例を示す図である。ＦＡ（Fractional Anisotropy)を各点について算出し、濃淡値や色で表示する方法で脳の断層画像を示した図である。Ｔ２強調像により脳の断層画像を示した図である。本発明の方法を実施するための構造強調濃淡画像表示装置の実施の形態の一例の構成を示すブロックである。図７の実施の形態をコンピュータを用いて実現する場合に用いるソフトウエアのアルゴリズムを示すフローチャートである。固有値・固有ベクトル画像生成処理ステップＳＴ１の詳細を示すフローチャートである。拡散テンソル画像生成処理ステップＳＴ２の詳細を示すフローチャートである。濃淡画像生成処理ステップＳＴ３の詳細を示すローチャートである。拡散シミュレーション処理ステップＳＴ６の詳細を示すフローチャートである。脳の神経領域をランダムドットとした脳の断層画像の初期の濃淡画像（Ｔ＝０の画像）の例を示す図である。Ｔ＝１０（反復回数１０回）、ｗ＝（１，１，１）としたときの脳の断層画像の濃淡画像である。Ｔ＝１０（反復回数１０回）、ｗ＝（１，０，０）の脳の断層画像の濃淡画像である。Ｔ＝１０（反復回数１０回）、ｗ＝（１，Ｃ_ｐ，０）の脳の断層画像の濃淡画像である。

符号の説明

１一次の拡散テンソル画像出力装置
２画像表示手段
３構造強調濃淡画像表示装置
４固有値分解手段
５固有値分解保存手段
６重み決定手段
７二次の拡散テンソル計算手段
８二次の拡散テンソル画像保存手段
９濃淡画像生成手段
１０時間微分算出手段
１１更新濃淡画像生成手段
１２構造強調濃淡画像作成手段

Claims

画像を表示する各点に拡散テンソルが与えられて構成される一次の拡散テンソル画像に基づいてその構造を強調した構造強調濃淡画像を形成してそれを画像表示手段に表示する拡散テンソル場の構造強調濃淡画像表示方法であって、
前記一次の拡散テンソル画像を構成する前記各点の一次の拡散テンソルについて固有値分解を行って得た二次の拡散テンソルに基づいて二次の拡散テンソル画像を生成して保存する第１のステップと、
前記二次の拡散テンソル画像と同じ解像度の濃淡画像を生成する第２のステップと、
前記濃淡画像の各点において、濃淡値の時間微分を前記濃淡画像と前記各点の二次の拡散テンソルを用いて算出する第３のステップと、
前記時間微分を設定された時間幅で積分して得た積分値を前記濃淡画像の対応する前記各点の濃淡値に加算して更新された濃淡画像を生成する第４のステップと、
前記更新された濃淡画像に対して、前記第３及び第４のステップを繰り返して実行することにより前記構造強調濃淡画像を作成する第５のステップとからなる拡散テンソル場の構造強調濃淡画像表示方法。
前記第１のステップは、
前記一次の拡散テンソル画像を構成する前記各点の一次の拡散テンソルについて固有値分解を行って前記各点について固有値と固有ベクトルの画像を保存するステップと、
異方性や解像度を調整するための制御値と前記固有値とから前記各点の固有値の重みを決定して、前記重みと、前記固有値と前記固有ベクトルとから二次の拡散テンソルを計算するステップと、
前記二次の拡散テンソルに基づく前記二次の拡散テンソル画像を保存するステップとからなる請求項１に記載の拡散テンソル場の構造強調濃淡画像表示方法。
前記各点の固有値の重みが、その絶対値が最大の固有値に対して１となり、その絶対値が前記最大の固有値以外の値に対して０である請求項２に記載の拡散テンソル場の構造強調濃淡画像表示方法。
前記各点の固有値の重みが、前記二次の拡散テンソルの一部を構成する３つの前記固有値の関数である請求項２または３に記載の拡散テンソル場の構造強調濃淡画像表示方法。
前記第２のステップでは、二次の拡散テンソル画像の各点について、濃淡値を調整するための制御値と前記固有値分解により得た固有値とから濃淡値を決定して、前記二次の拡散テンソル画像と同じ解像度になる前記濃淡画像を生成することを特徴とする請求項１に記載の拡散テンソル場の構造強調濃淡画像表示方法。
前記濃淡画像の点各の前記濃淡値は、構造強調を行う前記濃淡画像の点においては乱数であり、それ以外の前記濃淡画像の点では一定値であることを特徴とする請求項５に記載の拡散テンソル場の構造強調濃淡画像表示方法。
前記二次の拡散テンソル画像の前記解像度が、前記一次の拡散テンソル画像の解像度よりも大きいことを特徴とする請求項１乃至６のいずれか１つに記載の拡散テンソル場の構造強調濃淡画像表示方法。
前記二次の拡散テンソル画像の前記解像度が、前記一次の拡散テンソル画像の解像度と同じであることを特徴とする請求項１乃至６のいずれか１つに記載の拡散テンソル場の構造強調濃淡画像表示方法。
前記一次の拡散テンソル画像が、ＭＲＩ計測装置から出力される計測データである請求項１に記載の拡散テンソル場の構造強調濃淡画像表示方法。
画像を表示する各点に拡散テンソルが与えられて構成される一次の拡散テンソル画像に基づいてその構造を強調した構造強調濃淡画像を形成してそれを画像表示手段に表示することをコンピュータを利用して実現するために用いるプログラムであって、
前記一次の拡散テンソル画像を構成する前記各点の一次の拡散テンソルについて固有値分解を行って得た前記各点の二次の拡散テンソルに基づいて前記二次の拡散テンソル画像を生成して保存する第1の機能と、
前記二次の拡散テンソル画像と同じ解像度の濃淡画像を生成する第２の機能と、
前記濃淡画像の各点において、濃淡値の時間微分を前記濃淡画像と前記各点の二次の拡散テンソルとを用いて算出する第３の機能と、
前記時間微分を設定された時間幅で積分して得た積分値を前記濃淡画像の対応する前記各点の濃淡値に加算して更新された濃淡画像を生成する第４の機能と、
前記更新された濃淡画像に対して、前記第３及び第４の機能を繰り返して実行することにより前記構造強調濃淡画像を作成する第５の機能とを前記コンピュータに実現させる構成を備えたことを特徴とする拡散テンソル場の構造強調濃淡画像表示用プログラム。
画像を表示する各点に拡散テンソルが与えられて構成される一次の拡散テンソル画像に基づいてその構造を強調した構造強調濃淡画像を形成してそれを画像表示手段に表示する拡散テンソル場の構造強調濃淡画像表示装置であって、
前記一次の拡散テンソル画像を構成する前記各点の一次の拡散テンソルについて固有値分解を行う固有値分解手段と、
前記固有値分解により得た前記各点についての固有値と固有ベクトルの画像を保存する固有値分解保存手段と、
制御値と前記固有値とから前記各点の固有値の重みを決定する重み決定手段と、
前記重みと、前記固有値と前記固有ベクトルとから二次の拡散テンソルを計算する二次の拡散テンソル計算手段と、
前記二次の拡散テンソルに基づく前記二次の拡散テンソル画像を保存する二次の拡散テンソル画像保存手段と、
前記二次の拡散テンソル画像と同じ解像度の濃淡画像を生成する濃淡画像生成手段と、
前記濃淡画像の各点において、濃淡値の時間微分を前記濃淡画像と前記各点の二次の拡散テンソルを用いて算出する時間微分算出手段と、
前記時間微分を設定された時間幅で積分して得た積分値を前記濃淡画像に加算して更新された濃淡画像を生成する更新濃淡画像生成手段と、
前記更新された濃淡画像の各点について、繰り返し前記時間微分算出手段及び前記更新濃淡画像生成手段による処理を行わせて前記構造強調濃淡画像を作成する構造強調濃淡画像作成手段とからなる拡散テンソル場の構造強調濃淡画像表示装置。