WO2012137488A1

WO2012137488A1 - 温度推定方法、温度推定装置及びプログラム

Info

Publication number: WO2012137488A1
Application number: PCT/JP2012/002326
Authority: WO
Inventors: リラン; シュウフェンファン; コクセンチョン; アンツアントラン; 近藤　敏志
Original assignee: パナソニック株式会社
Priority date: 2011-04-07
Filing date: 2012-04-03
Publication date: 2012-10-11
Also published as: JPWO2012137488A1; US20130066584A1

Abstract

　超音波信号によって対象領域をスキャンすることにより得られるスキャン信号を受信し、スキャン信号に基づいて、超音波信号が前記対象領域を通過するのに要する時間が対象領域の温度に依存して変化する変化量であるエコーシフトを算出するエコーシフト算出ステップと、算出されたエコーシフトに基づいて、超音波信号が対象領域を通過する際の移動距離が対象領域の温度に依存して見かけ上変化する変化割合であるストレインを算出するストレイン算出ステップと、算出されたストレインに基づいて、ストレインの時間的な変化割合であるストレインレートを算出するストレインレート算出ステップと、予め定められたストレインとストレインレートと温度との関係を用いて、算出されたストレイン及びストレインレートに対応する対象領域の温度を推定する温度推定ステップと、を含む。

Description

温度推定方法、温度推定装置及びプログラム

　本発明は、超音波信号を用いることによって、例えば生体内部の所定部位の温度を推定する温度推定方法、温度推定装置及びプログラムに関する。

　近年、高集束超音波（ＨＩＦＵ：Ｈｉｇｈ　Ｉｎｔｅｎｓｉｔｙ　Ｆｏｃｕｓｅｄ　Ｕｌｔｒａｓｏｕｎｄ）を利用したＨＩＦＵ治療が行われている。ＨＩＦＵ治療とは、強力な超音波信号を患者の体表を介して腫瘍等の治療対象部位に集束させることにより、治療対象部位の温度を所定温度（例えば８０～９０℃程度）まで上昇させて治療対象部位を壊死させる治療法である。このＨＩＦＵ治療においては、確実な治療を行う観点から、治療対象部位が所定温度まで上昇したことを確認するために、治療対象部位の温度を推定する必要がある。また、安全性の確保の観点から、治療対象部位の周辺の部位の温度が必要以上に上昇していないことを確認するために、治療対象部位の周辺の温度を推定する必要がある。

　このように治療対象部位等の温度を推定する方法として、例えば、超音波信号を用いた温度推定方法（例えば、特許文献１～３並びに非特許文献１及び２参照）が行われている。

米国特許第４４５２０８１号明細書米国特許第７２１１０４４号明細書米国特許出願公開第２００７／０１０６１５７号明細書

Ｍｉｌｌｅｒ，Ｎ．Ｒ．，Ｊ．Ｃ．Ｂａｍｂｅｒ，ａｎｄ　Ｐ．Ｍ．　Ｍｅａｎｅｙ，Ｆｕｎｄａｍｅｎｔａｌ　ｌｉｍｉｔａｔｉｏｎｓ　ｏｆ　ｎｏｎｉｎｖａｓｉｖｅ　ｔｅｍｐｅｒａｔｕｒｅ　ｉｍａｇｉｎｇ　ｂｙ　ｍｅａｎｓ　ｏｆ　ｕｌｔｒａｓｏｕｎｄ　ｅｃｈｏ　ｓｔｒａｉｎ　ｅｓｔｉｍａｔｉｏｎ．Ｕｌｔｒａｓｏｕｎｄ　ｉｎ　ｍｅｄｉｃｉｎｅ　＆　ｂｉｏｌｏｇｙ，２００２．２８（１０）：ｐ．１３１９．Ｖａｌｖａｎｏ，　Ｊ．Ｗ．，　Ｂｉｏｈｅａｔ　ｔｒａｎｓｆｅｒ，ｉｎ　Ｅｎｃｙｃｌｏｐｅｄｉａ　ｏｆ　Ｍｅｄｉｃａｌ　Ｄｅｖｉｃｅｓ　ａｎｄ　Ｉｎｓｔｒｕｍｅｎｔａｔｉｏｎ．２００５，Ｗｉｌｅｙ．

　しかしながら、従来の超音波信号を用いた温度推定方法では、精度良く温度を推定することができず、温度を推定可能な温度範囲が狭いという問題がある。

　本発明は、上記従来の課題を解決するものであり、その目的は、高精度に温度を推定することができるとともに、広い温度範囲で温度を推定することができる温度推定方法、温度推定装置及びプログラムを提供することである。

　上記目的を達成するために、本発明の一態様に係る温度推定方法は、超音波信号を用いて対象領域の温度を推定する温度推定方法であって、超音波信号によって前記対象領域をスキャンすることにより得られるスキャン信号を受信し、前記スキャン信号に基づいて、超音波信号が前記対象領域を通過するのに要する時間が前記対象領域の温度に依存して変化する変化量であるエコーシフトを算出するエコーシフト算出ステップと、算出されたエコーシフトに基づいて、超音波信号が前記対象領域を通過する際の移動距離が前記対象領域の温度に依存して見かけ上変化する変化割合であるストレインを算出するストレイン算出ステップと、算出されたストレインに基づいて、ストレインの時間的な変化割合であるストレインレートを算出するストレインレート算出ステップと、予め定められたストレインとストレインレートと温度との関係を用いて、算出されたストレイン及びストレインレートに対応する前記対象領域の温度を推定する温度推定ステップと、を含む。

　なお、本発明は、方法として実現できるだけでなく、方法を構成するステップを処理手段とする装置として実現したり、それらステップをコンピュータに実行させるプログラムとして実現したり、そのプログラムを記録したコンピュータ読み取り可能なＣＤ－ＲＯＭ等の記録媒体として実現したり、そのプログラムを示す情報、データ又は信号として実現することもできる。それらプログラム、情報、データ及び信号は、インターネット等の通信ネットワークを介して配信してもよい。

　本発明の温度推定方法では、予め定められたストレインとストレインレートと温度との関係を用いて、算出されたストレイン及びストレインレートに対応する対象領域の温度を推定することにより、高精度に温度を推定することができるとともに、広い温度範囲で温度を推定することができる。

図１は、本発明の実施の形態１に係る温度推定装置の構成を示すブロック図である。図２は、スキャン信号を説明するための図である。図３は、図１の信号処理部の構成を具体的に示すブロック図である。図４は、前処理部の構成を具体的に示すブロック図である。図５は、スキャン信号に生じた破損フレームを示す図である。図６は、エコーシフト算出部の構成を具体的に示すブロック図である。図７は、エコーシフト計算ステップを説明するための図である。図８は、異常値修正ステップを説明するための図である。図９は、ストレイン算出部の構成を具体的に示すブロック図である。図１０は、ストレインレート算出部の構成を具体的に示すブロック図である。図１１は、温度推定部の構成を具体的に示すブロック図である。図１２は、本発明の実施の形態１に係る温度推定方法による温度の推定の流れを示すフローチャートである。図１３は、温度推定ステップを説明するための図である。図１４は、本発明の温度推定方法と従来の温度推定方法とを比較した実験結果を示す図である。図１５は、本発明の温度推定方法と従来の温度推定方法とを比較した実験結果を示す図である。図１６は、本発明の実施の形態２に係る温度推定装置の構成を示すブロック図である。図１７は、本発明の実施の形態３に係る温度推定装置の構成を示すブロック図である。図１８は、従来の温度推定方法における、モデル化された温度とストレインとの関係を示す図である。

　（本発明の基礎となった知見）
　まず、本発明の実施の形態を説明する前に、本発明者が見出した、従来の温度推定方法において生じる不具合について説明する。

　（従来の温度推定方法１）
　図１８は、従来の温度推定方法における、モデル化された温度とストレインとの関係を示す図である。従来の温度推定方法１では、超音波信号の速度が温度に依存して変化することを利用して、治療対象部位等の温度を推定する。例えば、治療対象部位に加熱用の超音波信号を集束させてＨＩＦＵ治療を行う前後において、治療対象部位の温度が変化することに伴って、温度推定用の超音波信号の速度が変化する。このように超音波信号の速度が変化することにより、超音波信号が治療対象部位を通過する際の移動距離が見かけ上変化する。この移動距離の変化割合をストレイン（ｓｔｒａｉｎ）という。

　従来の温度推定方法１では、治療前における超音波信号と治療中における超音波信号との差に基づいてストレインを計算により求めた後に、温度とストレインとの関係を示すモデル式に基づいて温度を推定する。図１８中の実線のグラフ４０１で示すように、このモデル式は、２次関数で近似されている。

　しかしながら、上述したモデル式では、図１８に示すように、温度変化に対するストレイン変化が小さい領域であるフラットレンジ（ｆｌａｔ　ｒａｎｇｅ）４０２が存在し、このフラットレンジ４０２においては、ノイズに対する影響が大きくなる。そのため、従来の温度推定方法１では、精度良く温度を推定することができないという問題がある。

　（従来の温度推定方法２）
　従来の温度推定方法２では、上述と同様に、例えばＨＩＦＵ治療前における超音波信号とＨＩＦＵ治療中における超音波信号との差に基づいてストレインを計算により求めた後に、温度とストレインとの関係を示すモデル式に基づいて温度を推定する。図１８中の破線のグラフ４０３で示すように、このモデル式は、１次関数で近似されている。

　しかしながら、上述したモデル式では、比較的低い温度領域（３０～５０℃程度）でしか温度を推定することができない。ＨＩＦＵ治療では、比較的高い温度領域（７０～９０℃程度）で温度を推定する必要があるため、従来の温度推定方法２をＨＩＦＵ治療における温度の推定に用いることが難しいという問題がある。

　（従来の温度推定方法３）
　従来の温度推定方法３では、モデル式として生体伝熱方程式（Ｂｉｏ－Ｈｅａｔ　Ｔｒａｎｓｆｅｒ　Ｅｑｕａｔｉｏｎ）が用いられる。まず、モデル式のパラメータをキャリブレーション時に求め、さらに、温度とストレインとの関係を示す関係式を求めておく。例えば、ＨＩＦＵ治療時には、モデル式に基づいて治療対象部位の温度を推定する。その後、推定した温度からストレインを上記関係式に基づいて計算により求め、計算により求めたストレインとＲＦ（Ｒａｄｉｏ　Ｆｒｅｑｕｅｎｃｙ）信号から得られたストレインとの誤差が小さくなるように、モデル式のパラメータを逐次修正する。

　しかしながら、従来の温度推定方法３では、モデル式のパラメータを求めるのが困難であるため、実用化が難しいという問題がある。

　上述のような問題を解決するために、本発明の一態様に係る温度推定方法は、超音波信号を用いて対象領域の温度を推定する温度推定方法であって、超音波信号によって前記対象領域をスキャンすることにより得られるスキャン信号を受信し、前記スキャン信号に基づいて、超音波信号が前記対象領域を通過するのに要する時間が前記対象領域の温度に依存して変化する変化量であるエコーシフトを算出するエコーシフト算出ステップと、算出されたエコーシフトに基づいて、超音波信号が前記対象領域を通過する際の移動距離が前記対象領域の温度に依存して見かけ上変化する変化割合であるストレインを算出するストレイン算出ステップと、算出されたストレインに基づいて、ストレインの時間的な変化割合であるストレインレートを算出するストレインレート算出ステップと、予め定められたストレインとストレインレートと温度との関係を用いて、算出されたストレイン及びストレインレートに対応する前記対象領域の温度を推定する温度推定ステップと、を含む。

　本態様によれば、高精度に温度を推定することができるとともに、広い温度範囲で温度を推定することができる。

　例えば、本発明の一態様に係る温度推定方法において、前記温度推定ステップにおいては、ストレインとストレインレートと温度との関係を示す線形モデル式を用いて、算出されたストレイン及びストレインレートに対応する前記対象領域の温度を推定するように構成してもよい。

　本態様によれば、ストレインとストレインレートと温度との関係を示す線形モデル式を用いて、温度を予測することができる。

　例えば、本発明の一態様に係る温度推定方法において、さらに、外乱によって前記スキャン信号に生じた破損フレームを除去する前処理ステップを含み、前記エコーシフト算出ステップでは、前記前処理ステップにおいて破損フレームが除去されたスキャン信号に基づいて、エコーシフトを算出するように構成してもよい。

　本態様によれば、エコーシフトを算出する前の前処理として、スキャン信号に生じた破損フレームを除去することができる。

　例えば、本発明の一態様に係る温度推定方法において、前記エコーシフト算出ステップは、受信された前記スキャン信号に基づいて生のエコーシフトを計算するエコーシフト計算ステップと、計算された前記生のエコーシフトにおいて温度変化に起因しない異常値を修正することにより、温度変化に起因するエコーシフトを得る異常値修正ステップと、前記温度変化に起因するエコーシフトに含まれるノイズをノイズフィルタにより減衰させるノイズ減衰ステップと、を含むように構成してもよい。

　本態様によれば、計算された生のエコーシフトにおいて温度変化に起因しない異常値を修正することにより、温度変化に起因するエコーシフトを得ることができる。

　例えば、本発明の一態様に係る温度推定方法において、前記異常値修正ステップは、前記エコーシフト計算ステップにおいて計算された前記生のエコーシフトに基づいて、予測エコーシフトを計算するためのエコーシフト予測モデルのパラメータを調整するパラメータ調整ステップと、最適化された前記エコーシフト予測モデルに基づいて、前記予測エコーシフトを計算する予測エコーシフト計算ステップと、前記予測エコーシフト及び前記対象領域をスキャンすることにより得られるＲＦ信号に基づいて、閾値を設定する閾値設定ステップと、前記生のエコーシフトの複数のデータポイントのうち、前記閾値を超えた前記異常値であるデータポイントを、前記予測エコーシフトの対応するデータポイントに基づいて修正するデータポイント修正ステップと、を含むように構成してもよい。

　本態様によれば、異常値の検出を行うことができる。

　例えば、本発明の一態様に係る温度推定方法において、前記エコーシフト予測モデルは、誤差関数又は相補誤差関数であるように構成してもよい。

　本態様によれば、エコーシフト予測モデルを誤差関数又は相補誤差関数で構成することができる。

　例えば、本発明の一態様に係る温度推定方法において、前記閾値設定ステップは、前記閾値としての上側の閾値と下側の閾値との間の間隔を設定する間隔設定ステップと、前記ＲＦ信号の強度を重みとして用いることにより、前記上側の閾値及び前記下側の閾値を調節する閾値調節ステップと、を含むように構成してもよい。

　本態様によれば、異常値を検出するための閾値を設定することができる。

　例えば、本発明の一態様に係る温度推定方法において、前記ストレインレート算出ステップにおいては、深さ及び時間の積と深さとに基づいてエコーシフトを予測するためのモデル式を用いて、前記スキャン信号に基づいて算出されたエコーシフトと前記モデル式に基づいて予測されたエコーシフトとの誤差が最小となるように、ストレインレートを算出するように構成してもよい。

　本態様によれば、ストレインレートを算出することができる。

　例えば、本発明の一態様に係る温度推定方法において、前記ストレインレート算出ステップにおいて用いられる前記モデル式は、線形モデル式であるように構成してもよい。

　本態様によれば、ストレインレート算出ステップにおいて用いられるモデル式を線形モデル式で構成することができる。

　例えば、本発明の一態様に係る温度推定方法において、前記温度推定ステップにおいて用いられる前記線形モデル式のアルゴリズムパラメータは、加熱すべき前記対象領域に対して熱電対の位置を位置合わせする位置合わせステップと、前記熱電対によって前記対象領域の温度を測定する温度測定ステップと、前記対象領域をスキャンすることにより得られる前記スキャン信号に基づいて、前記対象領域の推定温度を導出する推定温度導出ステップと、前記温度測定ステップにおいて測定された温度と前記推定温度導出ステップにおいて導出された推定温度との誤差が最小となるように、前記アルゴリズムパラメータを同定するアルゴリズムパラメータ同定ステップと、を実行することにより同定されるように構成してもよい。

　本態様によれば、温度推定ステップにおいて用いられる線形モデル式のアルゴリズムパラメータを同定することができる。

　例えば、本発明の一態様に係る温度推定方法において、前記位置合わせステップは、Ｂモード画像によって、加熱すべき前記対象領域に対して前記熱電対の位置を位置合わせするステップと、前記対象領域の温度を変化させるために、超音波信号を出力する加熱源を作動させるステップと、前記加熱源を空間内で変位させながら、前記熱電対によって前記対象領域の温度を測定するステップと、前記熱電対により測定された温度が最小となったときに前記加熱源の変位を停止させるステップと、を含むように構成してもよい。

　本態様によれば、加熱すべき対象領域に対して熱電対の位置を位置合わせすることができる。

　例えば、本発明の一態様に係る温度推定方法において、さらに、温度の空間的な連続性を考慮した目的関数に基づいて、前記温度推定ステップにより推定された温度を補正する温度補正ステップを含むように構成してもよい。

　本態様によれば、温度の空間的な連続性を考慮して、温度推定ステップにより推定された温度を補正することができる。

　例えば、本発明の一態様に係る温度推定方法において、さらに、加熱用の超音波信号により加熱された熱源位置と前記対象領域との間の距離を考慮した目的関数に基づいて、前記温度推定ステップにより推定された温度を補正する温度補正ステップを含むように構成してもよい。

　本態様によれば、加熱用の超音波信号（例えばＨＩＦＵ）により加熱された熱源位置と対象領域との間の距離を考慮して、温度推定ステップにより推定された温度を補正することができる。

　本発明の一態様に係る温度推定装置は、超音波信号を用いて対象領域の温度を推定する温度推定装置であって、超音波信号によって前記対象領域をスキャンすることにより得られるスキャン信号を受信し、前記スキャン信号に基づいて、超音波信号が前記対象領域を通過するのに要する時間が前記対象領域の温度に依存して変化する変化量であるエコーシフトを算出するエコーシフト算出部と、算出されたエコーシフトに基づいて、超音波信号が前記対象領域を通過する際の移動距離が前記対象領域の温度に依存して見かけ上変化する変化割合であるストレインを算出するストレイン算出部と、算出されたストレインに基づいて、ストレインの時間的な変化割合であるストレインレートを算出するストレインレート算出部と、予め定められたストレインとストレインレートと温度との関係を用いて、算出されたストレイン及びストレインレートに対応する前記対象領域の温度を推定する温度推定部と、を備える。

　本発明の一態様に係るプログラムは、超音波信号を用いて対象領域の温度を推定するためのプログラムであって、超音波信号によって前記対象領域をスキャンすることにより得られるスキャン信号を受信し、前記スキャン信号に基づいて、超音波信号が前記対象領域を通過するのに要する時間が前記対象領域の温度に依存して変化する変化量であるエコーシフトを算出するエコーシフト算出ステップと、算出されたエコーシフトに基づいて、超音波信号が前記対象領域を通過する際の移動距離が前記対象領域の温度に依存して見かけ上変化する変化割合であるストレインを算出するストレイン算出ステップと、算出されたストレインに基づいて、ストレインの時間的な変化割合であるストレインレートを算出するストレインレート算出ステップと、予め定められたストレインとストレインレートと温度との関係を用いて、算出されたストレイン及びストレインレートに対応する前記対象領域の温度を推定する温度推定ステップと、をコンピュータに実行させる。

　なお、これらの全般的又は具体的な態様は、システム、方法、集積回路、コンピュータプログラム又は記録媒体で実現されてもよく、システム、方法、集積回路、コンピュータプログラム又は記録媒体の任意な組み合わせで実現されてもよい。

　（実施の形態）
　以下、本発明の実施の形態について、図面を用いて詳細に説明する。なお、以下で説明する実施の形態は、いずれも本発明の一具体例を示すものである。以下の実施の形態で示される数値、形状、材料、構成要素、構成要素の配置位置、ステップ及びステップの順序等は、一例であり、本発明を限定する主旨ではない。また、以下の実施の形態における構成要素のうち、最上位概念を示す独立請求項に記載されていない構成要素については、任意の構成要素として説明される。

　（実施の形態１）
　（温度推定装置の全体構成）
　図１は、本発明の実施の形態１に係る温度推定装置の構成を示すブロック図である。図１に示すように、本実施の形態の温度推定装置１０は、本体部１０１及び超音波振動子１０２を備えている。

　超音波振動子１０２は、生体５０の内部の対象領域５１（例えば、ＨＩＦＵ治療における治療対象部位）に対して超音波信号を送信するとともに、対象領域５１をスキャンして対象領域５１で反射された超音波信号を受信する。

　本体部１０１は、送信部１０３、受信部１０４、信号処理部１０５、表示部１０６及び制御部１０７を備えている。送信部１０３は、制御部１０７の制御に従って、超音波振動子１０２に対して電気信号（送信信号）を供給することにより、超音波振動子１０２に超音波信号を発生させる。受信部１０４は、制御部１０７の制御に従って、超音波振動子１０２からの電気信号（受信信号）を受信する。

　信号処理部１０５は、受信部１０４から供給された電気信号（スキャン信号）を処理することにより、対象領域５１の温度を推定する。信号処理部１０５による温度推定方法については後述する。表示部１０６は、制御部１０７の制御に従って、信号処理部１０５により推定された温度を表示する。表示部１０６は、例えば液晶ディスプレイ等で構成される。

　制御部１０７は、送信部１０３、受信部１０４、信号処理部１０５及び表示部１０６をそれぞれ制御することにより、温度推定装置１０の全体制御を行う。制御部１０７は、例えばマイクロプロセッサ等で構成される。

　（温度推定のモデル）
　図２は、スキャン信号を説明するための図である。図２に示すように、スキャン信号は、３つの次元、即ち、深さ方向（ｄｅｐｔｈ　ｄｉｒｅｃｔｉｏｎ）（ｄ）、ライン方向（ｌｉｎｅ　ｄｉｒｅｃｔｉｏｎ）（ｌ）及びフレーム方向（ｆｒａｍｅ　ｄｉｒｅｃｔｉｏｎ）（ｎ）を有している。深さ方向は、生体５０の体表から生体５０の内部に向かう方向であり、ライン方向は、深さ方向に対して略垂直な方向である。フレーム方向は、フレームが切り替わる時間方向である。一般に、深さ方向における時間をファストタイム（ｆａｓｔ　ｔｉｍｅ）、フレーム方向における時間をスロータイム（ｓｌｏｗ　ｔｉｍｅ）という。

　対象領域５１の温度が変化することに伴って、超音波信号の速度、即ち、音速が変化する。このように超音波信号の速度が変化することにより、超音波信号が対象領域５１を通過するのに要する時間が変化する。このように、超音波信号が対象領域５１を通過するのに要する時間が対象領域５１の温度に依存して変化する変化量をエコーシフト（ｅｃｈｏ　ｓｈｉｆｔ）という。例えば、図２に示すように、フレーム１におけるＲＦ信号（スキャン信号）と、フレーム１の次のフレームであるフレーム２におけるＲＦ信号との間に生じる時間的なずれがエコーシフトである。

　図１に示すように、超音波信号は、生体５０の体表からその内部に向けて深さ方向に沿って走行する。生体５０の体表から対象領域５１の第１の部位５１ａまでの距離（深さ）をｄ_１、生体５０の体表から対象領域５１の第２の部位５１ｂまでの距離をｄ_２、対象領域５１の深さ方向における大きさをｄ（＝ｄ_２－ｄ_１）、対象領域５１の温度がＴ_０℃である場合における音速をｃ_０とする。このとき、超音波信号が対象領域５１を通過するのに要する時間（ファストタイム）ｔ_ｆ０は、

である。なお、式１以降において、ｔ_ｆはファストタイムを表し、ｎはスロータイムを表す。対象領域５１の温度がＴ_０℃からＴ_１℃に変化した際には、超音波信号が対象領域５１を通過するのに要する時間（ファストタイム）ｔ_ｆ１は、

である。なお、式２において、ｃ_１は、対象領域５１の温度がＴ_１℃である場合における音速である。

　従って、対象領域５１の温度がＴ_０℃からＴ_１℃に変化した際のエコーシフトψは、

で表される。

　対象領域５１の温度がＴ_０℃からＴ_１℃に変化した際における、超音波信号が対象領域５１を通過する際の移動距離の見かけ上の変化量Δｄｉｓは、

である。式４において、ｃ_{ｓｙｓｔｅｍ}は、超音波画像システムにおいて音速として用いられる定数であり、一般に、ｃ_{ｓｙｓｔｅｍ}＝１５４０ｍ／ｓｅｃである。

　ここで、超音波信号が対象領域５１を通過する際の移動距離が対象領域５１の温度に依存して見かけ上変化する変化割合であるストレインεは、

で表される。式５において、ｃ_{ｓｙｓｔｅｍ}は定数であるので、ストレインεは、

で表される。

　ここで、音速ｃ_１及び温度Ｔは、下記の二次多項式、

を満たすことが知られている。ａ_１，ａ_２及びａ_３はそれぞれ、一定の温度係数である。一般に、対象領域５１の温度に依存して音速が変化する変化量は、音速ｃ_０よりも十分に小さいことが知られている。即ち、

である。これにより、式５と式７とを組み合わせることにより、

が得られる。式９において、ａ_４，ａ_５及びａ_６はそれぞれ、

であり、定数である。なお、式９のモデル式は、上述した従来の温度推定方法１で用いられるモデル式と同様である。

　本実施の形態の温度推定方法では、式９のモデル式とＰｅｎｎｅｓの生体伝熱方程式とを組み合わせることにより導出される式１６のモデル式（後述する）を用いて、対象領域５１の温度を推定する。Ｐｅｎｎｅｓの生体伝熱方程式は、生体の内部の熱伝導を表す微分方程式であり、

で表される。式１１において、Ｑ_ｍｅｔ：代謝発熱、ΔＴ：環境上の温度上昇、ρ：組織密度、υ：組織の比熱、ρ_ｂ：血液密度、υ_ｂ：血液の比熱、ω：血液のかん流比率、κ：組織の熱伝導率である。

　ここで、１）フレーム方向における対象領域５１の温度変化のみが考慮され、２）組織及び血液の各々の特性は時間とともに変化せず、３）組織の温度分布は測定前において均一であり、４）組織の代謝発熱は時間とともに変化しない、と仮定する。これにより、ΔＴ＝Ｔ－Ｔ_０及び温度の基準を示すＴ_０はそれぞれ定数である。また、Ｑ_ｍｅｔ，ρ，υ，ρ_ｂ，υ_ｂ，ω及びκはそれぞれ定数である。さらに、∇^２Ｔは、時間変化しないので定数である。従って、式１１は、

と書き換えられる。式１２において、ａ_７，ａ_８はそれぞれ、

を満たす定数である。

　次に、ストレイン及び温度の関数はいずれもスロータイムに沿って連続的であるので、上記式９の両辺について偏導関数を計算することにより、ストレインレート（ｓｔｒａｉｎ－ｒａｔｅ）が得られる。

　式１３において、ξはストレインレートである。ストレインレートξは、ストレインεの時間的な変化割合である。換言すると、ストレインレートξは、深さ方向の深さｄ及び時間（スロータイム）ｎに沿ったエコーシフトψの混合偏導関数、

である。

　Ｐｅｎｎｅｓの生体伝熱方程式により推定された温度が実際の温度と等しいと仮定すると、式１２と式１３とを組み合わせることにより、次式１５が得られる。

　式９と式１５とを組み合わせて二次の項（Ｔ^２）を消去すると、

が得られる。式１６において、ａ_１２，ａ_１３及びａ_１４はそれぞれ、

を満たす定数である。

　式１６のモデル式は、予め定められたストレインεとストレインレートξと温度Ｔとの関係を示す線形モデル式である。本実施の形態の温度推定方法では、このモデル式を用いて、後述するようにして対象領域５１の温度を推定する。

　（信号処理部の構成）
　図３は、図１の信号処理部の構成を具体的に示すブロック図である。図３に示すように、信号処理部１０５は、前処理部２０１、エコーシフト算出部２０２、ストレイン算出部２０３、ストレインレート算出部２０４及び温度推定部２０５を有している。

　図４は、前処理部の構成を具体的に示すブロック図である。図４に示すように、前処理部２０１は、破損フレーム除去部３０１及びノイズフィルタ３０２を有している。破損フレーム除去部３０１は、外乱（例えば、環境ノイズ、外部の処理信号及びＨＩＦＵ治療で使用される高集束超音波等）によってスキャン信号に生じた破損フレームをスキャン信号から除去する（前処理ステップ）。なお、図３以降における「スキャン信号（ｄ，ｌ，ｎ）」等では、ｄは深さ方向の次元、ｌはライン方向の次元、ｎはフレーム方向の次元を表している。

　図５は、スキャン信号に生じた破損フレームを示す図である。図５に示すように、破損フレームにおいては、スキャン信号（ＲＦ信号）の強度の振幅が異常に大きくなる。

　破損フレーム除去部３０１は、スキャン信号の強度のグラフにおける上側の包絡線及び下側の包絡線をそれぞれ検出した後に、上側の包絡線と下側の包絡線との差分をフレーム毎に計算により求める。その後、破損フレーム除去部３０１は、２つの包絡線の差分が所定の閾値以上であるフレームを上述した破損フレームと判定し、スキャン信号から破損フレームを除去する。なお、２つの包絡線の差分が上記所定の閾値以上であるフレームは、例えばＨＩＦＵ治療において高集束超音波が照射される加熱段階に対応する。２つの包絡線の差分が上記所定の閾値よりも小さいフレームは、例えばＨＩＦＵ治療において高集束超音波が照射された後の冷却段階に対応する。これにより、ＨＩＦＵ治療の冷却段階におけるスキャン信号のフレームは、後述する温度推定部２０５における温度の推定に使用され、ＨＩＦＵ治療の加熱段階におけるスキャン信号のフレーム（即ち、破損フレーム）は、温度推定部２０５における温度の推定には使用されない。

　破損フレーム除去部３０１から出力されたスキャン信号は、ノイズフィルタ３０２においてノイズが減衰された後に、ノイズフィルタ３０２からフィルタリング信号（即ち、フィルタリングされたスキャン信号）として出力される。

　なお、前処理部２０１においては、体の動きの補償及び血液の流れの補償を考慮した処理を行うこともできる。体の動きは、動き推定技術における超音波信号（又は画像）を用いることにより、或いは、加速度センサ等を用いることにより検出される。血液の流れは、超音波ドプラ信号処理を用いることにより検出される。

　図６は、エコーシフト算出部の構成を具体的に示すブロック図である。図６に示すように、エコーシフト算出部２０２は、エコーシフト計算部３０３、異常値修正部３０４及びエコーシフトノイズフィルタ３０５を有している。エコーシフト算出部２０２は、後述するようにして、フィルタリング信号（スキャン信号）に基づいてエコーシフトを算出する（エコーシフト算出ステップ）。

　エコーシフト計算部３０３は、前処理部２０１から出力されたフィルタリング信号に基づいて、生のエコーシフトを計算する（エコーシフト計算ステップ）。生のエコーシフトを計算する方法として、例えば、自己相関方法及び相互相関方法がある。

　自己相関方法では、例えばＳＤｏｐｐ推定方法が用いられる。このＳＤｏｐｐ推定方法では、空間的に接触された多数のＩＱデータを用いることにより、エコーシフトが計算される。ＳＤｏｐｐ推定方法により計算されるエコーシフトは、

で表される。式１７において、ψ（ｄ，ｎ）は、深さｄ^ｔｈ及び時間（スロータイム）ｎ^ｔｈでのエコーシフトである。式１７において、ｋは、エコーシフトが計算された深さ方向におけるサンプルの数であり、ｙは、エコーシフトが計算されたフレームの数である。式１７において、Ｉ（ｄ，ｎ）は、深さｄ^ｔｈ及び時間（スロータイム）ｎ^ｔｈから選択されたＩＱセグメントである。ＳＤｏｐｐ推定方法において、（ｋ^＊ｙ）の窓の大きさは、一定値としてもよく、或いは、スキャン信号の特性の変化（例えば、振幅又はエネルギーの変化）に応じて可変することもできる。

　相互相関方法は、例えば、数学的に次式１８によって表される。

　式１８において、Ｓ_ｎは、フレームｎにおけるＲＦセグメントであり、Ｓ_ｎ+１は、フレームｎと隣接するフレームｎ+１におけるＲＦセグメントである。式１８において、Ｓ_ｎ ^ｍｅａｎ及びＳ_ｎ ^ＲＭＳはそれぞれ、セグメントの平均及び二乗平均平方根である。式１８において、ｋは、セグメントの窓の長さであり、ｑは、フレームのＲＦラインにおける探索範囲であり、γは、相互相関係数である。βは、γ（β）が最大値であるときの２つのセグメント間におけるエコーシフトと等しい。

　図７は、エコーシフト計算ステップを説明するための図である。図７に示すように、エコーシフト計算部３０３は、まず、フレームｎ+１における探索範囲ｑの上端に窓（Ｓ_ｎ+１）を変位させた状態で、セグメントＳ_ｎに窓（Ｓ_ｎ+１）を乗じる。この計算のステップは、データポイント毎に行われ、合計値は変数で蓄えられる。その後、エコーシフト計算部３０３は、窓（Ｓ_ｎ+１）が探索範囲ｑの下端に到達するまで、サンプル毎に窓（Ｓ_ｎ+１）を変位させながら、上述した計算を繰り返し行う。エコーシフト計算部３０３は、γ（β）が最大値に到達したときのβを、生のエコーシフトとして決定する。

　異常値修正部３０４は、エコーシフト計算部３０３により計算された生のエコーシフトにおいて、温度変化に起因しない異常値を修正する（異常値修正ステップ）。生のエコーシフトには、温度変化に起因するエコーシフト及び温度変化に起因しない（例えば、振動によって生じた）異常値が含まれる。この異常値は、温度推定の正確性を著しく低下させる。この異常値を修正することにより、温度変化に起因するエコーシフトを得ることができる。異常値修正部３０４は、例えば、エコーシフト予測モデルに基づいて異常値を除去する。このエコーシフト予測モデルは、誤差関数又は相補誤差関数で構成される。

　エコーシフト予測モデルとして誤差関数を用いた場合には、上述した異常値修正ステップは、次のようにして行われる。図８は、異常値修正ステップを説明するための図である。まず、エコーシフト計算ステップにおいて計算された生のエコーシフトに基づいて、予測エコーシフトを計算するためのエコーシフト予測モデルのパラメータを調整する（パラメータ調整ステップ）。次に、最適化されたエコーシフト予測モデルに基づいて、予測エコーシフトを計算する（予測エコーシフト計算ステップ）。その後、予測エコーシフト及び対象領域５１をスキャンすることにより得られるＲＦ信号に基づいて、上側の閾値及び下側の閾値を設定する（閾値設定ステップ）。この閾値設定ステップにおいては、上側の閾値と下側の閾値との間の間隔を予備的実験等に基づいて設定し（間隔設定ステップ）、ＲＦ信号の強度を重みとして用いることにより、上側の閾値及び下側の閾値を調節する（閾値調節ステップ）。その後、生のエコーシフトの複数のデータポイントのうち、上側の閾値及び下側の閾値を超えた異常値であるデータポイントを、予測エコーシフトの対応するデータポイントに修正する（データポイント修正ステップ）。

　エコーシフトノイズフィルタ３０５は、温度変化に起因するエコーシフトのＳＮ比を増大させることにより、エコーシフトに含まれるノイズを減衰させる（ノイズ減衰ステップ）。エコーシフトノイズフィルタ３０５は、例えば、ローパスフィルタ又はバンドパスフィルタで構成される。エコーシフトノイズフィルタ３０５によりノイズが低減されたエコーシフトは、エコーシフトノイズフィルタ３０５より出力される。

　図９は、ストレイン算出部の構成を具体的に示すブロック図である。ストレイン算出部２０３は、ストレイン計算部３０６及びストレインノイズフィルタ３０７を有している。ストレイン計算部３０６は、上記式６に示すように、深さｄに沿ったエコーシフトψの偏導関数、即ち、ストレインεを計算する（ストレイン算出ステップ）。

　ストレイン計算部３０６は、例えば、加重最小二乗アルゴリズム（ｗｅｉｇｈｔｅｄ　ｌｅａｓｔ　ｓｑｕａｒｅ　ａｌｇｏｒｉｔｈｍ）を用いてストレインを計算する。加重最小二乗アルゴリズムは、各データポイントに関連する重みを適した基準に導くことによる最小二乗法である。ストレイン計算部３０６は、まず、加重最小二乗アルゴリズムを用いて、３個又は３個以上のサンプルからエコーシフトを用いてストレインを計算する。これにより、設定されたデータにより一致させるために一次関数のパラメータを調整する、即ち、モデルとデータとの間の平均平方誤差を最小化する。加重最小二乗アルゴリズムを用いることにより、ストレインεは、

によって決定される。式１９において、Ｎはサンプルの数であり、ψはサンプル点におけるエコーシフトであり、ｄはサンプル点の深さ指数である。

　ＳＮ比を増大させるために、加重線形回帰が導入される。具体的には、ストレイン値を計算するために、各エコーシフトのサンプルには、ＲＦポイントの強度の割合で重みが付けられる。これにより、より高いＳＮ比でサンプル点を強調することができる。従って、式１９のストレインεは、

として修正される。式２０において、Ｉは、サンプルの変位を計算するためのＲＦポイントの強度である。

　ストレインノイズフィルタ３０７は、小さな環境外乱によってストレインに生じたノイズを減衰させる。ストレインノイズフィルタ３０７は、例えば、深さ及び時間に沿った２次元メジアンフィルタである。ストレインノイズフィルタ３０７によりノイズが低減されたストレインは、ストレインノイズフィルタ３０７から出力される。

　図１０は、ストレインレート算出部の構成を具体的に示すブロック図である。図１０に示すように、ストレインレート算出部２０４は、ストレインレート計算部３０８及びストレインレートノイズフィルタ３０９を有している。ストレインレート算出部２０４は、算出されたストレインに基づいて、ストレインレートを算出する（ストレインレート算出ステップ）。上述したように、ストレインレートξは、深さｄ及び時間（スロータイム）ｎに沿ったエコーシフトψの混合偏導関数である（式１４参照）。

　ストレインレート計算部３０８は、エコーシフトから間接的に又は直接的にストレインレートを計算する。間接的にストレインレートを計算する方法としては、例えば、ストレイン（深さに沿ったエコーシフトの偏導関数）又は速度（時間（スロータイム）に沿ったエコーシフトの偏導関数）に基づいてストレインレートを計算する方法があるが、これに限定されない。

　直接的にストレインレートを計算する方法としては、例えば、数値微分法（ｎｕｍｅｒｉｃ　ｄｉｆｆｅｒｅｎｔｉａｔｉｏｎ　ｍｅｔｈｏｄ）がある。エコーシフトψは、深さｄ及び時間（スロータイム）ｎを独立変数として、

で表される。

　エコーシフトの混合偏導関数は、ｄ－ｎ平面で連続的であるので、

と表される。

　式２２の混合偏導関数を数値解析することにより、

が得られる。式２３において、ａ_ｈ＝ｄ_ｉ＋１―ｄ_ｉ、ａ_ｋ＝ｎ_ｊ＋１―ｎ_ｊ＝ｎ_ｊ－ｎ_ｊ－１であり、ｉ，ｊはそれぞれ、データポイントの添え字である。

　直接的にストレインレートを計算する方法の他の例として、最小二乗法（ｌｅａｓｔ　ｓｑｕａｒｅｓ　ｍｅｔｈｏｄ）がある。エコーシフトは、下記の曲線に一致すると仮定される。

　式２４において、

は、ポイント（ｄ_ｉ，ｎ_ｊ）で一致されたエコーシフトの値である。

　ストレインレートξは、

として定義される。

　式２４及び式２５を組み合わせることにより、深さ及び時間（スロータイム）の積と深さとに基づいてエコーシフトを予測するための線形モデル式が得られる。この線形モデル式は、

で表される。一組のデータ（ｄ_１，ｎ_１，ψ_１，１），（ｄ_１，ｎ_２，ψ_１，２），・・・（ｄ_ｎ，ｎ_ｍ，ψ_ｎ，ｍ）について、次式２７のように、スキャン信号に基づいて算出されたエコーシフトψ_ｉ，ｊと、式２６の線形モデル式に基づいて予測されたエコーシフトとの誤差が最小になるようなξの値を求める。

　式２７において、ｋ，ｙはそれぞれ、ストレインレートが算出及び予測された深さ方向におけるサンプル数及びフレーム数である。全てのｄ_ｉ，ｎ_ｉ及びψ_ｉ，ｊは既知である一方、ξ，ｂ及びκは未知係数である。最小二乗誤差を得るために、未知係数ξ，ｂ及びκについて、一次導関数を生成する。

　式２８を解くことにより、ストレインレートξの値を得ることができる。

　ストレインレートノイズフィルタ３０９は、ストレインレートに生じたノイズを減衰させる。ストレインレートノイズフィルタ３０９によりノイズが低減されたストレインレートは、ストレインレートノイズフィルタ３０９から出力される。

　図１１は、温度推定部の構成を具体的に示すブロック図である。図１１に示すように、温度推定部２０５は、パラメータ計算部３１０及び温度計算部３１１を有している。温度推定部２０５は、算出されたストレイン及びストレインレートに基づいて、温度を推定する（温度推定ステップ）。

　パラメータ計算部３１０は、算出されたストレイン及びストレインレートに基づいて、式１６の線形モデル式のアルゴリズムパラメータａ_１２，ａ_１３及びａ_１４を計算する。これらのアルゴリズムパラメータは、次のようなステップを実行することにより同定される。

　まず、加熱すべき対象領域５１に対して熱電対（図示せず）の位置を位置合わせする（位置合わせステップ）。その後、熱電対によって対象領域５１の温度を測定する（温度測定ステップ）。熱電対によって測定された温度は、

であると仮定される。その後、対象領域５１をスキャンすることにより得られるスキャン信号に基づいて、対象領域５１の推定温度を導出する（推定温度導出ステップ）。ストレイン及びストレインレートに基づいて推定される温度は、次の一組のデータ、

におけるストレイン及びストレインレートを式１６に代入することにより求められる。温度測定ステップで測定された温度と推定温度導出ステップで導出された推定温度との誤差は、残差平方和として表現される。

　式２９における関数ｅｒｒを最小化するために、未知係数ａ_１２，ａ_１３及びａ_１４について一次導関数を生成する。

　式３０を解くことにより、アルゴリズムパラメータａ_１２，ａ_１３及びａ_１４が同定される（アルゴリズムパラメータ同定ステップ）。

　温度計算部３１１は、計算されたアルゴリズムパラメータを式１６に適用し、算出されたストレイン及びストレインレートに基づいて温度を計算することにより、温度を推定する。

　なお、上述した位置合わせステップにおいては、温度変化の全範囲をカバーするために、加熱すべき対象領域５１に対して熱電対の位置を正確に位置合わせする必要がある。位置合わせステップにおいて、まず、Ｂモード画像によって、加熱すべき対象領域５１に対して熱電対の位置を位置合わせする。その後、対象領域５１の温度を変化させるために、超音波信号を出力する加熱源（図示せず）を作動させ、加熱源を空間内で変位させながら、熱電対によって対象領域５１の温度を測定する。その後、熱電対により測定された温度が最小となったときに、加熱源の変位を停止させる。このようにして、熱電対と加熱源との相対的な位置関係が決定される。

　（温度推定方法）
　図１２は、本発明の実施の形態１に係る温度推定方法による温度の推定の流れを示すフローチャートである。図１３は、温度推定ステップを説明するための図である。

　まず、超音波信号によって対象領域５１をスキャンすることにより得られるスキャン信号を受信し、外乱によってスキャン信号に生じた破損フレームを除去する（前処理ステップ）（Ｓ１１）。その後、破損フレームが除去されたスキャン信号に基づいてエコーシフトを算出する（エコーシフト算出ステップ）（Ｓ１２）。その後、算出されたエコーシフトに基づいてストレインを算出し（ストレイン算出ステップ）（Ｓ１３）、さらに、算出されたストレインに基づいてストレインレートを算出する（ストレインレート算出ステップ）（Ｓ１４）。その後、図１３に示すように、算出されたストレイン及びストレインレートに基づいて、式１６の線形モデル式を用いて対象領域５１の温度を推定する（温度推定ステップ）（Ｓ１５）。

　以上のような温度推定方法によって、従来よりも高精度に温度を推定することができる。図１４は、本発明の温度推定方法と従来の温度推定方法とを比較した実験結果を示す図である。図１４において、実線のグラフは、実施の形態１の温度推定方法により推定された温度の時間変化を示し、一点鎖線のグラフは、上述した従来の温度推定方法１により推定された温度の時間変化を示している。破線のグラフは、記録された実際の温度の時間変化を示している。図１４から明らかなように、実施の形態１の温度推定方法により推定された温度は、従来の温度推定方法１により推定された温度と比べて、実際の温度との誤差が小さい。この実験結果から、本実施の形態の温度推定方法では、高精度に温度を推定することができることが理解できる。

　図１５は、本発明の温度推定方法と従来の温度推定方法とを比較した実験結果を示す図である。図１５において、上側のグラフは、従来の温度推定方法１により推定された温度と記録された実際の温度との誤差（具体的には、二乗平均平方根誤差）を求める実験を９回行った結果を示している。下側のグラフは、実施の形態１の温度推定方法により推定された温度と記録された実際の温度との誤差を求める実験を９回行った結果を示している。図１５から明らかなように、実施の形態１の温度推定方法により推定された温度の標準偏差は、従来の温度推定方法１により推定された温度の標準偏差よりも小さくなる。この実験結果からも、本実施の形態の温度推定方法では、高精度に温度を推定することができることが理解できる。

　（実施の形態２）
　図１６は、本発明の実施の形態２に係る温度推定装置の構成を示すブロック図である。本実施の形態の温度推定装置では、信号処理部１０５Ａは、さらに、メモリ３２１及び温度補正部３２２を有している。メモリ３２１には、温度推定部２０５により推定された１フレーム分（又は複数フレーム分）の温度が蓄積されている。温度補正部３２２は、温度の空間的な連続性を考慮した目的関数ｆ_１に基づいて、メモリ３２１に蓄積された１フレーム分（又は複数フレーム分）の温度を補正する（温度補正ステップ）。なお、温度補正ステップは、上述した温度推定ステップの後に実行される。

　上述した目的関数ｆ_１は、

で表される。式３１において、

は、温度補正部３２２により補正された温度であり、Ｔ_ｄ，ｌは、温度推定部２０５により推定された温度である。λ_１，ｇはそれぞれ、重み付け関数である。なお、重み付け関数ｇは、例えばＨｕｂｅｒ関数で構成される。ｄは深さ方向における空間位置を表し、ｌはライン方向における空間位置を表している。

　式３１において、第１の項は、温度の補正値を温度の推定値に近付けるための項であり、第２の項は、温度の推定値を空間微分するための項である。温度補正部３２２は、例えば共役勾配法等の数値計算を行うことにより、式３１の目的関数ｆ_１を最小にするような温度の補正値を求める。これにより、温度推定部２０５により推定された温度を、空間的な連続性が考慮された温度に補正することができる。

　（実施の形態３）
　図１７は、本発明の実施の形態３に係る温度推定装置の構成を示すブロック図である。本実施の形態の温度推定装置では、信号処理部１０５Ｂは、上記実施の形態２の温度補正部３２２に代えて、温度補正部３２３を有している。温度補正部３２３は、超音波信号により加熱された熱源位置と対象領域５１との間の距離を考慮した目的関数ｆ_２に基づいて、メモリ３２１に蓄積された１フレーム分（又は複数フレーム分）の温度を補正する（温度補正ステップ）。なお、温度補正ステップは、上述した温度推定ステップの後に実行される。

　上述した目的関数ｆ_２は、

で表される。式３２において、λ_２，αはそれぞれ、重み付け関数である。

　式３２において、第１の項は、温度の補正値を温度の推定値に近付けるための項であり、第２の項は、熱源位置（ｄ_０，ｌ_０）における温度と補正値を求めるべき位置（ｄ，ｌ）における温度との差分であり、第３の項は、熱源位置からの距離の３乗に反比例する項である。温度補正部３２３は、例えば共役勾配法等の数値計算を行うことにより、式３２の目的関数ｆ_２を最小にするような温度の補正値を求める。これにより、温度推定部２０５により推定された温度を、熱源位置と対象領域との間の距離が考慮された温度に補正することができる。

　なお、上記各実施の形態において、各構成要素は、専用のハードウェアで構成されるか、各構成要素に適したソフトウェアプログラムを実行することによって実現されてもよい。各構成要素は、ＣＰＵ又はプロセッサ等のプログラム実行部が、ハードディスク又は半導体メモリ等の記録媒体に記録されたソフトウェアプログラムを読み出して実行することによって実現されてもよい。ここで、上記各実施の形態の温度推定方法等を実現するソフトウェアは、次のようなプログラムである。

　即ち、このプログラムは、超音波信号によって前記対象領域をスキャンすることにより得られるスキャン信号を受信し、前記スキャン信号に基づいて、超音波信号が前記対象領域を通過するのに要する時間が前記対象領域の温度に依存して変化する変化量であるエコーシフトを算出するエコーシフト算出ステップと、算出されたエコーシフトに基づいて、超音波信号が前記対象領域を通過する際の移動距離が前記対象領域の温度に依存して見かけ上変化する変化割合であるストレインを算出するストレイン算出ステップと、算出されたストレインに基づいて、ストレインの時間的な変化割合であるストレインレートを算出するストレインレート算出ステップと、予め定められたストレインとストレインレートと温度との関係を用いて、算出されたストレイン及びストレインレートに対応する前記対象領域の温度を推定する温度推定ステップをコンピュータに実行させる。

　また、上記プログラムをコンピュータ読み取り可能な不揮発性の記録媒体、例えば、フレキシブルディスク、ハードディスク、ＣＤ－ＲＯＭ、ＭＯ、ＤＶＤ、ＤＶＤ－ＲＯＭ、ＤＶＤ－ＲＡＭ、ＢＤ（Ｂｌｕ－ｒａｙ　Ｄｉｓｃ（登録商標））及び半導体メモリ等に記録したものでもよい。

　以上、本発明の一つ又は複数の態様に係る温度推定方法、温度推定装置及びプログラムについて、実施の形態１～３について説明したが、本発明は、これら実施の形態１～３に限定されるものではない。本発明の趣旨を逸脱しない限り、当業者が思い付く各種変形を実施の形態１～３に施したものや、異なる実施の形態における構成要素を組み合わせて構築される形態も、本発明の一つ又は複数の態様の範囲内に含まれてもよい。

　本発明は、高精度に温度を推定することができるとともに、広い温度範囲で温度を推定することができる温度推定方法、温度推定装置及びプログラムとして適用することができる。

１０　温度推定装置
５０　生体
５１　対象領域
１０１　本体部
１０２　超音波振動子
１０３　送信部
１０４　受信部
１０５，１０５Ａ，１０５Ｂ　信号処理部
１０６　表示部
１０７　制御部
２０１　前処理部
２０２　エコーシフト算出部
２０３　ストレイン算出部
２０４　ストレインレート算出部
２０５　温度推定部
３０１　破損フレーム除去部
３０２　ノイズフィルタ
３０３　エコーシフト計算部
３０４　異常値修正部
３０５　エコーシフトノイズフィルタ
３０６　ストレイン計算部
３０７　ストレインノイズフィルタ
３０８　ストレインレート計算部
３０９　ストレインレートノイズフィルタ
３１０　パラメータ計算部
３１１　温度計算部
３２１　メモリ
３２２，３２３　温度補正部

Claims

　超音波信号を用いて対象領域の温度を推定する温度推定方法であって、
　超音波信号によって前記対象領域をスキャンすることにより得られるスキャン信号を受信し、前記スキャン信号に基づいて、超音波信号が前記対象領域を通過するのに要する時間が前記対象領域の温度に依存して変化する変化量であるエコーシフトを算出するエコーシフト算出ステップと、
　算出されたエコーシフトに基づいて、超音波信号が前記対象領域を通過する際の移動距離が前記対象領域の温度に依存して見かけ上変化する変化割合であるストレインを算出するストレイン算出ステップと、
　算出されたストレインに基づいて、ストレインの時間的な変化割合であるストレインレートを算出するストレインレート算出ステップと、
　予め定められたストレインとストレインレートと温度との関係を用いて、算出されたストレイン及びストレインレートに対応する前記対象領域の温度を推定する温度推定ステップと、を含む
　温度推定方法。
　前記温度推定ステップにおいては、ストレインとストレインレートと温度との関係を示す線形モデル式を用いて、算出されたストレイン及びストレインレートに対応する前記対象領域の温度を推定する
　請求項１に記載の温度推定方法。
　さらに、外乱によって前記スキャン信号に生じた破損フレームを除去する前処理ステップを含み、
　前記エコーシフト算出ステップでは、前記前処理ステップにおいて破損フレームが除去されたスキャン信号に基づいて、エコーシフトを算出する
　請求項１又は２に記載の温度推定方法。
　前記エコーシフト算出ステップは、
　受信された前記スキャン信号に基づいて生のエコーシフトを計算するエコーシフト計算ステップと、
　計算された前記生のエコーシフトにおいて温度変化に起因しない異常値を修正することにより、温度変化に起因するエコーシフトを得る異常値修正ステップと、
　前記温度変化に起因するエコーシフトに含まれるノイズをノイズフィルタにより減衰させるノイズ減衰ステップと、を含む
　請求項１～３のいずれか１項に記載の温度推定方法。
　前記異常値修正ステップは、
　前記エコーシフト計算ステップにおいて計算された前記生のエコーシフトに基づいて、予測エコーシフトを計算するためのエコーシフト予測モデルのパラメータを調整するパラメータ調整ステップと、
　最適化された前記エコーシフト予測モデルに基づいて、前記予測エコーシフトを計算する予測エコーシフト計算ステップと、
　前記予測エコーシフト及び前記対象領域をスキャンすることにより得られるＲＦ信号に基づいて、閾値を設定する閾値設定ステップと、
　前記生のエコーシフトの複数のデータポイントのうち、前記閾値を超えた前記異常値であるデータポイントを、前記予測エコーシフトの対応するデータポイントに基づいて修正するデータポイント修正ステップと、を含む
　請求項４に記載の温度推定方法。
　前記エコーシフト予測モデルは、誤差関数又は相補誤差関数である
　請求項５に記載の温度推定方法。
　前記閾値設定ステップは、
　前記閾値としての上側の閾値と下側の閾値との間の間隔を設定する間隔設定ステップと、
　前記ＲＦ信号の強度を重みとして用いることにより、前記上側の閾値及び前記下側の閾値を調節する閾値調節ステップと、を含む
　請求項５又は６に記載の温度推定方法。
　前記ストレインレート算出ステップにおいては、
　深さ及び時間の積と深さとに基づいてエコーシフトを予測するためのモデル式を用いて、前記スキャン信号に基づいて算出されたエコーシフトと前記モデル式に基づいて予測されたエコーシフトとの誤差が最小となるように、ストレインレートを算出する
　請求項１～７のいずれか１項に記載の温度推定方法。
　前記ストレインレート算出ステップにおいて用いられる前記モデル式は、線形モデル式である
　請求項８に記載の温度推定方法。
　前記温度推定ステップにおいて用いられる前記線形モデル式のアルゴリズムパラメータは、
　加熱すべき前記対象領域に対して熱電対の位置を位置合わせする位置合わせステップと、
　前記熱電対によって前記対象領域の温度を測定する温度測定ステップと、
　前記対象領域をスキャンすることにより得られる前記スキャン信号に基づいて、前記対象領域の推定温度を導出する推定温度導出ステップと、
　前記温度測定ステップにおいて測定された温度と前記推定温度導出ステップにおいて導出された推定温度との誤差が最小となるように、前記アルゴリズムパラメータを同定するアルゴリズムパラメータ同定ステップと、を実行することにより同定される
　請求項２に記載の温度推定方法。
　前記位置合わせステップは、
　Ｂモード画像によって、加熱すべき前記対象領域に対して前記熱電対の位置を位置合わせするステップと、
　前記対象領域の温度を変化させるために、超音波信号を出力する加熱源を作動させるステップと、
　前記加熱源を空間内で変位させながら、前記熱電対によって前記対象領域の温度を測定するステップと、
　前記熱電対により測定された温度が最小となったときに前記加熱源の変位を停止させるステップと、を含む
　請求項１０に記載の温度推定方法。
　さらに、温度の空間的な連続性を考慮した目的関数に基づいて、前記温度推定ステップにより推定された温度を補正する温度補正ステップを含む
　請求項１～１１のいずれか１項に記載の温度推定方法。
　さらに、加熱用の超音波信号により加熱された熱源位置と前記対象領域との間の距離を考慮した目的関数に基づいて、前記温度推定ステップにより推定された温度を補正する温度補正ステップを含む
　請求項１～１１のいずれか１項に記載の温度推定方法。
　超音波信号を用いて対象領域の温度を推定する温度推定装置であって、
　超音波信号によって前記対象領域をスキャンすることにより得られるスキャン信号を受信し、前記スキャン信号に基づいて、超音波信号が前記対象領域を通過するのに要する時間が前記対象領域の温度に依存して変化する変化量であるエコーシフトを算出するエコーシフト算出部と、
　算出されたエコーシフトに基づいて、超音波信号が前記対象領域を通過する際の移動距離が前記対象領域の温度に依存して見かけ上変化する変化割合であるストレインを算出するストレイン算出部と、
　算出されたストレインに基づいて、ストレインの時間的な変化割合であるストレインレートを算出するストレインレート算出部と、
　予め定められたストレインとストレインレートと温度との関係を用いて、算出されたストレイン及びストレインレートに対応する前記対象領域の温度を推定する温度推定部と、を備える
　温度推定装置。
　超音波信号を用いて対象領域の温度を推定するためのプログラムであって、
　超音波信号によって前記対象領域をスキャンすることにより得られるスキャン信号を受信し、前記スキャン信号に基づいて、超音波信号が前記対象領域を通過するのに要する時間が前記対象領域の温度に依存して変化する変化量であるエコーシフトを算出するエコーシフト算出ステップと、
　算出されたエコーシフトに基づいて、超音波信号が前記対象領域を通過する際の移動距離が前記対象領域の温度に依存して見かけ上変化する変化割合であるストレインを算出するストレイン算出ステップと、
　算出されたストレインに基づいて、ストレインの時間的な変化割合であるストレインレートを算出するストレインレート算出ステップと、
　予め定められたストレインとストレインレートと温度との関係を用いて、算出されたストレイン及びストレインレートに対応する前記対象領域の温度を推定する温度推定ステップと、をコンピュータに実行させる
　プログラム。