超音波造影装置
【課題】超音波造影剤を用いた画像診断において、超音波造影剤を含む生体組織のコントラスト分解能を高め、的確な診断に寄与することができる超音波造影装置を提供する。
【解決手段】超音波造影剤を含む生体組織に対して超音波パルスを送信し、任意の組織部位から発せられる散乱波形を受信するエコー受信部６と、散乱波形のパワースペクトルを算出するパワースペクトル算出手段７、８；１０、１１と、パワースペクトルの分散、標準偏差又は半値幅で定義されるスペクトル幅を特徴量として算出するスペクトル幅算出手段９、１２と、対象とする生体組織内全てにおいて求めて特徴量を求めその空間分布を得るように各手段を制御する制御手段と、該制御手段の制御により、対象とする生体組織内全てにおいて求めた空間分布を出力する出力手段を備え、超音波造影剤を含む生体組織を識別ないし可視化することを特徴とする超音波造影装置。
【解決手段】超音波造影剤を含む生体組織に対して超音波パルスを送信し、任意の組織部位から発せられる散乱波形を受信するエコー受信部６と、散乱波形のパワースペクトルを算出するパワースペクトル算出手段７、８；１０、１１と、パワースペクトルの分散、標準偏差又は半値幅で定義されるスペクトル幅を特徴量として算出するスペクトル幅算出手段９、１２と、対象とする生体組織内全てにおいて求めて特徴量を求めその空間分布を得るように各手段を制御する制御手段と、該制御手段の制御により、対象とする生体組織内全てにおいて求めた空間分布を出力する出力手段を備え、超音波造影剤を含む生体組織を識別ないし可視化することを特徴とする超音波造影装置。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、超音波造影剤を用いた画像診断において、超音波造影剤を含む生体組織のコントラスト分解能を高め、的確な診断に寄与する超音波造影装置に関するものである。
【背景技術】
【０００２】
従来の超音波造影剤又はコントラスト剤に係る超音波受信信号の処理法は数多く提案されており、フィルタ法、サブトラクション法、パルスインバージョン法の３つに分類することができる。特許文献１は、いわゆるフィルタ法について述べたものであり、超音波造影剤又はコントラスト剤の第２高調波発生が、周囲組織の第２高調波発生よりも高いことを利用し、超音波受信信号に対して基本周波数の２倍の帯域での信号成分をバンドパスフィルタにより抽出してこれを画像化する技術である。
【０００３】
特許文献２は、いわゆるサブトラクション法について述べたものであり、特に超音波音圧によって容易に崩壊する超音波造影剤を対象としており、超音波を照射する前後の画像の差分をとることにより、超音波造影剤の崩壊によって急激に低下輝度領域を逆に強調して画像化する技術である。
【０００４】
特許文献３は、いわゆるパルスインバージョン法について述べたものであり、考え方はフィルタ法と同様であるが、第２高調波成分を抽出するために、位相反転した２発以上の超音波を照射し、これらの差分をとることで実時間的に超音波造影剤を画像化する技術である。実際の超音波診断装置において、最も使用されている信号処理方式でもある。
【特許文献１】特開２００２−２２４１１０号公報「超音波造影描画装置」
【特許文献２】特開２００２−２８１６０「超音波診断装置」
【特許文献３】特開２００３−１３５４６７「超音波診断装置」
【発明の開示】
【発明が解決しようとする課題】
【０００５】
従来用いられてきた超音波造影剤の描画手法は、超音波造影剤の非線形振動による生じる高調波成分の検出（フィルタ法およびパルスインバージョン法）又は超音波造影剤の崩壊現象を利用するものであった（サブトラクション法）。
【０００６】
非線形振動は超音波造影剤のみならず正常組織又は周囲軟組織においても生じる。フィルタ法は、周波数領域において、高調波成分のみをバンドパスフィルタにより抽出し、この成分を画像化するものであるが、高調波成分には、少なからず基本波成分が重畳しているため、フィルタ法により完全に高調波成分のみを抽出し、かつ超音波造影剤のコントラストを著しく向上させることには限界があった。
【０００７】
パルスインバージョン法はこうした欠点が低減されていると考えられているが、周波数領域で考えると本質的にフィルタ法と同等な問題を抱えているものと推察できる。さらに、パルスインバージョン法は、原理的に流れが存在する場合のコントラスト向上には限界があることが指摘されており、超音波造影剤は主として血流中にて用いられることを鑑みれば、更なる改善を要することは自明である。
【０００８】
またサブトラクション法は、音圧による超音波造影剤の崩壊現象を利用した画像化法であるが、これまで主流であった超音波造影剤であるレボビストにおいては、ある程度有効であったが、近年登場した長い造影効果を有するソナゾイドにおいては、音圧に対して頑健であり、崩壊しにくいために、今後、発展が見込める超音波造影手法とは必ずしも言えない。
【０００９】
上記のように、現在、超音波造影剤の描画手法としては、超音波造影剤が発する高調波の検出と、超音波造影剤の崩壊の検出が主流であるが、各手法による検出能ないしコントラスト分解能の向上は次第に飽和状態に近づきつつあり、現状の延長上では著しいコントラスト分解能の向上を期待し難い。従って、臨床診断の場では、超音波造影剤のさらなるコントラスト分解能の向上が期待されており、高調波発生及び崩壊現象以外の物理現象を利用した超音波造影剤の検出及び画像化手法が求められている。
【００１０】
本発明は、このような従来技術の実状に鑑みてなされたもので、超音波造影剤を用いた画像診断において、超音波造影剤を含む生体組織のコントラスト分解能を高め、的確な診断に寄与することができる超音波造影装置を提供することを課題とする。
【課題を解決するための手段】
【００１１】
本発明は、上記課題を解決するために、第１には、超音波造影剤を含む生体組織に対して超音波パルスを送受信する超音波プローブを有し、該超音波プローブにより超音波を送信した後、任意の組織部位から発せられる散乱波形を、該プローブを介して受信する散乱波形取得手段と、該散乱波形取得手段が取得した散乱波形のパワースペクトルを算出するパワースペクトル算出手段と、該パワースペクトル算出手段が算出したパワースペクトルの分散、標準偏差又は半値幅で定義されるスペクトル幅を特徴量として算出するスペクトル幅算出手段と、該スペクトル幅算出手段が算出した特徴量を、対象とする生体組織内全てにおいて求めてその空間分布を得るように各手段を制御する制御手段と、該制御手段の制御により、対象とする生体組織内全てにおいて求めた空間分布を出力する出力手段を備え、超音波造影剤を含む生体組織を識別ないし可視化することを特徴とする超音波造影装置を提供する。
【００１２】
また、第２には、上記第１の発明において、該パワースペクトル算出手段は、任意の組織部位から発せられる散乱波形を、その組織部位の深さに対応する時間位置及び時間区間にのみ非零の値を持つステップ状ゲートと受信された散乱波形との乗算を行うことにより抽出する局所波形ゲート抽出部と、該局所波形ゲート抽出部により抽出された散乱波形に対してフーリエ変換を行ってパワースペクトルを算出する抽出波形スペクトル計算部を備え、該スペクトル幅算出手段は、超音波造影剤を含む生体組織内における該スペクトル幅を、該抽出波形スペクトル計算部により算出されたパワースペクトルに対して重み積分処理ないし該パワースペクトルのピークの半値幅検出を行うことにより算出することを特徴とする超音波造影装置を提供する。
【００１３】
また、第３には、上記第１の発明において、該パワースペクトル算出手段は、散乱波形の中心周波数成分を含む、所定個数の特定周波数成分を持つ正弦波及び余弦波を各々散乱波形に乗算する特定周波数成分乗算部と、該特定周波数成分乗算部により得られた所定個数の乗算波形に対して平均化処理を行ってパワースペクトルを離散化した所定個数の点を求める特定周波数パワー抽出部を備え、該スペクトル幅算出手段は、超音波造影剤を含む生体組織内における該スペクトル幅を、該特定周波数パワー抽出部が求めた所定個数のパワースペクトル離散点から、曲線近似ないし多項式近似を行うことにより推定されたパワースペクトルに対して重み積分処理ないし該パワースペクトルのピークの半値幅検出を行うことにより算出することを特徴とする超音波造影装置を提供する。
【００１４】
また、第４には、上記第３の発明において、該スペクトル幅算出手段は、超音波造影剤を含む生体組織内における該スペクトル幅を、該特定周波数パワー抽出部が求めた所定個数のパワースペクトル離散点を用いて台形公式、シンプソン公式及びガウス積分公式のうちのいずれかの数値積分法を適用し、時間領域において直接算出することを特徴とする超音波造影装置を提供する。
【００１５】
また、第５には、上記第１の発明において、該制御手段が、超音波造影剤を含む生体組織内における該スペクトル幅から、その逆数を計算することにより超音波減衰振動における時定数を算出し、その空間分布を呈示ないし画像化することにより、超音波造影剤を含む生体組織を識別ないし可視化することを特徴とする超音波造影装置を提供する。
【００１６】
また、第６には、上記第１の発明において、該制御手段が、粒径が均一な超音波造影剤を含む生体組織内における該スペクトル幅を、気泡挙動の方程式により予め理論的に得たスペクトル幅対粒径分布分散又は標準偏差の特性曲線の関数に代入することで粒径分布の分散又は標準偏差を算出し、その空間分布を呈示ないし画像化することにより、超音波造影剤を含む生体組織を識別ないし可視化することを特徴とする超音波造影装置を提供する。
【００１７】
さらに、第７には、上記第１、５又は６の発明において、該スペクトル幅算出手段が、超音波造影剤を含む生体組織から発せられる散乱波形より得られるスペクトル幅と、周辺の生体組織から発せられる散乱波形より得られるスペクトル幅を各々検出し、該制御手段が、それらの対数、あるいは時定数空間分布の対数、あるいは粒径分布の分散又は標準偏差の対数を算出し、その空間分布を呈示ないし画像化することにより、超音波造影剤を含む生体組織を識別ないし可視化することを特徴とする超音波造影装置を提供する。
【発明の効果】
【００１８】
本発明では、従来の高調波検出や崩壊現象の検出を利用した手法とは異なり、スペクトル幅を利用した超音波造影剤の検出ないしイメージングという物理現象の検出を利用しており、従来の超音波造影装置に比べて、周囲組織との高いコントラスト比が達成できる。これにより、本発明によれば、超音波造影剤を用いた診断に対して、近年登場した長い造影効果を有するソナゾイドにも十分対応でき、従来法よりも高いコントラスト強調を伴う超音波造影装置を提供することが可能となる。
【発明を実施するための最良の形態】
【００１９】
以下に、本発明の一実施形態を説明する。
【００２０】
図１は、本発明の一実施形態に係る超音波造影装置における超音波受信信号の処理部分構成を示す全体図である。同図において、１は超音波プローブ、２は周辺の生体組織、３は超音波造影剤を含む生体組織、４は超音波送受信伝搬経路、５は走査線を示す。超音波プローブ１から超音波伝搬経路４に沿って超音波を送受信し、Ｎ回の走査を繰り返しながら、生体組織断面からの散乱エコーを、エコー受信部６にて受信するようになっている。受信した散乱エコーは、その局所的なスペクトル幅を算出するために、図１に示すような標準ルーチン又は高速ルーチンのシステムで処理が行われる。ここでスペクトル幅とは、パワースペクトルの分散、標準偏差又は半値幅で定義される特徴量である。
【００２１】
標準ルーチンは、局所的な高速フーリエ変換（ＦＦＴ）を行ってパワースペクトルを求め、そのパワースペクトルを用いてスペクトル幅を算出するルーチンであり、計算時間はある程度かかるが、より正確なスペクトル幅計算が可能である。一方、高速ルーチンは、標準ルーチンにおけるスペクトル幅計算を簡素化して近似値を求めるルーチンであり、計算時間を著しく短縮することが可能である。
【００２２】
標準ルーチンと高速ルーチンは、用途により使い分けがなされる。例えば、リアルタイム処理において、あまり高いフレームレートが必要でなく、かつ正確なスペクトル幅計算が必要な場合は、得られた散乱エコーに対して、標準ルーチンを用いて処理することができる。
【００２３】
一方、リアルタイム処理において、高いフレームレートが必要で、高速にスペクトル幅を算出する必要がある場合は、得られた散乱エコーに対して、高速ルーチンを用いて処理することができる。
【００２４】
標準ルーチンと高速ルーチンは、このように用途により使い分けられるが、適用はどちらか一方を選択しなければならないわけではなく、両者を組み合わせて用いてもよい。
【００２５】
次に、標準ルーチンと高速ルーチンを実行するための各構成を説明する。まず標準ルーチンは、局所波形ゲート抽出部７、抽出波形スペクトル計算部８、スペクトル幅計算部９からなる。
【００２６】
局所波形ゲート抽出部７の詳細は図２に示すものとなる。ここでは、１走査線における受信エコー１５に対し、ゲート１６を各々乗算する。図２の横方向は時間に相当する。ゲート１６は、超音波伝搬経路４上における所望の深さからの散乱エコーを抽出するためのもので、抽出する深さに相当する時間位置及び時間区間にのみ非零（例えば、振幅１）を持ち、その他の振幅が０であり、散乱エコーと乗算することにより、抽出波形１７が得られる。ゲート１６において（１）〜（M）に示すようにゲートの時間位置を変えながら、この処理を全ての深さでM回繰り返すことにより、各深さからの散乱エコーが抽出される。ここで、パワースペクトルの形状は、ゲート１６における時間区間幅の影響を受ける。全ての深さで相対化して均等な条件下で後のスペクトル幅の算出が可能となるよう、ゲート１６の時間区間幅は一定値に固定する。
【００２７】
抽出波形スペクトル計算部８、スペクトル幅計算部９の詳細は、図３に示すものとなる。まず抽出波形スペクトル計算部８では、ゲート（ｉ）（ここでｉ＝１〜M）を用いて抽出された散乱エコー（図２の１７）に対し、高速フーリエ変換（ＦＦＴ）を適用する。そのため抽出波形スペクトル計算部８はM個のＦＦＴ演算部１８を有している。ＦＦＴ演算部１８のＦＦＴ出力をスペクトル幅計算部９にてパワースペクトルに変換し、パワースペクトルに対して重み積分処理を行って、スペクトル幅を算出する。さらに、深さ毎に算出されたスペクトル幅の値を、スペクトル幅プロファイル構成部１３おいて集合させ、合成すれば、超音波伝搬経路４上におけるスペクトル幅のプロファイルが得られ、これを図１におけるＮ本の走査線５において繰り返せば、空間分布が得られる。この空間分布を、最終的に、空間分布を呈示もしくは画像化させ、例えば表示部１４にそれを表示させることにより、超音波造影剤を含む生体組織を識別ないし可視化することが可能となる。
【００２８】
この一連の過程は図４のように数式表記される。すなわち、散乱エコーをｘ（ｔ）としたとき、ゲートを乗算し、深さ毎に算出されるパワースペクトルは、次のように表される。
【００２９】
【数１】
【００３０】
ここで、Ｐ（ω，ｔ）は局所的に抽出された波形のパワースペクトル、Ｔは局所抽出時間区間（ゲート時間区間幅）である。式（１）により、深さ毎にパワースペクトルが算出され、以下の重み積分処理を行って、スペクトル幅のプロファイルｂｗ（ｔ）が算出される。
【００３１】
【数２】
【００３２】
式（２）に示すプロファイルは、１走査線において得られるものであるため、画像を構成する全ての走査線５においてこの処理を行えば、スペクトル幅の空間分布が得られることになる。図３のパワースペクトル幅集合により構成されたプロファイル２０は、中間領域に超音波造影剤を含む生体組織が位置し、この部分のスペクトル幅が、周辺組織のスペクトル幅よりも小さくなることで、コントラスト分解能が向上する例示であり、この逆となることも有り得る。また、ここでのスペクトル幅の算出において、式（２）に示す重み積分処理は一例であり、別の方法としては、パワースペクトルのピーク値を求め、その半値となる周波数幅を求めてもよい。
【００３３】
次に、高速ルーチンの構成を説明する。高速ルーチンは、標準ルーチンを簡便化して、より高速にスペクトル幅の算出を行うためのルーチンである。高速ルーチンは、特定周波数成分乗算部１０、特定周波数パワー抽出部１１、スペクトル幅計算部（近似）１２からなる。
【００３４】
特定周波数成分乗算部１０の詳細は、図５に示すものとなる。ここでは、１走査線における受信エコー１５に対し、特定周波数成分２１を各々乗算する。ここで、特定周波数成分とは、受信エコー１５の中心角周波数ω０を起点とし、その前後、すなわち、ω０−ｎΔω〜ω０＋ｎΔωの範囲の（２ｎ＋１）個の周波数成分を指す。ｎは整数、Δωは角周波数の刻み幅である。これらの周波数成分を持つ正弦波及び余弦波をそれぞれ受信エコー１５に乗算することにより、受信エコー１５を周波数変調した複素数波形２２が得られる。
【００３５】
次に、図６の特定周波数パワー抽出部１１では、複素数波形２２を重み１の有限幅を持つ窓内で平滑化（すなわち局所的な積分を実行）し、得られた複素数成分のパワー（複素数の絶対値。実部と虚部の二乗和）を図６の特定周波数パワー抽出部１１に示すように周波数軸上にプロットする。すわわち、ω０−ｎΔωの周波数成分を持つ正弦波及び余弦波を乗算して得られた複素数波形から得られたパワーをω０−ｎΔωの位置にプロットし、これを（２ｎ＋１）回繰り返すと、パワースペクトルの離散点が得られる。その後、スペクトル幅計算部（近似）１２では、パワースペクトルの離散点を基に、曲線近似又は多項式近似によりパワースペクトル全体の分布形を推定する。パワースペクトルの分布形を推定した後は、標準ルーチンと同様、重み積分処理又はパワースペクトルピーク値の半値幅検出により、スペクトル幅を算出することができる。さらに、深さ毎に算出されたスペクトル幅の値を、スペクトル幅プロファイル構成部１３において集合させ、合成すれば、超音波伝搬経路４上におけるスペクトル幅のプロファイルが得られ、これを図１におけるＮ本の走査線５において繰り返せば、標準ルーチンと同様のスペクトル幅の空間分布が得られる。このスペクトル幅の空間分布を最終的に、空間分布を呈示もしくは画像化させ、表示部１４に表示させることにより、超音波造影剤を含む生体組織を識別ないし可視化することが可能となる。
【００３６】
この一連の過程は図７のように数式表記される。すなわち、散乱エコーをｘ（ｔ）としたとき、特定周波数ω０＋ｉΔω（ただし、ｉ＝−ｎ・・・０・・・ｎ，ｎ＞０）を持つ正弦波及び余弦波を乗算し、時間幅２Ｔを持つ窓内で単純平滑化（重み関数の係数が１の平滑化。局所的な積分に相当）すると、以下の式が得られる。
【００３７】
【数３】
【００３８】
この式は、式（１）を離散化したものに相当し、結局、式（３）は式（１）と等価な演算となる。式（３）を用いて計算されたパワースペクトルの離散点Ｐｉ（ｔ）を用いたパワースペクトル分布形の推定における曲線近似の一例として、J次多項式近似を行う場合、パワースペクトルの分布形は、
【００３９】
【数４】
【００４０】
のようにモデル化ないし近似されるので、得られた離散点データ（ωｉ，Ｐｉ（ｔ））を用いた最小二乗法により、係数行列{ａｊ}を求めれば、パワースペクトル分布形の推定が可能となり、スペクトル幅の算出が行われる。
【００４１】
ただし、パワースペクトル分布形の推定において、式（４）のモデルは一例であり、例えば、超音波パルス送信時のパワースペクトル分布関数モデルや、ガウス分布関数モデルなどパワースペクトル分布の近似関数モデルを利用することもできる。
【００４２】
上記では、式（３）により計算されたパワースペクトルの離散点Ｐｉ（ｔ）を用いたパワースペクトル分布形を、例えば式（４）の多項式近似モデルを適用して求め、その後、式（２）に示す重み積分処理を実行して、スペクトル幅の算出を精度よく行うための方法を述べたが、計算負荷の増大が懸念されるケースも生じ得る。このような場合への対処、すなわちスペクトル幅の算出において計算精度を多少犠牲にしても計算速度向上が望まれる場合、上記パワースペクトルの分布形計算処理を省略して、式（３）により計算されたパワースペクトルの離散点Ｐｉ（ｔ）に対して数値積分法を適用して、式（２）の重み積分を直接計算することができる。ここでの数値積分には、台形公式及びシンプソン公式及びガウス積分公式を適用できる。また、このときパワースペクトルの離散点Ｐｉ（ｔ）は時間領域で既に得られているため、上記数値積分法も時間領域において適用可能であり、すなわち、時間領域にて直接スペクトル幅を算出できる。
【００４３】
算出されたスペクトル幅は、そのまま空間分布として画像呈示できる他、制御手段において、該スペクトル幅から超音波減衰振動における時定数、または超音波造影剤における粒径分布の広がり分散又は標準偏差を算出し、呈示することもできる。
【００４４】
スペクトル幅からの減衰振動における時定数を得るために、まず、受信波形ｘ（ｔ）を次式で表す。
【００４５】
【数５】
【００４６】
ここで、Ａはｔ＝０での振幅、ｔは時間、Tは時定数、ω0は受信波の中心角周波数を表す。この受信波のパワースペクトルＰ（ω）は、次式のように表される。
【００４７】
【数６】
【００４８】
このパワースペクトルの半値幅は、２／Ｔとなり、すなわち、時定数の逆数の２倍に等しくなる。実際の受信波は、多くの超音波造影剤からの散乱波を全て足し合わせた総和になるため、実際のパワースペクトルも、個々の散乱波から得られるパワースペクトルの総和として得られる。しかしながら粒径分布が均一に近ければ、上記パワースペクトルの総和もまた、式（６）のように近似することができる。従って、パワースペクトルのスペクトル幅における半値幅の逆数を計算することにより、減衰振動の時定数もまた算出することができる。
【００４９】
一方、超音波造影剤における粒径分布の分散又は標準偏差は、気泡挙動の方程式により予め理論的に得たスペクトル幅対粒径分布分散または標準偏差の特性曲線の関数にスペクトル幅を代入することで、算出することができる。ここでは、気泡挙動の方程式に基づくスペクトル幅対粒径分布分散または標準偏差の特性曲線の関数について述べる。まず、超音波造影剤の共振周波数はその径に大きく依存し、粒径のばらつきが大きければその分スペクトル幅も広がり、また逆に粒径のばらつきが小さければその分スペクトル幅も狭くなることを利用し、算出されたスペクトル幅から粒径分布の広がり程度、すなわち分散又は標準偏差を評価しようとすることが、ここでの基本的な概念である。しかしながら、スペクトル幅として半値幅、粒径分布広がりとして標準偏差に着目するとき、散乱波パワースペクトルの半値幅と粒径分布の標準偏差との関係を表す関数は、超音波造影剤の平均粒径に依存し、さらには線形関係にない。粒径分布をガウス分布であると仮定すると、その標準偏差をσは一般に次のように記述できる。
【００５０】
【数７】
【００５１】
ここで、ｆは非線形関数を意味し、またＲｍとＦＷＨＭ（full width half maximum）はそれぞれ平均粒径と散乱波パワースペクトルの半値幅を表す。平均粒径は、超音波造影剤の粒径として現実的な３〜４μｍの定数と見なせるので、このときのｆの関数形が分かれば、散乱波パワースペクトルの半値幅ＦＷＨＭを式（７）に代入することにより、粒径分布の標準偏差σが得られることになる。ここで、ｆの関数形を予め求めておくための気泡挙動の方程式としては、よく知られたRayleigh-Plesset方程式を用いることができる。まずこの方程式に基づいて粒径対共振周波数を得る。これは、反比例的な関係式となる。次に、粒径分布がガウス分布であるとし、この分布を、上記の粒径対共振周波数の関係式に代入すると、共振周波数成分の集合としてのパワースペクトルが得られる。このパワースペクトルの半値幅は、ここでの粒径分布であるガウス分布の標準偏差と１対１の対応関係にある。従って、粒径分布におけるガウス分布の標準偏差を順次少しずつ変えたときに、それに対応する散乱波パワースペクトル半値幅の変化を求めてゆくと、曲線適合により、式（７）におけるｆの関数形を定めることができる。我々の検討によると、ｆは、σ＝（ＦＷＨＭ）ｑのような関数として近似できる（ただし、ｑは累乗における指数であり、平均粒径に応じて決まり、１未満の実数をとる。）。そして、このｑは、上記の曲線適合により決定することができるため、パワースペクトルの半値幅を代入することにより、粒径分布の標準偏差を算出することができる。さらに上記の累乗関数は、パワースペクトルの半値幅として現実的な０〜２ＭＨｚ付近までのレンジでは、急峻な勾配特性を持つため、僅かなパワースペクトル半値幅の違いが増幅され、粒径分布の標準偏差として得られることになる。これにより結果的に、周辺組織と、超音波造影剤を含む生体組織とのコントラスト分解能が向上することになる。
【００５２】
ここで、超音波造影剤を含む生体組織と周辺の生体組織における減衰振動特性の相違、あるいは周辺の生体組織の周波数特性と超音波造影剤の均一性に起因する周波数特性の相違のために、超音波造影剤を含む生体組織から発せられる散乱波形より得られるスペクトル幅と周辺の生体組織から発せられる散乱波形より得られるスペクトル幅はそれぞれ異なる値を有する。得られたスペクトル幅の生体組織内における空間分布自体をそのまま画像化すれば、超音波造影剤を含む生体組織と周辺の生体組織のコントラスト差は視認することが可能になる。しかしながら、このコントラスト差が顕著であり、表示画像のコントラスト分解能（例えば２５６階調）を超える場合は、得られたスペクトル幅を対数化してダイナミックレンジを圧縮することで、欠落なく全ての対象領域を画像化できる。また上記の減衰振動における時定数または粒径分布の標準偏差についても、これらのコントラスト差が顕著であり、表示画像のコントラスト分解能（例えば２５６階調）を超える場合は、得られた時定数ないし粒径分布の標準偏差を対数化してダイナミックレンジを圧縮することで、欠落なく全ての対象領域を画像化できる。
【００５３】
（実験例）
スペクトル幅算出によるコントラスト分解能向上を評価するため、上記で述べた装置を用い、造影実験を実施した。ここでは、超音波造影剤として中空マイクロカプセル及び市販の超音波造影剤であるレボビストを使用し、また参照波としてこれら造影剤を含むＰＰ（ポリプロピレン）容器の壁からの反射波を採取した。スペクトル幅算出の妥当性検証、特に、粒径が均一の超音波造影剤に対する検証には、本実験に用いた中空マイクロカプセルはより均一な粒径で製造可能である利点を利用することができる。Rayleigh-Plesset方程式に基づくと、気泡の共振周波数は、その径に大きく依存することになる。従って、造影剤の集合体において粒径のばらつきが大きければその分多様な共振周波数が存在することとなり、すなわちスペクトルにおけるメインローブのスペクトル幅は広がると考えられる。逆に、造影剤の集合体において粒径のばらつきが小さければ線スペクトルに近づき、スペクトル幅は狭まると考えられる。周囲組織が広帯域となる場合、このようなスペクトル幅の違いを局所的に算出できれば、コントラスト分解能の向上には有用である。中空マイクロカプセルの共振周波数に近い１ＭＨｚを用いたとき、パワースペクトル分布を図８（ａ）に、得られたスペクトル幅を図８（ｂ）に示す。なお、散乱波の受信も１ＭＨｚのプローブで行われたため、明確な高調波はここでは観察されていない。図８（ｂ）に示すように、中空マイクロカプセルは、参照となるＰＰ容器壁、すなわち疑似周囲組織に比べて有意に低下し、さらに粒径のばらつきが比較的大きな従来の超音波造影剤であるレボビストよりも半値幅は小さくなっている。これは、本手法により粒径分布の広がりが妥当に検出され得ることを示しており、結果的に、このようなスペクトル幅の差を呈示することにより、従来の超音波造影剤イメージング法よりも高いコントラスト分解能にもまた到達できる。
【図面の簡単な説明】
【００５４】
【図１】本発明の一実施形態に係る超音波造影装置における超音波受信信号の処理部分の全体構成を示す図である。
【図２】局所波形ゲート抽出部７の処理内容を説明した図である。
【図３】抽出波形スペクトル計算部８、スペクトル幅計算部９、スペクトル幅プロファイル構成部１３に至る処理内容を説明した図である。
【図４】図１における標準ルーチンの処理流れを数式を用いて説明した図である。
【図５】特定周波数成分乗算部１０の処理内容を説明した図である。
【図６】特定周波数パワー抽出部１１、スペクトル幅計算部（近似）１２、スペクトル幅プロファイル構成部１３に至る処理内容を説明した図である。
【図７】図１における高速ルーチンの処理流れを数式を用いて説明した図である。
【図８】本発明の装置によりスペクトル幅を算出した実験例である。
【符号の説明】
【００５５】
１ 超音波プローブ
２ 周辺組織
３ 超音波造影剤
４ 超音波送受信伝搬経路
５ 走査線
６ エコー受信部
７ （標準ルーチン）局所波形ゲート抽出部
８ （標準ルーチン）局所波形スペクトル幅計算部
９ （標準ルーチン）スペクトル幅計算部
１０ （高速ルーチン）特定周波数成分乗算部
１１ （高速ルーチン）特定周波数パワー抽出部
１２ （高速ルーチン）スペクトル幅計算部（近似値）
１３ スペクトル幅プロファイル構成部
１４ 表示部
１５ 受信エコー／走査線
１６ ゲート
１７ 抽出波形
１８ ＦＦＴ演算部
１９ パワースペクトル幅
２０ パワースペクトル幅集合により構成されたプロファイル
２１ 特定周波数成分
２２ 乗算された波形
２３ 離散的な特定周波数パワー値（丸）と近似されるパワースペクトル分布
２４ パワースペクトル分布を曲線近似により推定し、スペクトル幅を計算するルーチン
【技術分野】
【０００１】
本発明は、超音波造影剤を用いた画像診断において、超音波造影剤を含む生体組織のコントラスト分解能を高め、的確な診断に寄与する超音波造影装置に関するものである。
【背景技術】
【０００２】
従来の超音波造影剤又はコントラスト剤に係る超音波受信信号の処理法は数多く提案されており、フィルタ法、サブトラクション法、パルスインバージョン法の３つに分類することができる。特許文献１は、いわゆるフィルタ法について述べたものであり、超音波造影剤又はコントラスト剤の第２高調波発生が、周囲組織の第２高調波発生よりも高いことを利用し、超音波受信信号に対して基本周波数の２倍の帯域での信号成分をバンドパスフィルタにより抽出してこれを画像化する技術である。
【０００３】
特許文献２は、いわゆるサブトラクション法について述べたものであり、特に超音波音圧によって容易に崩壊する超音波造影剤を対象としており、超音波を照射する前後の画像の差分をとることにより、超音波造影剤の崩壊によって急激に低下輝度領域を逆に強調して画像化する技術である。
【０００４】
特許文献３は、いわゆるパルスインバージョン法について述べたものであり、考え方はフィルタ法と同様であるが、第２高調波成分を抽出するために、位相反転した２発以上の超音波を照射し、これらの差分をとることで実時間的に超音波造影剤を画像化する技術である。実際の超音波診断装置において、最も使用されている信号処理方式でもある。
【特許文献１】特開２００２−２２４１１０号公報「超音波造影描画装置」
【特許文献２】特開２００２−２８１６０「超音波診断装置」
【特許文献３】特開２００３−１３５４６７「超音波診断装置」
【発明の開示】
【発明が解決しようとする課題】
【０００５】
従来用いられてきた超音波造影剤の描画手法は、超音波造影剤の非線形振動による生じる高調波成分の検出（フィルタ法およびパルスインバージョン法）又は超音波造影剤の崩壊現象を利用するものであった（サブトラクション法）。
【０００６】
非線形振動は超音波造影剤のみならず正常組織又は周囲軟組織においても生じる。フィルタ法は、周波数領域において、高調波成分のみをバンドパスフィルタにより抽出し、この成分を画像化するものであるが、高調波成分には、少なからず基本波成分が重畳しているため、フィルタ法により完全に高調波成分のみを抽出し、かつ超音波造影剤のコントラストを著しく向上させることには限界があった。
【０００７】
パルスインバージョン法はこうした欠点が低減されていると考えられているが、周波数領域で考えると本質的にフィルタ法と同等な問題を抱えているものと推察できる。さらに、パルスインバージョン法は、原理的に流れが存在する場合のコントラスト向上には限界があることが指摘されており、超音波造影剤は主として血流中にて用いられることを鑑みれば、更なる改善を要することは自明である。
【０００８】
またサブトラクション法は、音圧による超音波造影剤の崩壊現象を利用した画像化法であるが、これまで主流であった超音波造影剤であるレボビストにおいては、ある程度有効であったが、近年登場した長い造影効果を有するソナゾイドにおいては、音圧に対して頑健であり、崩壊しにくいために、今後、発展が見込める超音波造影手法とは必ずしも言えない。
【０００９】
上記のように、現在、超音波造影剤の描画手法としては、超音波造影剤が発する高調波の検出と、超音波造影剤の崩壊の検出が主流であるが、各手法による検出能ないしコントラスト分解能の向上は次第に飽和状態に近づきつつあり、現状の延長上では著しいコントラスト分解能の向上を期待し難い。従って、臨床診断の場では、超音波造影剤のさらなるコントラスト分解能の向上が期待されており、高調波発生及び崩壊現象以外の物理現象を利用した超音波造影剤の検出及び画像化手法が求められている。
【００１０】
本発明は、このような従来技術の実状に鑑みてなされたもので、超音波造影剤を用いた画像診断において、超音波造影剤を含む生体組織のコントラスト分解能を高め、的確な診断に寄与することができる超音波造影装置を提供することを課題とする。
【課題を解決するための手段】
【００１１】
本発明は、上記課題を解決するために、第１には、超音波造影剤を含む生体組織に対して超音波パルスを送受信する超音波プローブを有し、該超音波プローブにより超音波を送信した後、任意の組織部位から発せられる散乱波形を、該プローブを介して受信する散乱波形取得手段と、該散乱波形取得手段が取得した散乱波形のパワースペクトルを算出するパワースペクトル算出手段と、該パワースペクトル算出手段が算出したパワースペクトルの分散、標準偏差又は半値幅で定義されるスペクトル幅を特徴量として算出するスペクトル幅算出手段と、該スペクトル幅算出手段が算出した特徴量を、対象とする生体組織内全てにおいて求めてその空間分布を得るように各手段を制御する制御手段と、該制御手段の制御により、対象とする生体組織内全てにおいて求めた空間分布を出力する出力手段を備え、超音波造影剤を含む生体組織を識別ないし可視化することを特徴とする超音波造影装置を提供する。
【００１２】
また、第２には、上記第１の発明において、該パワースペクトル算出手段は、任意の組織部位から発せられる散乱波形を、その組織部位の深さに対応する時間位置及び時間区間にのみ非零の値を持つステップ状ゲートと受信された散乱波形との乗算を行うことにより抽出する局所波形ゲート抽出部と、該局所波形ゲート抽出部により抽出された散乱波形に対してフーリエ変換を行ってパワースペクトルを算出する抽出波形スペクトル計算部を備え、該スペクトル幅算出手段は、超音波造影剤を含む生体組織内における該スペクトル幅を、該抽出波形スペクトル計算部により算出されたパワースペクトルに対して重み積分処理ないし該パワースペクトルのピークの半値幅検出を行うことにより算出することを特徴とする超音波造影装置を提供する。
【００１３】
また、第３には、上記第１の発明において、該パワースペクトル算出手段は、散乱波形の中心周波数成分を含む、所定個数の特定周波数成分を持つ正弦波及び余弦波を各々散乱波形に乗算する特定周波数成分乗算部と、該特定周波数成分乗算部により得られた所定個数の乗算波形に対して平均化処理を行ってパワースペクトルを離散化した所定個数の点を求める特定周波数パワー抽出部を備え、該スペクトル幅算出手段は、超音波造影剤を含む生体組織内における該スペクトル幅を、該特定周波数パワー抽出部が求めた所定個数のパワースペクトル離散点から、曲線近似ないし多項式近似を行うことにより推定されたパワースペクトルに対して重み積分処理ないし該パワースペクトルのピークの半値幅検出を行うことにより算出することを特徴とする超音波造影装置を提供する。
【００１４】
また、第４には、上記第３の発明において、該スペクトル幅算出手段は、超音波造影剤を含む生体組織内における該スペクトル幅を、該特定周波数パワー抽出部が求めた所定個数のパワースペクトル離散点を用いて台形公式、シンプソン公式及びガウス積分公式のうちのいずれかの数値積分法を適用し、時間領域において直接算出することを特徴とする超音波造影装置を提供する。
【００１５】
また、第５には、上記第１の発明において、該制御手段が、超音波造影剤を含む生体組織内における該スペクトル幅から、その逆数を計算することにより超音波減衰振動における時定数を算出し、その空間分布を呈示ないし画像化することにより、超音波造影剤を含む生体組織を識別ないし可視化することを特徴とする超音波造影装置を提供する。
【００１６】
また、第６には、上記第１の発明において、該制御手段が、粒径が均一な超音波造影剤を含む生体組織内における該スペクトル幅を、気泡挙動の方程式により予め理論的に得たスペクトル幅対粒径分布分散又は標準偏差の特性曲線の関数に代入することで粒径分布の分散又は標準偏差を算出し、その空間分布を呈示ないし画像化することにより、超音波造影剤を含む生体組織を識別ないし可視化することを特徴とする超音波造影装置を提供する。
【００１７】
さらに、第７には、上記第１、５又は６の発明において、該スペクトル幅算出手段が、超音波造影剤を含む生体組織から発せられる散乱波形より得られるスペクトル幅と、周辺の生体組織から発せられる散乱波形より得られるスペクトル幅を各々検出し、該制御手段が、それらの対数、あるいは時定数空間分布の対数、あるいは粒径分布の分散又は標準偏差の対数を算出し、その空間分布を呈示ないし画像化することにより、超音波造影剤を含む生体組織を識別ないし可視化することを特徴とする超音波造影装置を提供する。
【発明の効果】
【００１８】
本発明では、従来の高調波検出や崩壊現象の検出を利用した手法とは異なり、スペクトル幅を利用した超音波造影剤の検出ないしイメージングという物理現象の検出を利用しており、従来の超音波造影装置に比べて、周囲組織との高いコントラスト比が達成できる。これにより、本発明によれば、超音波造影剤を用いた診断に対して、近年登場した長い造影効果を有するソナゾイドにも十分対応でき、従来法よりも高いコントラスト強調を伴う超音波造影装置を提供することが可能となる。
【発明を実施するための最良の形態】
【００１９】
以下に、本発明の一実施形態を説明する。
【００２０】
図１は、本発明の一実施形態に係る超音波造影装置における超音波受信信号の処理部分構成を示す全体図である。同図において、１は超音波プローブ、２は周辺の生体組織、３は超音波造影剤を含む生体組織、４は超音波送受信伝搬経路、５は走査線を示す。超音波プローブ１から超音波伝搬経路４に沿って超音波を送受信し、Ｎ回の走査を繰り返しながら、生体組織断面からの散乱エコーを、エコー受信部６にて受信するようになっている。受信した散乱エコーは、その局所的なスペクトル幅を算出するために、図１に示すような標準ルーチン又は高速ルーチンのシステムで処理が行われる。ここでスペクトル幅とは、パワースペクトルの分散、標準偏差又は半値幅で定義される特徴量である。
【００２１】
標準ルーチンは、局所的な高速フーリエ変換（ＦＦＴ）を行ってパワースペクトルを求め、そのパワースペクトルを用いてスペクトル幅を算出するルーチンであり、計算時間はある程度かかるが、より正確なスペクトル幅計算が可能である。一方、高速ルーチンは、標準ルーチンにおけるスペクトル幅計算を簡素化して近似値を求めるルーチンであり、計算時間を著しく短縮することが可能である。
【００２２】
標準ルーチンと高速ルーチンは、用途により使い分けがなされる。例えば、リアルタイム処理において、あまり高いフレームレートが必要でなく、かつ正確なスペクトル幅計算が必要な場合は、得られた散乱エコーに対して、標準ルーチンを用いて処理することができる。
【００２３】
一方、リアルタイム処理において、高いフレームレートが必要で、高速にスペクトル幅を算出する必要がある場合は、得られた散乱エコーに対して、高速ルーチンを用いて処理することができる。
【００２４】
標準ルーチンと高速ルーチンは、このように用途により使い分けられるが、適用はどちらか一方を選択しなければならないわけではなく、両者を組み合わせて用いてもよい。
【００２５】
次に、標準ルーチンと高速ルーチンを実行するための各構成を説明する。まず標準ルーチンは、局所波形ゲート抽出部７、抽出波形スペクトル計算部８、スペクトル幅計算部９からなる。
【００２６】
局所波形ゲート抽出部７の詳細は図２に示すものとなる。ここでは、１走査線における受信エコー１５に対し、ゲート１６を各々乗算する。図２の横方向は時間に相当する。ゲート１６は、超音波伝搬経路４上における所望の深さからの散乱エコーを抽出するためのもので、抽出する深さに相当する時間位置及び時間区間にのみ非零（例えば、振幅１）を持ち、その他の振幅が０であり、散乱エコーと乗算することにより、抽出波形１７が得られる。ゲート１６において（１）〜（M）に示すようにゲートの時間位置を変えながら、この処理を全ての深さでM回繰り返すことにより、各深さからの散乱エコーが抽出される。ここで、パワースペクトルの形状は、ゲート１６における時間区間幅の影響を受ける。全ての深さで相対化して均等な条件下で後のスペクトル幅の算出が可能となるよう、ゲート１６の時間区間幅は一定値に固定する。
【００２７】
抽出波形スペクトル計算部８、スペクトル幅計算部９の詳細は、図３に示すものとなる。まず抽出波形スペクトル計算部８では、ゲート（ｉ）（ここでｉ＝１〜M）を用いて抽出された散乱エコー（図２の１７）に対し、高速フーリエ変換（ＦＦＴ）を適用する。そのため抽出波形スペクトル計算部８はM個のＦＦＴ演算部１８を有している。ＦＦＴ演算部１８のＦＦＴ出力をスペクトル幅計算部９にてパワースペクトルに変換し、パワースペクトルに対して重み積分処理を行って、スペクトル幅を算出する。さらに、深さ毎に算出されたスペクトル幅の値を、スペクトル幅プロファイル構成部１３おいて集合させ、合成すれば、超音波伝搬経路４上におけるスペクトル幅のプロファイルが得られ、これを図１におけるＮ本の走査線５において繰り返せば、空間分布が得られる。この空間分布を、最終的に、空間分布を呈示もしくは画像化させ、例えば表示部１４にそれを表示させることにより、超音波造影剤を含む生体組織を識別ないし可視化することが可能となる。
【００２８】
この一連の過程は図４のように数式表記される。すなわち、散乱エコーをｘ（ｔ）としたとき、ゲートを乗算し、深さ毎に算出されるパワースペクトルは、次のように表される。
【００２９】
【数１】
【００３０】
ここで、Ｐ（ω，ｔ）は局所的に抽出された波形のパワースペクトル、Ｔは局所抽出時間区間（ゲート時間区間幅）である。式（１）により、深さ毎にパワースペクトルが算出され、以下の重み積分処理を行って、スペクトル幅のプロファイルｂｗ（ｔ）が算出される。
【００３１】
【数２】
【００３２】
式（２）に示すプロファイルは、１走査線において得られるものであるため、画像を構成する全ての走査線５においてこの処理を行えば、スペクトル幅の空間分布が得られることになる。図３のパワースペクトル幅集合により構成されたプロファイル２０は、中間領域に超音波造影剤を含む生体組織が位置し、この部分のスペクトル幅が、周辺組織のスペクトル幅よりも小さくなることで、コントラスト分解能が向上する例示であり、この逆となることも有り得る。また、ここでのスペクトル幅の算出において、式（２）に示す重み積分処理は一例であり、別の方法としては、パワースペクトルのピーク値を求め、その半値となる周波数幅を求めてもよい。
【００３３】
次に、高速ルーチンの構成を説明する。高速ルーチンは、標準ルーチンを簡便化して、より高速にスペクトル幅の算出を行うためのルーチンである。高速ルーチンは、特定周波数成分乗算部１０、特定周波数パワー抽出部１１、スペクトル幅計算部（近似）１２からなる。
【００３４】
特定周波数成分乗算部１０の詳細は、図５に示すものとなる。ここでは、１走査線における受信エコー１５に対し、特定周波数成分２１を各々乗算する。ここで、特定周波数成分とは、受信エコー１５の中心角周波数ω０を起点とし、その前後、すなわち、ω０−ｎΔω〜ω０＋ｎΔωの範囲の（２ｎ＋１）個の周波数成分を指す。ｎは整数、Δωは角周波数の刻み幅である。これらの周波数成分を持つ正弦波及び余弦波をそれぞれ受信エコー１５に乗算することにより、受信エコー１５を周波数変調した複素数波形２２が得られる。
【００３５】
次に、図６の特定周波数パワー抽出部１１では、複素数波形２２を重み１の有限幅を持つ窓内で平滑化（すなわち局所的な積分を実行）し、得られた複素数成分のパワー（複素数の絶対値。実部と虚部の二乗和）を図６の特定周波数パワー抽出部１１に示すように周波数軸上にプロットする。すわわち、ω０−ｎΔωの周波数成分を持つ正弦波及び余弦波を乗算して得られた複素数波形から得られたパワーをω０−ｎΔωの位置にプロットし、これを（２ｎ＋１）回繰り返すと、パワースペクトルの離散点が得られる。その後、スペクトル幅計算部（近似）１２では、パワースペクトルの離散点を基に、曲線近似又は多項式近似によりパワースペクトル全体の分布形を推定する。パワースペクトルの分布形を推定した後は、標準ルーチンと同様、重み積分処理又はパワースペクトルピーク値の半値幅検出により、スペクトル幅を算出することができる。さらに、深さ毎に算出されたスペクトル幅の値を、スペクトル幅プロファイル構成部１３において集合させ、合成すれば、超音波伝搬経路４上におけるスペクトル幅のプロファイルが得られ、これを図１におけるＮ本の走査線５において繰り返せば、標準ルーチンと同様のスペクトル幅の空間分布が得られる。このスペクトル幅の空間分布を最終的に、空間分布を呈示もしくは画像化させ、表示部１４に表示させることにより、超音波造影剤を含む生体組織を識別ないし可視化することが可能となる。
【００３６】
この一連の過程は図７のように数式表記される。すなわち、散乱エコーをｘ（ｔ）としたとき、特定周波数ω０＋ｉΔω（ただし、ｉ＝−ｎ・・・０・・・ｎ，ｎ＞０）を持つ正弦波及び余弦波を乗算し、時間幅２Ｔを持つ窓内で単純平滑化（重み関数の係数が１の平滑化。局所的な積分に相当）すると、以下の式が得られる。
【００３７】
【数３】
【００３８】
この式は、式（１）を離散化したものに相当し、結局、式（３）は式（１）と等価な演算となる。式（３）を用いて計算されたパワースペクトルの離散点Ｐｉ（ｔ）を用いたパワースペクトル分布形の推定における曲線近似の一例として、J次多項式近似を行う場合、パワースペクトルの分布形は、
【００３９】
【数４】
【００４０】
のようにモデル化ないし近似されるので、得られた離散点データ（ωｉ，Ｐｉ（ｔ））を用いた最小二乗法により、係数行列{ａｊ}を求めれば、パワースペクトル分布形の推定が可能となり、スペクトル幅の算出が行われる。
【００４１】
ただし、パワースペクトル分布形の推定において、式（４）のモデルは一例であり、例えば、超音波パルス送信時のパワースペクトル分布関数モデルや、ガウス分布関数モデルなどパワースペクトル分布の近似関数モデルを利用することもできる。
【００４２】
上記では、式（３）により計算されたパワースペクトルの離散点Ｐｉ（ｔ）を用いたパワースペクトル分布形を、例えば式（４）の多項式近似モデルを適用して求め、その後、式（２）に示す重み積分処理を実行して、スペクトル幅の算出を精度よく行うための方法を述べたが、計算負荷の増大が懸念されるケースも生じ得る。このような場合への対処、すなわちスペクトル幅の算出において計算精度を多少犠牲にしても計算速度向上が望まれる場合、上記パワースペクトルの分布形計算処理を省略して、式（３）により計算されたパワースペクトルの離散点Ｐｉ（ｔ）に対して数値積分法を適用して、式（２）の重み積分を直接計算することができる。ここでの数値積分には、台形公式及びシンプソン公式及びガウス積分公式を適用できる。また、このときパワースペクトルの離散点Ｐｉ（ｔ）は時間領域で既に得られているため、上記数値積分法も時間領域において適用可能であり、すなわち、時間領域にて直接スペクトル幅を算出できる。
【００４３】
算出されたスペクトル幅は、そのまま空間分布として画像呈示できる他、制御手段において、該スペクトル幅から超音波減衰振動における時定数、または超音波造影剤における粒径分布の広がり分散又は標準偏差を算出し、呈示することもできる。
【００４４】
スペクトル幅からの減衰振動における時定数を得るために、まず、受信波形ｘ（ｔ）を次式で表す。
【００４５】
【数５】
【００４６】
ここで、Ａはｔ＝０での振幅、ｔは時間、Tは時定数、ω0は受信波の中心角周波数を表す。この受信波のパワースペクトルＰ（ω）は、次式のように表される。
【００４７】
【数６】
【００４８】
このパワースペクトルの半値幅は、２／Ｔとなり、すなわち、時定数の逆数の２倍に等しくなる。実際の受信波は、多くの超音波造影剤からの散乱波を全て足し合わせた総和になるため、実際のパワースペクトルも、個々の散乱波から得られるパワースペクトルの総和として得られる。しかしながら粒径分布が均一に近ければ、上記パワースペクトルの総和もまた、式（６）のように近似することができる。従って、パワースペクトルのスペクトル幅における半値幅の逆数を計算することにより、減衰振動の時定数もまた算出することができる。
【００４９】
一方、超音波造影剤における粒径分布の分散又は標準偏差は、気泡挙動の方程式により予め理論的に得たスペクトル幅対粒径分布分散または標準偏差の特性曲線の関数にスペクトル幅を代入することで、算出することができる。ここでは、気泡挙動の方程式に基づくスペクトル幅対粒径分布分散または標準偏差の特性曲線の関数について述べる。まず、超音波造影剤の共振周波数はその径に大きく依存し、粒径のばらつきが大きければその分スペクトル幅も広がり、また逆に粒径のばらつきが小さければその分スペクトル幅も狭くなることを利用し、算出されたスペクトル幅から粒径分布の広がり程度、すなわち分散又は標準偏差を評価しようとすることが、ここでの基本的な概念である。しかしながら、スペクトル幅として半値幅、粒径分布広がりとして標準偏差に着目するとき、散乱波パワースペクトルの半値幅と粒径分布の標準偏差との関係を表す関数は、超音波造影剤の平均粒径に依存し、さらには線形関係にない。粒径分布をガウス分布であると仮定すると、その標準偏差をσは一般に次のように記述できる。
【００５０】
【数７】
【００５１】
ここで、ｆは非線形関数を意味し、またＲｍとＦＷＨＭ（full width half maximum）はそれぞれ平均粒径と散乱波パワースペクトルの半値幅を表す。平均粒径は、超音波造影剤の粒径として現実的な３〜４μｍの定数と見なせるので、このときのｆの関数形が分かれば、散乱波パワースペクトルの半値幅ＦＷＨＭを式（７）に代入することにより、粒径分布の標準偏差σが得られることになる。ここで、ｆの関数形を予め求めておくための気泡挙動の方程式としては、よく知られたRayleigh-Plesset方程式を用いることができる。まずこの方程式に基づいて粒径対共振周波数を得る。これは、反比例的な関係式となる。次に、粒径分布がガウス分布であるとし、この分布を、上記の粒径対共振周波数の関係式に代入すると、共振周波数成分の集合としてのパワースペクトルが得られる。このパワースペクトルの半値幅は、ここでの粒径分布であるガウス分布の標準偏差と１対１の対応関係にある。従って、粒径分布におけるガウス分布の標準偏差を順次少しずつ変えたときに、それに対応する散乱波パワースペクトル半値幅の変化を求めてゆくと、曲線適合により、式（７）におけるｆの関数形を定めることができる。我々の検討によると、ｆは、σ＝（ＦＷＨＭ）ｑのような関数として近似できる（ただし、ｑは累乗における指数であり、平均粒径に応じて決まり、１未満の実数をとる。）。そして、このｑは、上記の曲線適合により決定することができるため、パワースペクトルの半値幅を代入することにより、粒径分布の標準偏差を算出することができる。さらに上記の累乗関数は、パワースペクトルの半値幅として現実的な０〜２ＭＨｚ付近までのレンジでは、急峻な勾配特性を持つため、僅かなパワースペクトル半値幅の違いが増幅され、粒径分布の標準偏差として得られることになる。これにより結果的に、周辺組織と、超音波造影剤を含む生体組織とのコントラスト分解能が向上することになる。
【００５２】
ここで、超音波造影剤を含む生体組織と周辺の生体組織における減衰振動特性の相違、あるいは周辺の生体組織の周波数特性と超音波造影剤の均一性に起因する周波数特性の相違のために、超音波造影剤を含む生体組織から発せられる散乱波形より得られるスペクトル幅と周辺の生体組織から発せられる散乱波形より得られるスペクトル幅はそれぞれ異なる値を有する。得られたスペクトル幅の生体組織内における空間分布自体をそのまま画像化すれば、超音波造影剤を含む生体組織と周辺の生体組織のコントラスト差は視認することが可能になる。しかしながら、このコントラスト差が顕著であり、表示画像のコントラスト分解能（例えば２５６階調）を超える場合は、得られたスペクトル幅を対数化してダイナミックレンジを圧縮することで、欠落なく全ての対象領域を画像化できる。また上記の減衰振動における時定数または粒径分布の標準偏差についても、これらのコントラスト差が顕著であり、表示画像のコントラスト分解能（例えば２５６階調）を超える場合は、得られた時定数ないし粒径分布の標準偏差を対数化してダイナミックレンジを圧縮することで、欠落なく全ての対象領域を画像化できる。
【００５３】
（実験例）
スペクトル幅算出によるコントラスト分解能向上を評価するため、上記で述べた装置を用い、造影実験を実施した。ここでは、超音波造影剤として中空マイクロカプセル及び市販の超音波造影剤であるレボビストを使用し、また参照波としてこれら造影剤を含むＰＰ（ポリプロピレン）容器の壁からの反射波を採取した。スペクトル幅算出の妥当性検証、特に、粒径が均一の超音波造影剤に対する検証には、本実験に用いた中空マイクロカプセルはより均一な粒径で製造可能である利点を利用することができる。Rayleigh-Plesset方程式に基づくと、気泡の共振周波数は、その径に大きく依存することになる。従って、造影剤の集合体において粒径のばらつきが大きければその分多様な共振周波数が存在することとなり、すなわちスペクトルにおけるメインローブのスペクトル幅は広がると考えられる。逆に、造影剤の集合体において粒径のばらつきが小さければ線スペクトルに近づき、スペクトル幅は狭まると考えられる。周囲組織が広帯域となる場合、このようなスペクトル幅の違いを局所的に算出できれば、コントラスト分解能の向上には有用である。中空マイクロカプセルの共振周波数に近い１ＭＨｚを用いたとき、パワースペクトル分布を図８（ａ）に、得られたスペクトル幅を図８（ｂ）に示す。なお、散乱波の受信も１ＭＨｚのプローブで行われたため、明確な高調波はここでは観察されていない。図８（ｂ）に示すように、中空マイクロカプセルは、参照となるＰＰ容器壁、すなわち疑似周囲組織に比べて有意に低下し、さらに粒径のばらつきが比較的大きな従来の超音波造影剤であるレボビストよりも半値幅は小さくなっている。これは、本手法により粒径分布の広がりが妥当に検出され得ることを示しており、結果的に、このようなスペクトル幅の差を呈示することにより、従来の超音波造影剤イメージング法よりも高いコントラスト分解能にもまた到達できる。
【図面の簡単な説明】
【００５４】
【図１】本発明の一実施形態に係る超音波造影装置における超音波受信信号の処理部分の全体構成を示す図である。
【図２】局所波形ゲート抽出部７の処理内容を説明した図である。
【図３】抽出波形スペクトル計算部８、スペクトル幅計算部９、スペクトル幅プロファイル構成部１３に至る処理内容を説明した図である。
【図４】図１における標準ルーチンの処理流れを数式を用いて説明した図である。
【図５】特定周波数成分乗算部１０の処理内容を説明した図である。
【図６】特定周波数パワー抽出部１１、スペクトル幅計算部（近似）１２、スペクトル幅プロファイル構成部１３に至る処理内容を説明した図である。
【図７】図１における高速ルーチンの処理流れを数式を用いて説明した図である。
【図８】本発明の装置によりスペクトル幅を算出した実験例である。
【符号の説明】
【００５５】
１ 超音波プローブ
２ 周辺組織
３ 超音波造影剤
４ 超音波送受信伝搬経路
５ 走査線
６ エコー受信部
７ （標準ルーチン）局所波形ゲート抽出部
８ （標準ルーチン）局所波形スペクトル幅計算部
９ （標準ルーチン）スペクトル幅計算部
１０ （高速ルーチン）特定周波数成分乗算部
１１ （高速ルーチン）特定周波数パワー抽出部
１２ （高速ルーチン）スペクトル幅計算部（近似値）
１３ スペクトル幅プロファイル構成部
１４ 表示部
１５ 受信エコー／走査線
１６ ゲート
１７ 抽出波形
１８ ＦＦＴ演算部
１９ パワースペクトル幅
２０ パワースペクトル幅集合により構成されたプロファイル
２１ 特定周波数成分
２２ 乗算された波形
２３ 離散的な特定周波数パワー値（丸）と近似されるパワースペクトル分布
２４ パワースペクトル分布を曲線近似により推定し、スペクトル幅を計算するルーチン
【特許請求の範囲】
【請求項１】
超音波造影剤を含む生体組織に対して超音波パルスを送受信する超音波プローブを有し、該超音波プローブにより超音波を送信した後、任意の組織部位から発せられる散乱波形を、該プローブを介して受信する散乱波形取得手段と、
該散乱波形取得手段が取得した散乱波形のパワースペクトルを算出するパワースペクトル算出手段と、
該パワースペクトル算出手段が算出したパワースペクトルの分散、標準偏差又は半値幅で定義されるスペクトル幅を特徴量として算出するスペクトル幅算出手段と、
該スペクトル幅算出手段が算出した特徴量を、対象とする生体組織内全てにおいて求めてその空間分布を得るように各手段を制御する制御手段と、
該制御手段の制御により、対象とする生体組織内全てにおいて求めた空間分布を出力する出力手段を備え、
超音波造影剤を含む生体組織を識別ないし可視化することを特徴とする超音波造影装置。
【請求項２】
該パワースペクトル算出手段は、任意の組織部位から発せられる散乱波形を、その組織部位の深さに対応する時間位置及び時間区間にのみ非零の値を持つステップ状ゲートと受信された散乱波形との乗算を行うことにより抽出する局所波形ゲート抽出部と、該局所波形ゲート抽出部により抽出された散乱波形に対してフーリエ変換を行ってパワースペクトルを算出する抽出波形スペクトル計算部を備え、
該スペクトル幅算出手段は、超音波造影剤を含む生体組織内における該スペクトル幅を、該抽出波形スペクトル計算部により算出されたパワースペクトルに対して重み積分処理ないし該パワースペクトルのピークの半値幅検出を行うことにより算出することを特徴とする請求項１に記載の超音波造影装置。
【請求項３】
該パワースペクトル算出手段は、散乱波形の中心周波数成分を含む、所定個数の特定周波数成分を持つ正弦波及び余弦波を各々散乱波形に乗算する特定周波数成分乗算部と、該特定周波数成分乗算部により得られた所定個数の乗算波形に対して平均化処理を行ってパワースペクトルを離散化した所定個数の点を求める特定周波数パワー抽出部を備え、
該スペクトル幅算出手段は、超音波造影剤を含む生体組織内における該スペクトル幅を、該特定周波数パワー抽出部が求めた所定個数のパワースペクトル離散点から、曲線近似ないし多項式近似を行うことにより推定されたパワースペクトルに対して重み積分処理ないし該パワースペクトルのピークの半値幅検出を行うことにより算出することを特徴とする請求項１に記載の超音波造影装置。
【請求項４】
該スペクトル幅算出手段は、超音波造影剤を含む生体組織内における該スペクトル幅を、該特定周波数パワー抽出部が求めた所定個数のパワースペクトル離散点を用いて台形公式、シンプソン公式及びガウス積分公式のうちのいずれかの数値積分法を適用し、時間領域において直接算出することを特徴とする請求項３に記載の超音波造影装置。
【請求項５】
該制御手段が、超音波造影剤を含む生体組織内における該スペクトル幅から、その逆数を計算することにより超音波減衰振動における時定数を算出し、その空間分布を呈示ないし画像化することにより、超音波造影剤を含む生体組織を識別ないし可視化することを特徴とする請求項１に記載の超音波造影装置。
【請求項６】
該制御手段が、粒径が均一な超音波造影剤を含む生体組織内における該スペクトル幅を、気泡挙動の方程式により予め理論的に得たスペクトル幅対粒径分布分散又は標準偏差の特性曲線の関数に代入することで粒径分布の分散又は標準偏差を算出し、その空間分布を呈示ないし画像化することにより、超音波造影剤を含む生体組織を識別ないし可視化することを特徴とする請求項１に記載の超音波造影装置。
【請求項７】
該スペクトル幅算出手段が、超音波造影剤を含む生体組織から発せられる散乱波形より得られるスペクトル幅と、周辺の生体組織から発せられる散乱波形より得られるスペクトル幅を各々検出し、
該制御手段が、それらの対数、あるいは時定数空間分布の対数、あるいは粒径分布の分散又は標準偏差の対数を算出し、その空間分布を呈示ないし画像化することにより、超音波造影剤を含む生体組織を識別ないし可視化することを特徴とする請求項１、５又は６に記載の超音波造影装置。
【請求項１】
超音波造影剤を含む生体組織に対して超音波パルスを送受信する超音波プローブを有し、該超音波プローブにより超音波を送信した後、任意の組織部位から発せられる散乱波形を、該プローブを介して受信する散乱波形取得手段と、
該散乱波形取得手段が取得した散乱波形のパワースペクトルを算出するパワースペクトル算出手段と、
該パワースペクトル算出手段が算出したパワースペクトルの分散、標準偏差又は半値幅で定義されるスペクトル幅を特徴量として算出するスペクトル幅算出手段と、
該スペクトル幅算出手段が算出した特徴量を、対象とする生体組織内全てにおいて求めてその空間分布を得るように各手段を制御する制御手段と、
該制御手段の制御により、対象とする生体組織内全てにおいて求めた空間分布を出力する出力手段を備え、
超音波造影剤を含む生体組織を識別ないし可視化することを特徴とする超音波造影装置。
【請求項２】
該パワースペクトル算出手段は、任意の組織部位から発せられる散乱波形を、その組織部位の深さに対応する時間位置及び時間区間にのみ非零の値を持つステップ状ゲートと受信された散乱波形との乗算を行うことにより抽出する局所波形ゲート抽出部と、該局所波形ゲート抽出部により抽出された散乱波形に対してフーリエ変換を行ってパワースペクトルを算出する抽出波形スペクトル計算部を備え、
該スペクトル幅算出手段は、超音波造影剤を含む生体組織内における該スペクトル幅を、該抽出波形スペクトル計算部により算出されたパワースペクトルに対して重み積分処理ないし該パワースペクトルのピークの半値幅検出を行うことにより算出することを特徴とする請求項１に記載の超音波造影装置。
【請求項３】
該パワースペクトル算出手段は、散乱波形の中心周波数成分を含む、所定個数の特定周波数成分を持つ正弦波及び余弦波を各々散乱波形に乗算する特定周波数成分乗算部と、該特定周波数成分乗算部により得られた所定個数の乗算波形に対して平均化処理を行ってパワースペクトルを離散化した所定個数の点を求める特定周波数パワー抽出部を備え、
該スペクトル幅算出手段は、超音波造影剤を含む生体組織内における該スペクトル幅を、該特定周波数パワー抽出部が求めた所定個数のパワースペクトル離散点から、曲線近似ないし多項式近似を行うことにより推定されたパワースペクトルに対して重み積分処理ないし該パワースペクトルのピークの半値幅検出を行うことにより算出することを特徴とする請求項１に記載の超音波造影装置。
【請求項４】
該スペクトル幅算出手段は、超音波造影剤を含む生体組織内における該スペクトル幅を、該特定周波数パワー抽出部が求めた所定個数のパワースペクトル離散点を用いて台形公式、シンプソン公式及びガウス積分公式のうちのいずれかの数値積分法を適用し、時間領域において直接算出することを特徴とする請求項３に記載の超音波造影装置。
【請求項５】
該制御手段が、超音波造影剤を含む生体組織内における該スペクトル幅から、その逆数を計算することにより超音波減衰振動における時定数を算出し、その空間分布を呈示ないし画像化することにより、超音波造影剤を含む生体組織を識別ないし可視化することを特徴とする請求項１に記載の超音波造影装置。
【請求項６】
該制御手段が、粒径が均一な超音波造影剤を含む生体組織内における該スペクトル幅を、気泡挙動の方程式により予め理論的に得たスペクトル幅対粒径分布分散又は標準偏差の特性曲線の関数に代入することで粒径分布の分散又は標準偏差を算出し、その空間分布を呈示ないし画像化することにより、超音波造影剤を含む生体組織を識別ないし可視化することを特徴とする請求項１に記載の超音波造影装置。
【請求項７】
該スペクトル幅算出手段が、超音波造影剤を含む生体組織から発せられる散乱波形より得られるスペクトル幅と、周辺の生体組織から発せられる散乱波形より得られるスペクトル幅を各々検出し、
該制御手段が、それらの対数、あるいは時定数空間分布の対数、あるいは粒径分布の分散又は標準偏差の対数を算出し、その空間分布を呈示ないし画像化することにより、超音波造影剤を含む生体組織を識別ないし可視化することを特徴とする請求項１、５又は６に記載の超音波造影装置。
【図２】
【図４】
【図５】
【図７】
【図１】
【図３】
【図６】
【図８】
【図４】
【図５】
【図７】
【図１】
【図３】
【図６】
【図８】
【公開番号】特開２０１０−６９２４３（Ｐ２０１０−６９２４３Ａ）
【公開日】平成２２年４月２日（２０１０．４．２）
【国際特許分類】
【出願番号】特願２００８−２４３０００（Ｐ２００８−２４３０００）
【出願日】平成２０年９月２２日（２００８．９．２２）
【出願人】（３０１０２１５３３）独立行政法人産業技術総合研究所 (6,529)
【Ｆターム（参考）】
【公開日】平成２２年４月２日（２０１０．４．２）
【国際特許分類】
【出願日】平成２０年９月２２日（２００８．９．２２）
【出願人】（３０１０２１５３３）独立行政法人産業技術総合研究所 (6,529)
【Ｆターム（参考）】
[ Back to top ]