IP Force 特許公報掲載プロジェクト 2022.1.31 β版

知財求人 - 知財ポータルサイト「IP Force」

▶ スーパー ソニック イマジンの特許一覧

特許7237116超音波減衰パラメータを推定するための方法及びシステム
<>
  • 特許-超音波減衰パラメータを推定するための方法及びシステム 図1
  • 特許-超音波減衰パラメータを推定するための方法及びシステム 図2
  • 特許-超音波減衰パラメータを推定するための方法及びシステム 図3
  • 特許-超音波減衰パラメータを推定するための方法及びシステム 図4
  • 特許-超音波減衰パラメータを推定するための方法及びシステム 図5
  • 特許-超音波減衰パラメータを推定するための方法及びシステム 図6
< >
(19)【発行国】日本国特許庁(JP)
(12)【公報種別】特許公報(B2)
(11)【特許番号】
(24)【登録日】2023-03-02
(45)【発行日】2023-03-10
(54)【発明の名称】超音波減衰パラメータを推定するための方法及びシステム
(51)【国際特許分類】
   A61B 8/14 20060101AFI20230303BHJP
【FI】
A61B8/14
【請求項の数】 17
【外国語出願】
(21)【出願番号】P 2021114177
(22)【出願日】2021-07-09
(65)【公開番号】P2022016414
(43)【公開日】2022-01-21
【審査請求日】2021-11-09
(31)【優先権主張番号】20315348
(32)【優先日】2020-07-10
(33)【優先権主張国・地域又は機関】EP
(73)【特許権者】
【識別番号】508291928
【氏名又は名称】スーパー ソニック イマジン
【氏名又は名称原語表記】SUPER SONIC IMAGINE
(74)【代理人】
【識別番号】100145403
【弁理士】
【氏名又は名称】山尾 憲人
(74)【代理人】
【識別番号】100135703
【弁理士】
【氏名又は名称】岡部 英隆
(72)【発明者】
【氏名】フラスキーニ,クリストフ
【審査官】最首 祐樹
(56)【参考文献】
【文献】特開昭61-025536(JP,A)
【文献】特開昭60-153849(JP,A)
【文献】KIM Hyungsuk, et al.,Attenuation Estimation using Spectral Cross-Correlation,IEEE transactions on ultrasonics, ferroelectrics, and frequency control,米国,IEEE,2007年03月26日,Volume 54, Issue 3,p510-519,https://ieeexplorer.ieee.org/abstract/document/4139331
(58)【調査した分野】(Int.Cl.,DB名)
A61B 8/00-8/15
(57)【特許請求の範囲】
【請求項1】
媒体中の領域における超音波減衰パラメータを推定するための方法であって、上記方法は、少なくとも1つの超音波トランスデューサ(2)に関連付けられた処理装置(8)によって実施され、上記方法は、
(a)上記媒体において少なくとも1つのパルスがトランスデューサによって送信される送信ステップ(101)と、
(b)上記パルスに応答してデータがトランスデューサによって取得される受信ステップ(102)と、
(c)上記処理装置によって上記データが処理されて上記領域の後方散乱された取得データを提供する処理ステップ(103)と、
(d)時空間領域における深度の関数である、上記後方散乱された取得データの自己相関関数がゼロ遅延で決定される関数決定ステップ(104)と、
(e)上記自己相関関数に基づいて超音波減衰パラメータが推定される減衰推定ステップ(105)とを含む、
方法。
【請求項2】
上記処理ステップ(c)は、上記データがビームフォーミング処理によって処理されて上記領域のビームフォーミングされた取得データを提供するビームフォーミングステップ(103)を含む、
請求項1記載の方法。
【請求項3】
上記減衰推定ステップ(e)において、上記超音波減衰パラメータは、深度方向の所定区間(Δz)にわたって推定される、
請求項1又は2記載の方法。
【請求項4】
上記減衰推定ステップ(e)において、上記超音波減衰パラメータは、上記自己相関関数を解くことで決定される、
請求項1~3のうちの1つに記載の方法。
【請求項5】
上記後方散乱された取得データの自己相関関数は、深度及び時間の関数である、
請求項1~4のうちの1つに記載の方法。
【請求項6】
上記関数決定ステップ(d)において、上記後方散乱された取得データにおける残留回折パターンは、予め決められた回析基準データに基づいて補償される、
請求項1~5のうちの1つに記載の方法。
【請求項7】
上記後方散乱された取得データにおける残留回析パターンは、予め決められた基準超音波減衰パラメータを有する基準媒体に対してステップ(a)~(e)を含む方法を適用し、上記自己相関関数を比較することで補償される、
請求項記載の方法。
【請求項8】
上記自己相関関数を比較することは、
上記媒体及び上記基準媒体における2つの異なる深度 及びz に関する自己相関関数の間のデルタ関数ΔE(z,z)を計算することと、
上記デルタ関数ΔE(z,z)=を解いて上記媒体の超音波減衰パラメータを決定することとを含む、
請求項7記載の方法。
【請求項9】
上記予め決められた回析基準データは、基準媒体における少なくとも2つの異なる深度(z,z)に係る予め決められたパラメータを含む、
請求項6~8のうちの1つに記載の方法。
【請求項10】
上記パルスは、広帯域パルス及び/又は複数のパルスである、
請求項1~9のうちの1つに記載の方法。
【請求項11】
媒体中の領域における局所的超音波減衰パラメータを推定するための方法であって、上記方法は、請求項1~10のうちの1つに記載の方法を含み、
上記方法は、上記媒体内の深度方向(z)における複数の軸線のそれぞれについて、
上記自己相関関数を、駆動入力を有するノイズありの状態空間モデルとしてモデリングすることと、
上記自己相関関数を、チホノフ正則化法、リッジ回帰法、及び/又はベイジアン回帰法を含む正則化法によって正則化することと、
サビツキー・ゴーレイ(Savitzky-Golay)フィルタを用いて上記自己相関関数をフィルタリングすることと、モンテカルロ法と、特異スペクトル分析と、
圧縮された検出及び/又は疎サンプリング法と、
のうちの少なくとも1つを実行することにより、Eyを上記自己相関関数としかつzを上記深度とするとき、数式F(Ey(z),β)=0における局所的超音波減衰βを推定することを含む、
方法。
【請求項12】
上記自己相関関数が、ノイズありの状態空間モデルとしてモデリングされる場合、上記駆動入力は、基準媒体の予め決められた減衰パラメータを含む、
請求項11記載の方法。
【請求項13】
上記モデリングされた状態空間モデルは、カルマンフィルタ又は予め定義されたベイジアン法を用いて、局所的超音波減衰パラメータを決定する、
請求項11又は12記載の方法。
【請求項14】
上記媒体の反射率パターンは、予め定義された収束しきい値が達成されるまで、上記モデリングされた状態空間モデルに対してEMアルゴリズムを繰り返すことで補償され、
上記EMアルゴリズムは、反射率パラメータを出力しかつ駆動入力として超音波減衰パラメータを用いる第1のカルマンフィルタと、超音波減衰パラメータを出力しかつ駆動入力として反射率パラメータを用いる第2のカルマンフィルタとを使用し、
両方のフィルタは、駆動入力として他方の出力をそれぞれ用いるように接続される、
請求項11~13のうちの1つに記載の方法。
【請求項15】
媒体の領域内における局所的超音波減衰の画像生成を行うための方法であって、請求項11~14のうちの1つに記載の方法を含み、
局所的減衰パラメータに基づいて減衰画像が構築される、
方法。
【請求項16】
媒体の領域内における局所的超音波減衰の画像生成を行うための請求項15記載の方法であって、
少なくとも1つの予め定義されたカテゴリにそれぞれ関連付けられた1つ又は複数の関心対象領域に画像をセグメント化するようにトレーニングされた機械学習モデルに減衰画像を供給することを含む、
方法。
【請求項17】
媒体中の領域における超音波減衰パラメータを推定するためのシステムであって、上記システムは、少なくとも1つの超音波トランスデューサ(2)に関連付けられた処理装置(8)を備え、
上記システムは、
(a)上記媒体においてトランスデューサにより少なくとも1つのパルスを送信し、
(b)上記パルスに応答してトランスデューサによりデータを取得し、
(c)上記処理装置によりデータを処理して上記領域の後方散乱された取得データを提供し、
(d)時空間領域における深度の関数である、上記後方散乱された取得データの自己相関関数をゼロ遅延で決定し、
(e)上記自己相関関数に基づいて超音波減衰パラメータを推定するように構成された、
システム。
【発明の詳細な説明】
【技術分野】
【0001】
本発明は、特に医療用イメージング(画像生成)のための、イメージング方法及び上記方法を実現する装置に関する。
【0002】
本開示は、特に、媒体中の領域における超音波減衰パラメータを推定するための方法に関し、より正確には、媒体の内部の領域における局所的な超音波減衰パラメータを推定して画像生成するための方法に関する。本方法は、上記媒体に対して所定関係を有する複数のトランスデューサ(例えば、ライン又はアレイ)に関連付けられた処理装置によって実現されてもよい。
【背景技術】
【0003】
古典的な超音波イメージングは、所与の点に合焦する円筒波を用いて媒体に対して超音波照射(insonification)を行うことからなる。この単一の超音波照射の後方散乱エコーを用いて、画像のライン全体が動的受信ビームフォーミング処理を用いて計算される。画像全体を構築するために、この手順は、所与の深度において横方向のラインに沿って走査する、合焦された一組の波(焦点面と呼ぶ)を送信することで繰り返される。合焦された各波について、動的ビームフォーミングが実行され、画像全体が1ラインずつ構築されて取得される。動的ビームフォーミングは、受信モードでは一様な合焦を保証するが、送信モードでは、焦点は所与の深度に固定される。最終的な画像は、焦点面において、また、焦点軸長に対応する媒体の限られた領域において最適である。しかしながら、回析の法則によって課されるこのエリアの外部では、画像品質は、(合焦されたビームの近傍界及び遠方界における)他の深度では急速に劣化する。
【0004】
この制限を克服するために、古典的な解決方法は、多焦点イメージングを実行するというものであり、画像全体にわたって均質な品質を得るために、異なる送信焦点深度が使用される。所与の焦点深度における各送信は、軸焦点距離によって定められた範囲を有する領域における部分画像の実施を可能にする最終的な画像は、様々な深度に対応するこれらの部分画像を再結合したものを用いて取得される。最適な多焦点画像は、典型的には、数十個の焦点面を必要とする。このことは、典型的には10フレーム/秒未満のフレームレート制限をもたらし、この制限は超音波イメージングでは許容できない。画像品質とフレームレートとの間での良好な妥協点は、約4つの焦点深度画像である。
【0005】
合成動的送信合焦を実行することで、画像品質の改善を期待することができる。そのようなアプローチは、互いに異なる一組の超音波照射のビームフォーミングを行って、次いで結合することにより、動的送信合焦(すなわち、画像中の画素数と同じ個数の焦点深度)を再合成することにある。
【0006】
さらに、本願出願人によって出願された特許文献1によれば、改善された合成超音波イメージング方法が知られ、これは、例えば特許文献2に開示される従来技術の平面波合成超音波イメージング方法を改善することを可能にする。
【0007】
特許文献1は、少なくとも以下のステップを含む超音波イメージングのための方法を提案する。
a)画像生成される領域に対して複数の超音波が送信され、各超音波に応答してトランスデューサのアレイによって各組の生データが取得される送信ステップであって、上記複数の超音波が互いに異なる空間周波数内容を有し、領域における画像生成される複数の場所のそれぞれについて、各組の生データが、対応する超音波に応答してトランスデューサによって受信された時間信号を表す送信ステップ。
b)画像生成される領域における複数の仮想送信焦点ゾーンのそれぞれについて、複数組の生データから少なくとも一組のコヒーレントデータが合成されるコヒーレンス向上ステップ。
c)仮想送信焦点ゾーンのそれぞれに含まれる複数の場所のそれぞれについて、一組のコヒーレントデータを用いてビームフォーミングによって画像の画素が計算されるビームフォーミングステップ。
【0008】
これらの性質のおかげで、生データの空間的コヒーレンスは、ビームフォーミング前にステップb)において回復され、それにより、様々な超音波の送信から受信されたデータを正確に結合することを可能にする。空間的コヒーレンスを回復する必要性は、画像生成される領域を空間的に広く拡散された波動場で照射する場合、媒体から戻るエコーを、画像生成される領域においてランダムに分散されたインコヒーレントなソース(散乱器)から発生された波動場とみなすことができるという事実に起因し、したがって、波動場の空間的コヒーレンスは生データでは失われる(又は非常に低くなる)。次いで、ビームフォーミングは、コヒーレンス回復ステップから得られるコヒーレントデータに対して実行されてもよく、その結果、より正確な画像をもたらす。
【0009】
合成ビームフォーミングも呼ばれるこの技術に基づいて、特許文献3は、低減されたスペックルノイズを有して、媒体の内部の領域の画像を生成するもう1つのイメージング方法を開示する。このため、本方法は以下のステップを含む。
(a)媒体の内部においてトランスデューサによって第1の複数の波が送信される送信ステップ。
(b)波に応答して上記トランスデューサによって一組のデータが取得される受信ステップ。
(c)第2の複数のビームフォーミング処理によって一組のデータが処理されて、画像の少なくとも一部のビームフォーミングされた画素値を提供するビームフォーミングステップであって、各ビームフォーミング処理は、送信重み付けベクトルで生成された波に対応する一組のデータを用いるか、又は、ビームフォーミングされた画素値の計算において送信重み付けベクトルを用いるビームフォーミングステップ。
(d)上記第2の複数のものに係るビームフォーミングされた画素値が互いに結合されて、画像の内部の画素の画素値を提供する結合ステップであって、送信重み付けベクトルが互いに異なりかつ直交する結合ステップ。
【0010】
これらの特徴のおかげで、各送信重み付けベクトルは、無相関のスペックルノイズを生成し、重み付けられたデータの結合は、低減されたスペックルノイズを有する領域の画像を計算することを可能にする。
【先行技術文献】
【特許文献】
【0011】
【文献】欧州特許出願公開第2101191号公報
【文献】米国特許第6551246号明細書
【文献】国際公開第2017/093778号
【非特許文献】
【0012】
【文献】LYMBERIS ET AL: "Estimation of frequency-dependent attenuation based on parametric spectral analysis and correlation lags of the demodulated echo signal", ULTRASONIC IMAGING, DYNAMEDIA INC., SILVER SPRING, MD, US, vol. 13, no. 1, 1 January 1991 (1991-01-01), pages 1-26, XP026409816, ISSN: 0161-7346.
【文献】HYUNGSUK KIM ET AL: "Attenuation estimation using spectral cross-correlation", IEEE TRANSACTIONS ON ULTRASONICS, FERROELECTRICS AND FREQUENCY CONTROL, IEEE, US, vol. 54, no. 3, 1 March 2007 (2007-03-01), pages 510-519, XP011175817, ISSN: 0885-3010, DOI: 10.1109/TUFFC.2007.274.
【発明の概要】
【発明が解決しようとする課題】
【0013】
しかしながら、超音波イメージングにおいて望ましくは考慮されるべき別の現象は、検査される媒体の内部における超音波減衰である。超音波減衰は、送信された超音波への応答に対して直接的に影響する。それによって、減衰は、敏感な、周波数及び深度に依存する現象を構成する。したがって、組織減衰を考慮して例えば時間利得補償によって従来行われるように、結果的に得られる計算された画像に対する減衰のいかなる影響も補償することが望ましい。
【0014】
一方、超音波減衰の信頼できる推定は、他の超音波診断目的で使用されてもよい。例えば、検査された媒体について、例えばヒトの肝臓又は筋肉について推定された大域的減衰パラメータは、その脂肪含有量を決定するために使用されてもよい。さらに、例えば画像を計算するために処理された、局所的な減衰の分布は、例えばヒトの乳房における、癌を検出するために使用されてもよい。
【0015】
しかしながら、超音波減衰を推定するための従来の方法は、通常、周波数(又はスペクトル)分析を使用する。しかしながら、そのような分析はいくつかの欠点を意味する。
【0016】
減衰が、敏感な、周波数及び深度に依存する現象であるので、正確な評価には、十分に大きな空間及び周波数の観測範囲を必要とする。このことは、大きな周波数及び空間的範囲の使用が望まれるかもしれないことを意味する。しかしながら、減衰を理由として、空間的範囲が大きくなり、また、周波数内容の変動が大きくなると、バイアスが大きくなる。一方、空間的範囲が小さくなると、スペクトル解析が不十分になる。さらに、空間的範囲が大きくなると、回析の影響が大きくなる。さらに、スペックルノイズは、減衰推定の変動を増大させる。
【0017】
従って、この必要なトレードオフを考慮すると、周波数(又はスペクトル)分析によって超音波減衰を推定することは、不正確な結果をもたらし、それと同時に、特に、必要とされるフーリエ変換に起因して、非常に大きな計算上のコストを必要とする。結果として、可能な結果の変動が増大し、従って、推定の信頼性が低下する、すなわち再現可能性が低下するという問題がある。従って、取得時間と精度のレベルとの間のトレードオフが考慮されなければならない。A.Lymberis他は、周波数に依存する減衰を推定する方法を提案する。特に、提案された方法2において、平均周波数推定量は、復調信号の利用可能なすべての自己相関遅延(ACn)から導出された。非特許文献1を参照。さらに、Hyungsuk他は、スペクトル相互相関を用いた減衰推定方法を開示する。非特許文献2を参照。
【0018】
現在、前述の問題を克服し、特に、媒体中の領域における超音波減衰パラメータを高い信頼性で推定するための方法及びシステムであって、より高速であり、計算上のコストがより低く、したがって、必要な処理パワーがより小さいという利点を有する方法及びシステムを提供することが引き続き望まれている。
【課題を解決するための手段】
【0019】
従って、本開示の実施形態によれば、媒体中の領域における超音波減衰パラメータを推定するための方法が提供される。上記方法は、少なくとも1つの超音波トランスデューサ(これは上記媒体と所定関係にあってもよい)に関連付けられた処理装置によって実施される。本方法は、
(a)上記媒体において少なくとも1つのパルスがトランスデューサによって送信される送信ステップと、
(b)上記パルスに応答してデータがトランスデューサによって取得される受信ステップと、
(c)上記処理装置によって上記データが処理されて上記領域の後方散乱された取得データを提供する処理ステップと、
(d)時空間領域における深度の関数である、上記後方散乱された取得データの自己相関関数がゼロ(0)遅延で決定される関数決定ステップと、
(e)上記自己相関関数に基づいて超音波減衰パラメータが推定される減衰推定ステップとを含む。
【0020】
上記方法は、例えば、検査された領域全体に係る、例えば超音波画像全体に係る1つの減衰パラメータの大域的減衰推定方法と呼ばれてもよい。
【0021】
そのような方法を提供することによって、周波数(又はスペクトル)領域の代わりに時空間領域において分析を実行することで超音波減衰パラメータを推定することが可能になる。従って、推定は、さらにずっと正確になるという利点があり、より少ない計算コストをともない、これにより、例えばリアルタイムの、又は少なくとも準リアルタイムの計算モードを改善する。さらに、増大した精度に起因して、変動の減少を達成することができ、したがって、増大した再現性を達成することができる。
【0022】
処理ステップ(c)は、データがビームフォーミング処理によって処理されて領域のビームフォーミングされた取得データを提供するビームフォーミングステップを含んでもよい。上記ビームフォーミングされた取得データ又は後方散乱された取得データは、例えば画素値又はボクセル値の形式で提示されてもよい。
【0023】
ビームフォーミング処理に起因して、取得データにおける回析パターンを低減することが可能になる。
【0024】
ビームフォーミング処理は、例えば、合成ビームフォーミング処理であってもよい。これは、回析パターンをさらに低減することを可能にする。
【0025】
さらに、スペックルノイズをさらに低減するようにステップ(c)のビームフォーミング処理を構成することが可能である。例えば、ステップ(c)は、
(c1)複数のビームフォーミング処理によってデータが処理されて、少なくとも領域のビームフォーミングされた取得データを提供するビームフォーミングステップであって、各ビームフォーミング処理は、送信重み付けベクトルで生成されたパルスに対応するデータを用いるか、又は、ビームフォーミングされた取得データの計算において送信重み付けベクトルを用いるビームフォーミングステップと、
(c2)上記複数のものに係るビームフォーミングされた取得データが互いに結合されて、生成された超音波画像の内部の画素の画素値を提供する結合ステップであって、送信重み付けベクトルは互いに異なりかつ直交する結合ステップとを含む。
【0026】
しかしながら、後方散乱された取得データを提供するために、上述したビームフォーミング処理の代替処理が使用されてもよい。例えば、音響レンズを介して送信パルス(及びオプションで受信信号)を音響的にビーム整形することが使用されてもよく、又は、機械的な掃引/走査が使用されてもよい。そのような場合、処理ステップ(c)においてビームフォーミング処理は省略されてもよい。
【0027】
減衰推定ステップ(e)において、超音波減衰パラメータは、深度方向の所定区間Δzにわたって推定されてもよい。例えば、上記区間Δzは、媒体における指定された第1の深度及び第2の深度の間で予め定義されてもよい。媒体における深度は、プローブの表面と関心対象点(空間的場所)との間の距離として定義されてもよいことに注意すべきである。
【0028】
減衰決定ステップにおいて、超音波減衰パラメータは、例えば区間Δzに係る、自己相関関数を解くことで決定されてもよい。
【0029】
後方散乱された取得データの自己相関関数は、深度及び時間の関数であってもよい。
【0030】
関数決定ステップ(d)において、後方散乱された取得データにおける残留回折パターンは、予め決められた回析基準データに基づいて補償されてもよい。
【0031】
後方散乱された取得データにおける残留回析パターンは、予め決められた基準超音波減衰パラメータを有する基準媒体に対してステップ(a’)~(e’)を含む方法を適用し、自己相関関数を比較することで補償される。
【0032】
言いかえれば、推定方法は、予め決められた、すなわち既知の超音波減衰パラメータを用いて、そのような基準媒体、例えばファントムに適用されてもよい。従って、決定された自己相関関数は基準自己相関関数と呼ばれてもよく、この自己相関関数を解くことで、後方散乱された取得データにおける残留回析パターンを決定してもよく、及び/又は、検査される関心対象領域を含む媒体の自己相関関数においてそれを補償してもよい。
【0033】
自己相関関数を比較することは、両方の関数の間のデルタ関数を計算してデルタ関数を解き、媒体の超音波減衰パラメータを決定することを含んでもよい。
【0034】
例えば、基準超音波減衰パラメータは、検査される媒体における区間Δzに類似した基準媒体における区間Δz’に係るパラメータであってもよい。
【0035】
予め決められた回析基準データは、基準媒体における異なる深度に係る予め決められたパラメータを含んでもよい。
【0036】
パルスは、広帯域パルス及び/又は複数のパルスであってもよい。例えば、1つ又は複数の広帯域パルスを使用することは、空間的な限定、すなわちより正確な区間決定を可能にするので、推定の精度を増大させうる。
【0037】
本開示はさらに、媒体中の領域における局所的超音波減衰パラメータを推定するための方法に関する。上記方法は上述した大域的減衰推定方法(すなわちステップ(a)~(e))を含んでもよい。媒体内の深度方向zにおける複数の軸線(例えば、ステップ(a)~(c)において走査されたライン)のそれぞれについて、
・上記自己相関関数を、駆動入力を有するノイズありの状態空間モデルとしてモデリングすることと、
・上記自己相関関数を、チホノフ正則化法、リッジ回帰法、及び/又はベイジアン回帰法を含む正則化法によって正則化することと、
・サビツキー・ゴーレイ(Savitzky-Golay)フィルタを用いて上記自己相関関数をフィルタリングすることと、モンテカルロ法と、
特異スペクトル分析と、
圧縮された検出及び/又は疎サンプリング法と、
のうちの少なくとも1つによって、Eyを自己相関関数としかつzを深度とするとき、数式F(Ey(z),B)=0における局所的超音波減衰Bを推定することを含む。
【0038】
上記方法は、例えば、検査された領域全体に係る、例えば1つの超音波画像に係る複数の減衰パラメータの局所的減衰推定方法と呼ばれてもよい。上記パラメータは、超音波画像における画素又はボクセルに割り当てられてもよい。さらに、それらは、検査された領域の追加の減衰マップの超音波画像におけるセグメント又はスーパー画素又はクラスタを形成してもよい。
【0039】
この方法のおかげで、関心対象領域に係る局所的超音波減衰パラメータを推定することが可能になり、領域にわたる減衰特性のマップ及び/又は1つ又は複数の画像を生成できるようになる。
【0040】
数式F(Ey(z),B)を解くために使用されうる従来のインバージョン法と比較して、提案した方法は、より速く収束するので、より少ない計算負荷を意味する。同時に、それらは、例えば、取得した局所推定値のより高い空間分解能について、より正確な結果をもたらす。
【0041】
自己相関関数は、局所的減衰関数を取得するために、深度方向(z)において離散化されてもよい。
【0042】
駆動入力は、基準媒体の予め決められた減衰パラメータを含んでもよい。
【0043】
モデリングされた状態空間モデルは、カルマンフィルタ又は予め定義されたベイジアン法を用いて、局所的超音波減衰パラメータを決定してもよい。
【0044】
カルマンフィルタを使用することは、最小平均二乗誤差基準を最適化するという利点を有しうる。
【0045】
媒体の反射率パターンは、予め定義された収束しきい値が達成されるまで、モデリングされた状態空間モデルに対してEMアルゴリズムを繰り返すことで補償されてもよい。例えば、EMアルゴリズムは、反射率パラメータを出力する第1のカルマンフィルタを使用してもよく、駆動入力として超音波減衰パラメータを使用してもよい。第2のカルマンフィルタは、超音波減衰パラメータを出力してもよく、駆動入力として反射率パラメータを使用してもよい。両方のフィルタは、例えば、予め定義された収束しきい値が達成されるまで、駆動入力として他方の出力をそれぞれ用いるように接続されてもよい。
【0046】
本開示はさらに、媒体の領域内における局所的超音波減衰の画像生成を行うための方法に関する。本方法は、前述した局所的減衰推定方法を含んでもよい。減衰画像は、例えば各軸線について、局所的減衰パラメータに基づいて構築されてもよい。
【0047】
媒体の領域内における局所的超音波減衰の画像生成を行うための方法は、少なくとも1つの予め定義されたカテゴリにそれぞれ関連付けられた1つ又は複数の関心対象領域に画像をセグメント化するようにトレーニング(学習)された機械学習モデルに減衰画像を供給することをさらに含んでもよい。
【0048】
上記機械学習モデルは、例えば人工ニューラルネットワーク、例えば畳み込みニューラルネットワークを含んでもよい(モデルは例えばコンピュータで実装されてもよい)。モデルは、トレーニングデータとして例えばそれぞれ注釈された画像を用いる教師ありトレーニング方法によって、又は、教師なしトレーニング方法によってトレーニングされてもよい。
【0049】
最後に、本開示は、(例えば上記媒体と所定関係にある)少なくとも1つの超音波トランスデューサ(2)に関連付けられているか接続されている処理装置を備える、媒体中の領域における超音波減衰パラメータを推定するためのシステムに関する。本システム又は処理装置は、
(a)上記媒体においてトランスデューサにより少なくとも1つのパルスを送信し、
(b)上記パルスに応答してトランスデューサによりデータを取得し、
(c)上記処理装置によりデータを処理して上記領域の後方散乱取得データを提供し、
(d)時空間領域における深度の関数である、上記後方散乱された取得データの自己相関関数をゼロ(0)遅延で決定し、
(e)上記自己相関関数に基づいて超音波減衰パラメータを推定するように構成される。
【0050】
本システムは、オプションで、少なくとも1つのトランスデューサを備えてもよい。
【0051】
少なくとも1つのトランスデューサは、パルスを送信して組織の応答を受信するように構成された単一のトランスデューサであってもよい。例えば、合焦されたトランスデューサは、例えば、凹面形式を有してもよく、又は、各レンズを有してもよい。また、単一のトランスデューサを掃引することも可能である。
【0052】
複数のトランスデューサ及び/又はトランスデューサアレイ2を用いることも可能である。例えば、典型的には、軸X(水平又はアレイ方向X)に沿って並置された数十個のトランスデューサ(例えば100~300)を含む、リニアアレイが提供されてもよい。本開示の実施例では、3Dプローブが使用されてもよい。
【0053】
パルスを送信して応答を受信するために1つ又は複数の同じトランスデューサが使用されてもよく、又は、送信及び受信のために異なるトランスデューサが使用されてもよい。
【0054】
本開示はさらに、コンピュータによって実行されるとき、上述した方法のうちの少なくとも1つに係るステップを実行するための命令を含むコンピュータプログラムに関してもよい。
【0055】
最後に、本開示はさらに、コンピュータによって読み取り可能な記録媒体であって、コンピュータによって実行されるとき、上述した方法のうちの少なくとも1つに係るステップを実行するための命令を含むコンピュータプログラムをその上に記録した記録媒体に関してもよい。
【0056】
矛盾しない限り、上述した構成要素及び本明細書内で説明したものが組み合わされもよいことが意図される。
【0057】
先述の概要説明及び後述の詳細説明の両方は例示及び説明にすぎず、実例の提示を目的とし、請求項に記載された開示の限定ではないことは理解されるべきである。
【0058】
添付図面は、本明細書に組み込まれ、本明細書の一部を構成し、その説明とともに本開示の実施形態を示し、その原理をサポート及び例証として提示される。
【図面の簡単な説明】
【0059】
図1】本開示の実施形態に係る超音波装置を示す概略図である。
図2図1の装置の一部を示すブロック図である。
図3】本開示に係る超音波減衰パラメータを推定するための方法であって、図1の装置で実装された方法を示す図である。
図4】本開示に係る方法の第1の例示的な実施形態(大域的減衰推定)を示す図である。
図5図3又は図4の方法によって取得された超音波画像の例を示す。
図6】本開示に係る方法の第2の例示的な実施形態(局所的減衰推定)を示す図である。
【発明を実施するための形態】
【0060】
ここで、本開示の例示的な実施形態が詳細に参照され、その実施例は添付図面に図示されている。図面の全体にわたって可能であればどこでも、同じ部分又は同様の部分を示すために同じ参照番号が使用される。
【0061】
図1に示す装置は、媒体の領域1、例えば、生体組織及び特に患者のヒト組織の画像生成を行うように適応化される。本装置は、例えば、以下の構成要素を含んでもよい。
【0062】
・少なくとも1つのトランスデューサは、例えば、パルスを送信して組織の応答を受信するように構成された単一のトランスデューサ。また、複数のトランスデューサ及び/又はトランスデューサアレイ2を用いることも可能である。例えば、典型的には、通常のプローブで既知のように、軸X(水平又はアレイ方向X)に沿って並置された数十個のトランスデューサ(例えば100~300)を含む、リニアアレイが提供されてもよい。この場合、アレイ2は、領域1の二次元(2D)イメージングを実行するように適応化されるが、アレイ2は、領域1の3Dイメージングを実行するように適応化された二次元アレイであってもよい。トランスデューサアレイ2は、曲線に沿って整列された複数のトランスデューサを含む凸面のアレイであってもよい。パルスを送信して応答を受信するために1つ又は複数の同じトランスデューサが使用されてもよく、又は、送信及び受信のために異なるトランスデューサが使用されてもよい。
・トランスデューサアレイを制御してそこから信号を取得する電子的ベイ3。
・電子的ベイ3を制御し、例えば、電子的ベイから取得された画像を見るためのマイクロコンピュータ4(変形例では、単一の電子装置が、電子的ベイ3及びマイクロコンピュータ4のすべての機能を満たしてもよい)。
【0063】
図1における軸Zは、軸Xに垂直な軸であり、それは、通常、アレイのトランスデューサによって生成された超音波ビームの方向であり、例えば、検査された媒体の深度方向にある。この方向は、本開示では、垂直方向又は軸方向と呼ばれる。
【0064】
図2に示すように、電子的ベイ3は、例えば、以下の構成要素を含んでもよい。
【0065】
・トランスデューサアレイ2のL個のトランスデューサ(Tl-TL)に個々に接続されたL個のアナログ/ディジタル変換器5(A/Di-A/DL);
・n個のアナログ/ディジタル変換器5にそれぞれ接続されたL個のバッファメモリー6(Bi-Bn);
・バッファメモリー6及びマイクロコンピュータ4と通信する中央処理装置8(CPU);
・中央処理装置8に接続されたメモリ9(MEM);
・中央処理装置8に接続されたディジタル信号プロセッサ10(DSP)。
【0066】
本願で開示した装置は超音波イメージングのための装置であり、トランスデューサは超音波トランスデューサであり、実現された方法は、領域1に係る超音波減衰パラメータを推定し、また、オプションで、領域1の超音波画像を生成してもよい。
【0067】
しかしながら、本装置は、超音波とは異なる波(超音波波長とは異なる波長を有する波)を用いる任意のイメージング装置であってもよく、この場合、トランスデューサ及び電子的ベイ構成要素は上記波に適応化される。
【0068】
図3は、本開示に係る超音波減衰パラメータを推定するための方法であって、図1の装置で実装された方法を示す図である。
【0069】
方法ステップは、主として中央処理装置8によって、最終的には、ディジタル信号プロセッサ10の支援を受けて、又は他の任意の手段によって制御される。本方法は、
(a)上記媒体において少なくとも1つのパルスがトランスデューサによって送信される送信ステップ(101)と、
(b)上記パルスに応答してデータがトランスデューサによって取得される受信ステップ(102)と、
(c)上記処理装置によって上記データが処理されて上記領域の後方散乱取得データを提供する処理ステップ(103)と、
(d)空間-時間領域における深度の関数である、上記後方散乱取得データの自己相関関数が決定される関数決定ステップ(104)と、
(e)上記自己相関関数に基づいて超音波減衰パラメータが推定される減衰推定ステップ(105)とを含む。
【0070】
オプションで、超音波画像は、ステップ103の後方散乱された取得データに基づいて生成されてもよい。
【0071】
例えば、推定された超音波減衰パラメータは、以下の目的で使用されてもよい。
【0072】
・それが大域的推定量を提供する場合、肝臓脂肪症を評価すること、
・それが大域的推定量を提供する場合、デュシェンヌ型筋ジストロフィーキャリアを検出すること、及び/又は
・それが減衰のマップを提供する場合、1つ又は複数の乳房腫瘤を特徴づけること。
【0073】
それらの実施例は、本開示で説明する方法の実証的なアプリケーションとして提供されるが、他のアプリケーション、例えば動物の身体及び/又は材料分析に関連したものもまた同様に関心対象であってもよい。
【0074】
ステップ101~105はオプションでループを形成してもよく、例えば、ステップ105からステップ101に戻ることでループ107を形成してもよい。このことで、複数の減衰パラメータを推定できるようになる可能性があり、オプションで、各サイクルのステップ103において、さらに、超音波画像が生成されてもよい。追加的又は代替的に、ステップ105からステップ104に戻ることで、ステップ104及び105をループするループ108を追加することが可能である。このことは、図6のコンテキストでさらに詳細に説明するように、複数の局所的推定量(オプションで、関心対象領域の異なるエリアにおける推定量)を推定できるようになる可能性がある。
【0075】
いくつかの実施形態では、上記方法は、図4及び図5のコンテキストでさらに詳細に説明するように、大域的超音波減衰パラメータを推定するために使用されてもよく、又は、図6のコンテキストで説明するように、局所的超音波減衰パラメータを推定するために使用されてもよい。しかしながら、本開示の詳細な実施形態について詳しく説明する前に、以下、いくつかの基本について提示する。
【0076】
以下の説明では、次の用語を用いる。
【0077】
ACF:自己相関関数(Auto Correlation Function)
β:超音波減衰(dB/cm/MHz又はNp/cm/MHz)
c:媒体における音速
:パルスの中心周波数
λ:パルスの波長(λ=c/f
MMSE:最小平均二乗誤差(Minimum Mean Square Error)
PSD:パワースペクトル密度(Power Spectral Density)
ROI:関心対象領域(Region Of Interest)
RFデータ:プローブ上で取得されたディジタル化信号を格納する二次元アレイ(トランスデューサ対時間)
σ:ガウス形状スペクトラムの標準偏差。
【0078】
以下では、媒体が均質であり、複数の弱い散乱器を含むと仮定する。送信パルスの解析的信号は、ガウシアン包絡線及び中心周波数fを有すると仮定する。
【0079】
【数1】
【0080】
以下では、本開示で使用されうる一般的な信号モデリングについて説明する。それによって、信号モデリングは、モデルに対して超音波エコー(すなわちPSD)、超音波減衰、回析、反射率を次第に追加することで、次第に向上される。
【0081】
まず、画素において測定されたパワー信号密度(Power Signal Density:PSD)がどのようにモデリングされうるかについて説明する。この説明セクションでは、簡単化のため、超音波減衰は考慮しない。合焦された開口について、時刻tにおいて等時性の体積によって後方散乱された信号は、次式のようにモデリングされてもよい。
【0082】
【数2】
【0083】
ここで、n(t)は、微視的な自己相関関数、すなわち次式を満たす、ゼロ平均のガウス確率過程である。
【0084】
【数3】
【0085】
従って、後方散乱された信号の自己相関関数の表現は、次式のように定式化されてもよい。
【0086】
【数4】
【0087】
ウィーナー・ヒンチン定理を用いると、後方散乱された信号のPSDであるS(f)は、次式のように表すことができる。
【0088】
【数5】
【0089】
ビームフォーミング処理は、関心対象の画素からトランスデューサへの伝搬時間に対応する遅延を有してトランスデューサにおいて取得された複数の信号を蓄積することにある。従って、画素において測定されたPSDは、数式(5)によって与えられる。
【0090】
次に、超音波減衰がどのようにモデリングされうるかについて説明する。特に、超音波減衰は、次式の伝達関数を有するフィルタとしてモデリングされてもよい。
【0091】
【数6】
【0092】
係数2は往復の伝搬を表わす。軸方向の送信及び受信ビームを仮定して、(5)及び(6)の組み合わせは、深度zにおいてビームフォーミングされた画素のPSDの表現を導く。
【0093】
【数7】
【0094】
ここで、次式を用いる。
【0095】
【数8】
【0096】
数式(7)で与えられるPSDにおいてスペクトル領域から時空間領域へのウィーナー・ヒンチン定理を呼び出すことで、次式が得られる。
【0097】
【数9】
【0098】
遅延τ=0において得られたACFは、画素値の平均パワーを与える。数式(8)を用いて次式が得られる。
【0099】
【数10】
【0100】
(9)から(10)への変形は、(9)の右辺がガウス関数の-∞から+∞までの積分であり、従って、その積分値が1になるという事実によって正当化される。
【0101】
次に、画素信全体号がどのようにモデリングされうるかについて説明する。トランスデューサの感度、ビームフォーミングアンテナの利得、回析を考慮し、数式(10)と組み合わせることで、次式の表現を導く。
【0102】
【数11】
【0103】
数式(11)は、その文献で説明された方法の核心を表す。そのフレームワークにおいて、constant(z)が推定されなければならない。
【0104】
数式(11)に基づく超音波減衰推定方法が以下の利点を提示するということが本開示の発明者によって発見された。
【0105】
・以下の理由により、回析効果にはそれほど敏感でない:
-合成ビームフォーミングデータに対して動作すること;従って、任意の空間的場所において、送信及び受信の両方の焦点が保証されること。
-時空間領域において動作し、回析効果はパルスの周波数範囲にわたって平均化されること。
-送信パルスが広帯域であること。
・周波数に基づく方法とは異なり、時間(又は空間的)周波数の妥協を必要としない。
・ビームフォーミングされたデータにおける局所的エネルギー推定は簡単かつ頑健であり、その結果、減衰推定量はより頑健になる。
【0106】
図4を参照して、大域的超音波減衰パラメータを推定するための方法について説明する。従って、図4は、本開示に係る方法の第1の例示的な実施形態(大域的減衰推定)を示す図である。特に、図4の方法は、図3の方法の特定の実施形態であってもよい。
【0107】
図4の方法は、受信ステップb、ここでは参照数字’102及び102’で開始し、このとき、データは、パルスに応答してトランスデューサによって取得される。このデータは、ここでは、検査された媒体に係る「RFデータ集合」と呼ばれ、また、基準媒体、例えばファントムに係る「RFトレーニングデータ集合」と呼ばれる。ステップ102~104は、検査された媒体のデータの処理を示し、一方、同様のステップ102’~104’は、基準媒体のデータの処理を示す。
【0108】
次のステップ103及び103’「合成ビームフォーミング」はそれぞれ、本開示及び図3のステップ(c)、すなわち、処理装置によってデータが処理されて領域の後方散乱された取得データを提供する処理ステップに対応する。
【0109】
次のステップ104、104’「平均パワー推定」はそれぞれ、本開示及び図3のステップ(d)、すなわち、時空間領域における深度の関数である、後方散乱取得データ(特に、ステップ103のビームフォーミングされた取得データ)の自己相関関数が決定される関数決定ステップに対応する。
【0110】
次のステップ104aにおいて、2つの自己相関関数が互いに比較され、後方散乱された取得データにおける回析パターン及び/又は未知のトランスデューサ感度を補償する。特に、両方の関数の間のデルタ関数が決定されてもよい。上記比較ステップの例示的な実施形態を以下で説明する。
【0111】
ステップ105において、2つの自己相関関数の比較、すなわち、特にそれらのデルタ関数を解くことで、超音波減衰パラメータが推定される。より詳しくは、以下の通りである。
【0112】
以下の説明では、次の基本的仮定を用いる。
【0113】
・媒体における減衰が一定であること。
・数式(11)のconstant(z)はzに依存しない、すなわち、それは関心対象領域ROIにわたって一定であること。
【0114】
図3及び図4の方法の目標は、特定のROIにおいてβを推定することである。ROIは、図5に示すような肝臓の均質部分であってもよい。
【0115】
及びzは(それぞれ)、ROIの頂部及び底部の深度を(それぞれ)示すものとする。深度は、プローブの表面と関心対象点(空間的場所)との間の距離を意味する。
【0116】
深度z及びzにおける平均パワーの対数の差を評価することで、次式が得られる。
【0117】
【数12】
【0118】
数式(12)(これはステップ104の結果として理解されてもよい)は、2つの未知数、すなわち、constant(z,z)及びβを含む。
【0119】
βを推定するために、constant(z,z)を推定する必要がある。constant(z,z)はいくつかのファクターに依存するが、媒体には依存せず、従って、それは、既知のβに関連付けられたトレーニングデータ集合(すなわちステップ102’の基準データ)において学習されてもよい。このことは、減衰について較正されたファントム(すなわち基準媒体)において超音波データを取得することで達成可能である。δEref(z,z)は、較正されたファントムにおける深度z及びzで測定された平均パワーの対数の差を示し、βrefは、その較正された(既知の)減衰を示し、ΔE(z,z)は、δE(z,z)及びδEref(z,z)の差を示すものとする。次いで、それは次式で決定可能である。
【0120】
【数13】
【0121】
数式(13)は、比較ステップ104aの例示的な結果であってもよい。数式(13)(すなわちステップ105)の検査は、ΔE(z,z)=0を解くことがβの推定値を与えることを明らかにする。
【0122】
数式(13)は次式で近似されうることに注意する。
【0123】
【数14】
【0124】
このことは、数式(13)の二次の項が比較的小さく、したがって無視できるという状況に起因する。
【0125】
この方法の重要なステップは、深度z及びzにおける平均パワーを推定することにあるということに注意する。その推定処理の第1ステージは、ステップ102の取得された超音波RFデータに対して、ステップ103で合成ビームフォーミングを行うことであってもよく、その目標は、回析効果を最小化することにある。第2ステージは、特定の深度におけるこれらのビームフォーミングされたデータの平均パワーを推定することであってもよい。しかしながら、ビームフォーミングされたデータは、超音波イメージングにおける公知の問題、すなわち、それらがスペックルノイズによって損なわれるという問題の影響を受ける。スペックルノイズは、空間分解能を犠牲にして空間的平均化を行うことにより平滑化されうる。代替として、ビームフォーミングされた取得データを得るためにビームフォーミングを使用することが可能である。この場合、平均パワー推定値の変動を減少させるために、合成ビームフォーミングステージにおいて開口における直交アポディゼーションを使用することが可能になる。
【0126】
本開示の方法、特に、上述した例示的な方法は、次の利点を有する。
【0127】
・方程式ソルバ(方程式解計算器)は、(例えば、比較的小さな二次の項に起因して)数式(13)及び(14)の両方に基づくものであってもよい。
・スペクトル解析は実行されず、従って、時間/空間対周波数分解能のトレードオフを考慮する必要はない。
・回析を最小化するために、広帯域のパルスが使用されてもよい。
・送信及び受信の両方の合焦を保証するために、従って、回析を低減するために、合成ビームフォーミングが実行されてもよい。
・スペックルを平滑化して平均パワー変動推定値を減少させるために、従って、超音波減衰推定値の変動を低減するために、直交アポディゼーションが使用されてもよい。
・回析又はプローブ感度のような、媒体に依存しないパラメータが、較正されたファントムにおいて学習される。
【0128】
図6は、本開示に係る方法の第2の例示的な実施形態(局所的減衰推定)を示す図である。
【0129】
図6の方法は図4のもの、すなわち、ステップ102~104及びステップ102’~104’に対応する。しかしながら、(図4のステップ104a及び105の代わりに)後のステップ106において、局所的減衰パラメータを推定するために、両方の自己相関関数がカルマンフィルタを用いて比較される。
【0130】
図6の方法の目標は、1つ又は複数の局所的減衰パラメータを、例えば局所的推定値のマップの形式で取得することにある。このことは、基本的な仮定が図4の方法のものと比較してわずかに異なることを意味する。特に、それは、ここでは、超音波減衰パラメータが空間的場所に応じて変化するということを考慮しなければならない。
【0131】
ここで、深度方向zにおいて軸線に沿った局所的超音波減衰を推定し、ラインに沿って走査すること、すなわち、Bモードイメージングで通常行われるものと同様の方法で走査することが提案される。
【0132】
従って、この処理の主な目標は、特定の走査線における超音波減衰パラメータを、深度の関数、すなわちβ(z)として推定することにある。このため、数式(11)の逆関数を求めることを提案する。この数式から直接的に逆関数を求めることは容易ではない。したがって、ノイズを含む状態空間系を観測(又は測定)することによって、平均パワーの対数を深度の関数としてモデリングすることが提案される。
【0133】
深度の関数としての平均パワーの対数は離散化されてもよく、次式の線形確率差分方程式(15)によって支配される、離散時間制御された処理の測定としてモデリングされてもよい。
【0134】
【数15】
【0135】
ここで、
・Xkは(隠れた)状態ベクトルである。
・Aは遷移行列である。
・w及びvは処理及び測定のノイズである。
・Hは測定行列である。
・LEは、テスト対象媒体及び深度kΔzにおける既知の減衰を有する基準媒体のインデックスkにおける平均パワーの対数と、深度kΔzにおける基準媒体の平均パワーの対数との差を含むベクトルである。
・Δzは深度サンプリングステップである。
・Bは制御行列である。
・uは駆動入力である。
【0136】
E[k]は、深度kΔzにおいて平均パワー推定値の対数を示すものとする。数式(11)のテーラー展開を行うことにより、次式が得られる。
【0137】
【数16】
【0138】
ここで、次式のように近似されてもよい。
【0139】
【数17】
【0140】
ΔE[k]は、テスト対象媒体及び既知の減衰βrefを有する基準媒体、例えば較正されたファントムのインデックスkにおける平均パワーの対数の差を示すものとする。この場合、次式が得られる。
【0141】
【数18】
【0142】
において(constant(k)-constant(k-1))が含まれる場合、次式が得られる。
【0143】
【数19】
【0144】
ここで、
・Erefは、基準媒体についてインデックスkにおいて推定された平均パワーの対数である。
・n[k]、n[k]、n[k]は、状態ノイズ、すなわち、ゼロ平均の白色ガウス雑音である。
【0145】
次式を満たす場合、数式(17)は、数式(15)の形式にすることができる。
【0146】
【数20】
【0147】
カルマンフィルタが、MMSEの意味で、数式(15)における状態ベクトルxの最適な推定量を提供することがわかった。その結果、β[k]の最適な推定量(MMSEの意味で)を提供する。
【0148】
請求項を含む本開示の全体にわたって、用語「1つの~を備える(含む)」は、他の記載がなければ、「少なくとも1つの~を備える(含む)」と同義のものとして理解されるべきである。さらに、請求項を含む本開示で示した任意の範囲は、他の記載がなければ、その端の値を含むものとして理解されるべきである。説明した構成要素に係る特定の値は、当業者に既知である、許容された製造時又は産業上の公差内にあるものと理解されるべきであり、用語「実質的に」及び/又は「約」及び/又は「概して」のいかなる使用も、そのような許容された公差内にあることを意味するものと理解されるべきである。
【0149】
本開示は特定の実施形態を参照してここに説明したが、これらの実施形態は、本開示の原理及び応用の単なる例示であることが理解されるべきである。
【0150】
明細書及び実施例が例示的なものとしてのみ考慮され、本開示の真の範囲は添付の特許請求の範囲によって示されることが意図される。
【0151】
要約すると、上述した本開示に係る方法は、より正確な減衰推定を可能にし、より少ない計算上のコストを意味し、このことは、特に、リアルタイム計算モードを改善する。さらに、増大した精度に起因して、変動の減少を達成することができ、したがって、増大した再現性を達成することができる。
図1
図2
図3
図4
図5
図6