JP6253869B2 - 超音波診断装置、超音波診断装置の作動方法および超音波診断装置の作動プログラム - Google Patents

Description

本発明は、超音波を用いて観測対象の組織を観測する超音波診断装置、超音波診断装置の作動方法および超音波診断装置の作動プログラムに関する。

従来、超音波を用いた検体等の観測対象の組織性状を観測する技術として、観測対象で後方散乱された超音波エコーを超音波振動子で受信した超音波信号の周波数スペクトルの特徴量を画像化する技術が知られている（例えば、特許文献１を参照）。なお、散乱とは、音波が媒体中で粒子と衝突あるいは相互作用して方向を変えられることである。さらに、後方散乱とは、このうち音源の方向に戻ってくることである。この現象は一般に反射とも言われるが、本願では以下、後方散乱を用いる。このときの音源は超音波振動子である。この技術では、観測対象の組織性状を表す解析値として周波数スペクトルの特徴量を抽出した後、この特徴量に対応する視覚的な情報、例えば色情報を付与した特徴量画像を生成して表示する。医師等の術者は、表示された特徴量画像を見ることによって検体の組織性状を診断する。

国際公開第２０１２／０１１４１４号

ところで、組織性状を鑑別するにあたり、観測対象が有する散乱体の大きさや数密度を把握することが重要である。しかしながら、特許文献１に記載の技術は、周波数スペクトルの特徴量と、観測対象が有する散乱体の直径や数密度といった物理量が対応付いておらず、特徴量に応じて付与される色情報から組織性状を鑑別するには熟練が必要であった。このため、熟練を要さずに、容易に組織性状を鑑別することができる技術が望まれていた。

本発明は、上記に鑑みてなされたものであって、特徴量に基づく組織性状の鑑別を容易かつ正確に行うことができる超音波診断装置、超音波診断装置の作動方法および超音波診断装置の作動プログラムを提供することを目的とする。

上述した課題を解決し、目的を達成するために、本発明に係る超音波診断装置は、観測対象へ超音波を送信し、該観測対象で後方散乱された超音波を受信する超音波振動子を備えた超音波プローブが取得した超音波信号に基づいて超音波画像を生成する超音波診断装置であって、前記観測対象から受信した超音波信号に基づき、特徴量を算出する特徴量算出部と、含有する散乱体の物理量が既知の対照物体の前記物理量、および該対象物体から得られた特徴量に基づき導出された関係と、前記特徴量算出部が算出した前記観測対象の前記特徴量とを用いて、前記観測対象が含む散乱体の物理量を推定する推定部と、前記推定部による推定結果を含み、表示部に表示させるための情報を生成する物理量情報生成部と、を備えたことを特徴とする。

本発明に係る超音波診断装置は、上記発明において、前記対照物体が含有する物理量は、前記散乱体の数密度、前記散乱体の大きさ、および前記散乱体の散乱強度の少なくとも１つを含む、ことを特徴とする。

本発明に係る超音波診断装置は、上記発明において、前記物理量は、前記対照物体が含有する散乱体の散乱強度を含み、前記散乱強度は、前記散乱体と媒体の振幅反射率、エネルギー反射率、および／または、それらの関数であることを特徴とする。

本発明に係る超音波診断装置は、上記発明において、前記推定部は、前記特徴量算出部が算出した前記特徴量を、前記関係としての関係式に代入することにより前記観測対象が含む散乱体の物理量を推定することを特徴とする。

本発明に係る超音波診断装置は、上記発明において、前記関係式、該関係式の係数、該関係式の定数項、および前記関係を記述したテーブルのうちの少なくとも一つを記憶する関係情報記憶部をさらに備え、前記推定部は、前記関係情報記憶部を参照して前記観測対象が含む散乱体の物理量を推定することを特徴とする。

本発明に係る超音波診断装置は、上記発明において、前記関係式は、前記物理量および前記特徴量のうち少なくとも一部を重回帰分析することにより導出されることを特徴とする。

本発明に係る超音波診断装置は、上記発明において、前記既知の散乱体の物理量のうち、前記散乱体の数密度、および前記散乱体の大きさを含む場合、前記散乱体の数密度、および前記散乱体を非線形変換して前記重回帰分析を行うことを特徴とする。

本発明に係る超音波診断装置は、上記発明において、前記特徴量は、前記超音波信号をもとに算出される周波数特徴量を含むことを特徴とする。

本発明に係る超音波診断装置は、上記発明において、前記特徴量は、前記超音波信号をもとに算出される減衰率を含むことを特徴とする。

本発明に係る超音波診断装置は、上記発明において、前記特徴量は、前記超音波信号をもとに算出される音速を含むことを特徴とする。

本発明に係る超音波診断装置は、上記発明において、前記物理量情報生成部は、前記推定部が推定した前記物理量に応じて視覚情報を付与した画像データを生成することを特徴とする。

本発明に係る超音波診断装置は、上記発明において、前記推定部が互いに異なる複数の物理量を推定した場合、前記物理量情報生成部は、前記複数の物理量を同時、順次または異なるタイミングで前記表示部に表示される前記情報を生成することを特徴とする。

本発明に係る超音波診断装置は、上記発明において、前記非線形変換が施された前記散乱体の数密度、および／または前記散乱体の大きさに対してさらに非線形変換を施す変数変換部をさらに備えたことを特徴とする。

本発明に係る超音波診断装置は、上記発明において、前記観測対象から受信した超音波信号に基づき算出した特徴量において減衰の影響を補正する減衰補正部をさらに備え、前記推定部が、前記含有する散乱体の物理量が既知の前記対照物体の前記物理量と、前記対照物体から得られた特徴量において減衰の影響を補正した値とに基づき導出された関係と、前記減衰補正部が補正した前記観測対象の前記特徴量の値とを用いて、前記観測対象が含む散乱体の物理量を推定することを特徴とする。

本発明に係る超音波診断装置の作動方法は、観測対象へ超音波を送信し、該観測対象で後方散乱された超音波を受信する超音波振動子を備えた超音波プローブが取得した超音波信号に基づいて超音波画像を生成する超音波診断装置の作動方法であって、特徴量算出部が、前記観測対象から受信した超音波信号に基づき、特徴量を算出する特徴量算出ステップと、推定部が、含有する散乱体の物理量が既知の対照物体の前記物理量、および該対照物体から得られた特徴量に基づき導出された関係と、前記特徴量算出部が算出した前記観測対象の前記特徴量とを用いて、前記観測対象が含む散乱体の物理量を推定する推定ステップと、物理量情報生成部が、前記推定部による推定結果を含み、表示部に表示させるための情報を生成する物理量情報生成ステップと、を含むことを特徴とする。

本発明に係る超音波診断装置の作動プログラムは、観測対象へ超音波を送信し、該観測対象で後方散乱された超音波を受信する超音波振動子を備えた超音波プローブが取得した超音波信号に基づいて超音波画像を生成する超音波診断装置の作動プログラムであって、特徴量算出部が、前記観測対象から受信した超音波信号に基づき、特徴量を算出する特徴量算出手順と、推定部が、含有する散乱体の物理量が既知の対象物体の前記物理量、および該対照物体から得られた特徴量に基づき導出された関係と、前記特徴量算出部が算出した前記観測対象の前記特徴量とを用いて、前記観測対象が含む散乱体の物理量を推定する推定手順と、物理量情報生成部が、前記推定部による推定結果を含み、表示部に表示させるための情報を生成する物理量情報生成手順と、を前記超音波診断装置に実行させることを特徴とする。

本発明によれば、特徴量に基づく組織性状の鑑別を容易かつ正確に行うことができるという効果を奏する。

図１は、本発明の実施の形態１に係る超音波診断装置を備えた超音波診断システムの構成を示すブロック図である。 図２は、本発明の実施の形態１に係る超音波診断装置の信号増幅部が行う増幅処理における受信深度と増幅率との関係を示す図である。 図３は、超音波振動子の走査領域とＢモード用受信データとを模式的に示す図である。 図４は、超音波信号の１つの音線におけるデータ配列を模式的に示す図である。 図５は、本発明の実施の形態１に係る超音波診断装置の周波数解析部により算出されたスペクトルデータの例を示す図である。 図６は、本発明の実施の形態１に係る超音波診断装置の減衰補正部が補正した補正後特徴量をパラメータとして有する直線を示す図である。 図７は、ミッドバンドフィットと、散乱体の直径と、数密度との関係を示す図であって、回帰平面を説明する図である。 図８は、本発明の実施の形態１に係る超音波診断装置が行う処理の概要を示すフローチャートである。 図９は、本発明の実施の形態１に係る超音波診断装置の周波数解析部が実行する処理の概要を示すフローチャートである。 図１０は、本発明の実施の形態１に係る超音波診断装置の表示装置における特徴量画像の表示例を模式的に示す図である。 図１１は、本発明の実施の形態１の変形例に係る超音波診断装置を備えた超音波診断システムの構成を示すブロック図である。 図１２は、本発明の実施の形態１の変形例に係る超音波診断装置の表示装置における特徴量画像の表示例を模式的に示す図である。 図１３は、本発明の実施の形態３に係る超音波診断装置を備えた超音波診断システムの構成を示すブロック図である。 図１４は、本発明の実施の形態４に係る超音波診断装置を備えた超音波診断システムの構成を示すブロック図である。 図１５は、本発明の実施の形態４に係る超音波診断装置が記憶するルックアップテーブルを説明するための図である。

以下、添付図面を参照して、本発明を実施するための形態（以下、「実施の形態」という）を説明する。

（実施の形態１）
図１は、本発明の実施の形態１に係る超音波診断装置３を備えた超音波診断システム１の構成を示すブロック図である。同図に示す超音波診断システム１は、観測対象である観測対象へ超音波を送信し、該観測対象で後方散乱された超音波を受信する超音波内視鏡２と、超音波内視鏡２が取得した超音波信号に基づいて超音波画像を生成する超音波診断装置３と、超音波診断装置３が生成した超音波画像を表示する表示装置４と、を備える。超音波診断装置３は、超音波内視鏡２を一つ、または超音波内視鏡２の同種、異種に関わらず複数を同時に接続することができる。本実施の形態では、超音波内視鏡２が、超音波プローブとして作用する。なお、図１では、実線の矢印が画像にかかる電気信号の伝送を示し、破線の矢印が制御にかかる電気信号の伝送を示している。

超音波内視鏡２は、その先端部に、超音波診断装置３から受信した電気的なパルス信号を超音波パルス（音響パルス）に変換して観測対象へ照射するとともに、観測対象で後方散乱された超音波エコーを電圧変化で表現する電気的なエコー信号に変換する超音波振動子２１を有する。

超音波内視鏡２は、通常はその先端部に、さらに、撮像光学系および撮像素子を有しており、観測対象の消化管（食道、胃、十二指腸、大腸）、または呼吸器（気管、気管支）へ挿入され、消化管や呼吸器、その周囲臓器（膵臓、胆嚢、胆管、胆道、リンパ節、縦隔臓器、血管等）を撮像することが可能である。超音波内視鏡２は、観測対象への長尺の挿入部を有している。そして、挿入部は、通常は撮像時に観測対象へ照射する照明光を導くライトガイドを有する。このライトガイドは、その先端部が挿入部の先端まで達している一方、基端部が照明光を発生する光源装置に接続されている。

超音波診断装置３は、超音波内視鏡２と電気的に接続され、所定の波形および送信タイミングに基づいて高電圧パルスからなる送信信号（パルス信号）を超音波振動子２１へ送信するとともに、超音波振動子２１から電気的な受信信号であるエコー信号を受信して高周波（ＲＦ：Radio Frequency）信号のデジタルデータ（以下、ＲＦデータという）を生成、出力する送受信部３１と、送受信部３１から受信したＲＦデータをもとにデジタルのＢモード用受信データを生成する信号処理部３２と、送受信部３１から受信したＲＦデータに対して所定の演算を施す演算部３３と、各種画像データを生成する画像処理部３４と、キーボード、マウス、タッチパネル等のユーザインタフェースを用いて実現され、各種情報の入力を受け付ける入力部３５と、超音波診断システム１全体を制御する制御部３６と、超音波診断装置３の動作に必要な各種情報を記憶する記憶部３７と、を備える。

送受信部３１は、エコー信号を増幅する信号増幅部３１１を有する。信号増幅部３１１は、受信深度が大きいエコー信号ほど高い増幅率で増幅するＳＴＣ（Sensitivity Time Control）補正を行う。図２は、信号増幅部３１１が行う増幅処理における受信深度と増幅率との関係を示す図である。図２に示す受信深度ｚは、超音波の受信開始時点からの経過時間に基づいて算出される量である。図２に示すように、増幅率β（ｄＢ）は、受信深度ｚが閾値ｚ_thより小さい場合、受信深度ｚの増加に伴ってβ_０からβ_th（＞β₀）へ線型に増加する。また、増幅率β（ｄＢ）は、受信深度ｚが閾値ｚ_th以上である場合、一定値β_thをとる。閾値ｚ_thの値は、観測対象から受信する超音波信号がほとんど減衰してしまい、ノイズが支配的になるような値である。より一般に、増幅率βは、受信深度ｚが閾値ｚ_thより小さい場合、受信深度ｚの増加に伴って単調増加すればよい。なお、図２に示す関係は、予め記憶部３７に記憶されている。

送受信部３１は、信号増幅部３１１によって増幅されたエコー信号と、増幅されない元のエコー信号との双方に対してフィルタリング等の処理を施した後、適当なサンプリング周波数（例えば５０ＭＨｚ）でサンプリングして離散化（いわゆるＡ／Ｄ変換処理）する。こうして、送受信部３１は、増幅されないエコー信号から離散化されたＲＦデータＡ、増幅されないエコー信号から離散化されたＲＦデータＮの２つのＲＦデータを生成し、信号処理部３２および演算部３３へ出力する。ここで、Ａは、“Ａｍｐｌｉｆｉｅｄ”の略、Ｎは、“Ｎｏｒｍａｌ”の略である。なお、超音波内視鏡２が複数の素子をアレイ状に設けた超音波振動子２１を電子的に走査させる構成を有する場合、送受信部３１は、複数の素子に対応したビーム合成用の多チャンネル回路を有する。

送受信部３１が送信するパルス信号の周波数帯域は、超音波振動子２１がパルス信号を超音波パルスへ電気音響変換をする際の線型応答周波数帯域をほぼカバーする広帯域にする。また、信号増幅部３１１におけるエコー信号の各種処理周波数帯域は、超音波振動子２１が超音波エコーをエコー信号へ音響電気変換する際の線型応答周波数帯域をほぼカバーする広帯域にする。これらにより、後述する周波数スペクトルの近似処理を実行する際、精度のよい近似を行うことが可能となる。

送受信部３１は、制御部３６が出力する各種制御信号を超音波内視鏡２に対して送信するとともに、超音波内視鏡２から識別用のＩＤなどを含む各種情報を受信して制御部３６へ送信する機能も有する。

信号処理部３２は、ＲＦデータＡに対してバンドパスフィルタ、包絡線検波、対数変換など公知の処理を施し、デジタルのＢモード用受信データを生成する。対数変換では、ＲＦデータＡを基準電圧Ｖ_cで除した量の常用対数をとってデシベル値で表現する。このＢモード用受信データでは、超音波パルスの後方散乱の強さを示す受信信号の振幅または強度が、超音波パルスの送受信方向（深度方向）に沿って並んでいる。図３は、超音波振動子２１の走査領域（以下、単に走査領域ということもある）とＢモード用受信データとを模式的に示す図である。図３に示す走査領域Ｓは扇形をなしている。なお、図３では、超音波振動子２１が、超音波が往復する経路（音線）を直線で、Ｂモード用受信データを各音線上に並んだ点で表現している。図３では、後の説明の都合上、各音線に、走査開始（図３右）から順に、１、２、３・・・と番号を付している。そして、１番目の音線をＳＲ₁、２番目の音線をＳＲ₂、３番目の音線をＳＲ₃、・・・、ｋ番目の音線をＳＲ_kと定義する。図３は、超音波振動子２１がコンベックス振動子である場合に相当している。また、図３では、Ｂモード用受信データの受信深度をｚとして記載している。超音波振動子２１の表面から照射された超音波パルスが受信深度ｚにある物体内で後方散乱し、超音波エコーとして超音波振動子２１へ戻ってきた場合、その往復距離Ｌと受信深度ｚとの間には、ｚ＝Ｌ／２の関係がある。信号処理部３２は、生成したＢモード用受信データを、画像処理部３４のＢモード画像データ生成部３４１へ出力する。信号処理部３２は、ＣＰＵ（Central Processing Unit）等の汎用プロセッサ、またはＡＳＩＣ（Application Specific Integrated Circuit）もしくはＦＰＧＡ（Field Programmable Gate Array）等の特定の機能を実行する専用の集積回路等を用いて実現される。

演算部３３は、送受信部３１が生成したＲＦデータＮに高速フーリエ変換（ＦＦＴ：Fast Fourier Transform）を施して周波数解析を行うことによりスペクトルデータを算出する周波数解析部３３１と、周波数解析部３３１により算出されたスペクトルデータを用いて、単回帰分析によりスペクトルデータの特徴量を算出する単回帰分析部３３２（特徴量算出部）と、単回帰分析部３３２が算出した特徴量、および記憶部３７に記憶されている関係式であって、散乱体の数密度（以下、数密度をｎとする）および直径（以下、散乱体の直径をｄとする）から導かれる関係式を用いて、数密度ｎの対数（logｎ）および直径ｄの対数（logｄ）を推定する推定部３３３と、を有する。演算部３３は、ＣＰＵ（Central Processing Unit）や各種演算回路等を用いて実現される。なお、散乱体の数密度とは、単位体積あたりに含まれる散乱体の個数を意味し、本実施の形態での単位は［個／ｃｍ³］である。また、直径ｄの本実施の形態での単位は［μｍ］である。また、散乱体の数密度と直径の対数は常用対数とし、底を１０とする。従って、logｎやlogｄは、ｎやｄの１０進表現の桁数−１を表現することになる。

周波数解析部３３１は、送受信部３１が生成した各音線のＲＦデータＮ（ラインデータ）を所定の時間間隔で再びサンプリングし、サンプルデータを生成する。周波数解析部３３１は、サンプルデータ群にＦＦＴ処理を施すことにより、ＲＦデータ上の多数の箇所（データ位置）における周波数スペクトルを算出する。ここでいう「周波数スペクトル」とは、サンプルデータ群をＦＦＴ処理を施すことによって得られた「ある受信深度ｚにおける強度の周波数分布」を意味する。また、ここでいう「強度」とは、例えばエコー信号の電圧、エコー信号の電力、超音波エコーの音圧、超音波エコーの音響エネルギー等のパラメータ、これらパラメータの振幅や時間積分値やその組み合わせのいずれかを指す。

本実施の形態では、強度としてエコー信号の電圧を採用し、周波数解析部３３１が、電圧振幅の周波数成分Ｖ（ｆ，Ｌ）をもとに周波数スペクトルのデータ（以下、スペクトルデータともいう）を生成するものとして説明する。ｆは、周波数である。周波数解析部３３１は、電圧振幅の周波数成分Ｖ（ｆ，Ｌ）を基準電圧Ｖ_cで除し、常用対数（ｌｏｇ）をとってデシベル単位で表現する対数変換処理を施した後、適当な正の定数Ａを乗ずることにより、次式（１）で与えられるスペクトルデータＦ（ｆ，Ｌ）を生成する。
Ｆ（ｆ，Ｌ）＝Ａ・ｌｏｇ｛Ｖ（ｆ，Ｌ）／Ｖ_c｝ ・・・（１）
ここで、ｌｏｇは常用対数である（以下、同じ）。

以下、周波数解析部３３１での周波数解析により電圧振幅の周波数成分Ｖ（ｆ，Ｌ）を求める方法について説明する。一般に、エコー信号の周波数スペクトルは、観測対象が生体組織である場合、超音波が走査された生体組織の性状によって異なる傾向を示す。これは、周波数スペクトルが、超音波を散乱する散乱体の大きさ、数密度、音響インピーダンス等と相関を有しているためである。ここでいう「生体組織の性状」とは、例えば悪性腫瘍（癌）、良性腫瘍、内分泌腫瘍、粘液性腫瘍、正常組織、嚢胞、脈管などのことである。

図４は、超音波信号の１つの音線ＳＲ_kにおけるデータ配列を模式的に示す図である。音線ＳＲ_kにおける白または黒の長方形は、１つのサンプル点におけるデータを意味している。また、音線ＳＲ_kにおいて、右側に位置するデータほど、超音波振動子２１から音線ＳＲ_kに沿って計った場合の深い箇所からのサンプルデータである（図４の矢印を参照）。音線ＳＲ_kは、前述の通り、送受信部３１でのＡ／Ｄ変換処理によりエコー信号からサンプリングされ、離散化されたＲＦデータを、さらに、周波数解析部３３１によりサンプリングされたサンプルデータである。図４では、番号ｋの音線ＳＲ_kの８番目のデータ位置を受信深度ｚの方向の初期値Ｚ^(k) ₀として設定した場合を示しているが、初期値の位置は任意に設定することができる。周波数解析部３３１による算出結果は複素数で得られ、記憶部３７に格納される。

図４に示すデータ群Ｆ_j（ｊ＝１、２、・・・、Ｋ）は、ＦＦＴ処理の対象となるサンプルデータ群である。一般に、ＦＦＴ処理を行うためには、サンプルデータ群が２のべき乗のデータ数を有している必要がある。この意味で、サンプルデータ群Ｆ_j（ｊ＝１、２、・・・、Ｋ−１）はデータ数が１６（＝２⁴）で正常なデータ群である一方、サンプルデータ群Ｆ_Kは、データ数が１２であるため異常なデータ群である。異常なデータ群に対してＦＦＴ処理を行う際には、不足分だけゼロデータを挿入することにより、正常なサンプルデータ群を生成する処理を行う。この点については、周波数解析部３３１の処理を説明する際に詳述する（図９を参照）。この後、周波数解析部３３１は、前述の通り、ＦＦＴ処理を実行し、電圧振幅の周波数成分Ｖ（ｆ，Ｌ）を算出し、前述の式（１）に基づいてスペクトルデータＦ（ｆ，Ｌ）を算出し、単回帰分析部３３２へ出力する。

図５は、周波数解析部３３１により算出されたスペクトルデータの例を示す図である。図５では、横軸が周波数ｆである。また、図５では、縦軸が、上式（１）で与えられるスペクトルデータＦ（ｆ，Ｌ）である。図５に示す直線Ｌ₁₀については後述する。なお、本実施の形態において、曲線および直線は、離散的な点の集合からなる。

図５に示すスペクトルデータＣ₁において、以後の演算に使用する周波数帯域の下限周波数ｆ_Lおよび上限周波数ｆ_Hは、超音波振動子２１の周波数帯域、送受信部３１が送信するパルス信号の周波数帯域などをもとに決定されるパラメータである。以下、図５において、下限周波数ｆ_Lおよび上限周波数ｆ_Hによって定まる周波数帯域を「周波数帯域Ｕ」という。

単回帰分析部３３２は、周波数解析部３３１から出力された複数のスペクトルデータを直線で近似することによってスペクトルデータの特徴量（以下、補正前特徴量という）を算出する近似部３３２ａと、近似部３３２ａが算出した補正前特徴量に対して周波数に依存した減衰を補正することによって特徴量を算出する減衰補正部３３２ｂと、を有する。

近似部３３２ａは、所定周波数帯域におけるスペクトルデータの単回帰分析を行ってスペクトルデータを一次式（回帰直線）で近似することにより、この近似した一次式を特徴付ける補正前特徴量を算出する。単回帰分析とは、独立変数が１種類のみの場合の回帰分析である。本実施の形態での単回帰分析の独立変数は周波数ｆにあたる。例えば、スペクトルデータが図５に示すＣ₁の状態である場合、近似部３３２ａは、周波数帯域Ｕで単回帰分析を行いスペクトルデータＣ₁の回帰直線である直線Ｌ₁₀を得る。次に、近似部３３２ａは、直線Ｌ₁₀の傾きａ₀、切片ｂ₀、および周波数帯域Ｕの中心周波数（すなわち、「ミッドバンド」）ｆ_M＝（ｆ_L＋ｆ_H）／２の回帰直線上の値であるミッドバンドフィット（Mid-band fit）ｃ₀＝ａ₀ｆ_M＋ｂ₀を補正前特徴量として算出する。このように直線Ｌ₁₀を特徴付ける一次式のパラメータ（傾きａ₀、切片ｂ₀、ミッドバンドフィットｃ₀）でスペクトルデータＣ₁を表現することで、スペクトルデータＣ₁を一次式に近似したことになる。

３つの補正前特徴量のうち、傾きａ_0、切片ｂ₀は、超音波を散乱する散乱体の大きさ、散乱体の散乱強度、散乱体の数密度（濃度）等と相関を有していると考えられる。ミッドバンドフィットｃ₀は、有効な周波数帯域内の中心におけるスペクトルの強度を与える。このため、ミッドバンドフィットｃ₀は、散乱体の大きさ、散乱体の散乱強度、散乱体の数密度に加えて、Ｂモード画像の輝度とある程度の相関を有していると考えられる。この後、近似部３３２ａは、これら補正前特徴量ａ₀、ｂ₀、ｃ₀を減衰補正部３３２ｂへ出力する。なお、近似部３３２ａは、回帰分析によって二次以上の多項式でスペクトルデータを近似するようにしてもよい。

減衰補正部３３２ｂが行う補正について説明する。一般に、超音波の振幅は指数的に減衰する。従って、振幅を対数変換し、デシベル表現にした場合、超音波が受信深度０と受信深度ｚとの間を往復する間に生じる減衰量Ａ（ｆ，ｚ）は、往復する前後の線形の変化（デシベル表現での差）として表現できる。この振幅の減衰量Ａ（ｆ，ｚ）は、観測対象が生体である場合には周波数に依存し、高周波では減衰が大きく、低周波では減衰が小さいことが知られている。特に、一様な組織内では周波数に比例することが経験的に知られており、以下の式（２）で表現される。
Ａ（ｆ，ｚ）＝２ζｚｆ ・・・（２）
ここで、比例定数ζは減衰率と呼ばれる量である。また、ｚは超音波の受信深度であり、ｆは周波数である。減衰率ζの具体的な値は、観測対象が生体である場合、生体の部位や組織に応じて定まる。正常肝では概ね、０．５５ｄＢ／ｃｍ／ＭＨｚである。なお、本実施の形態において、減衰率ζの値は記憶部３７に予め記憶されており、減衰補正部３３２ｂは適宜、記憶部３７から減衰率ζの値を読み出して用いる。入力部３５が、超音波内視鏡２による超音波の送信の前に、予め、観測対象の部位名や組織名の入力を術者から受けた場合には、減衰補正部３３２ｂは、部位名や組織名に対応した減衰率ζの適当な値を読み出し、以下の減衰補正に用いる。さらに、入力部３５が、減衰率ζの値を術者から直接受けた場合には、減衰補正部３３２ｂは、その値を以下の減衰補正に用いる。入力部３５が、一切の入力を術者から受けなかった場合には、減衰補正部３３２ｂは、上記０．５５ｄＢ／ｃｍ／ＭＨｚを以下の減衰補正に用いる。

減衰補正部３３２ｂは、近似部３３２ａが抽出した補正前特徴量（傾きａ₀、切片ｂ₀、ミッドバンドフィットｃ₀）に対し、以下に示す式（３）〜（５）にしたがって減衰補正を行うことにより、補正後特徴量ａ、ｂ、ｃ（以下、単に特徴量と呼ぶ）を算出する。
ａ＝ａ₀＋２ζｚ ・・・（３）
ｂ＝ｂ₀ ・・・（４）
ｃ＝ｃ₀＋Ａ（ｆ_M，ｚ）＝ｃ₀＋２ζｚｆ_M（＝ａｆ_M＋ｂ） ・・・（５）
式（３）、（５）からも明らかなように、減衰補正部３３２ｂは、超音波の受信深度ｚが大きいほど、補正量が大きい補正を行う。また、式（４）によれば、切片に関する補正は恒等変換である。これは、切片が周波数０（Ｈｚ）に対応する周波数成分であって減衰の影響を受けないためである。

図６は、減衰補正部３３２ｂが算出した特徴量ａ、ｂ、ｃをパラメータとして有する直線を示す図である。Ｉを、図６の縦軸の値とすると、直線Ｌ₁の式は、
Ｉ＝ａｆ＋ｂ＝（ａ₀＋２ζｚ）ｆ＋ｂ₀ ・・・（６）
で表される。この式（６）からも明らかなように、直線Ｌ₁は、減衰補正前の直線Ｌ₁₀と比較して、傾きが大きく（ａ＞ａ₀）、かつ切片が同じ（ｂ＝ｂ₀）である。この後、単回帰分析部３３２は、これら減衰補正された特徴量ａ、ｂ、ｃを推定部３３３へ出力する。

推定部３３３は、単回帰分析部３３３が減衰補正を行った特徴量ａおよび特徴量ｃと、記憶部３７に記憶されている定数α、β、γ、α´、β´、γ´とを用いて、数密度ｎおよび直径ｄ、実際にはそれらのオーダーである常用対数logｎおよびlogｄを推定する。具体的には、特徴量ａ（傾き）および特徴量ｃ（ミッドバンドフィット）が、散乱体の数密度ｎや、散乱体の直径ｄのおよその値であるオーダー（桁数−１）、即ち、logｎおよびlogｄと線形関係を有するものと仮定し、推定部３３３は、以下の回帰平面の方程式（７）、（８）から導いた式（１０）をもとに、logｎおよびlogｄを推定する。推定部３３３は、推定したlogｎおよび／またはlogｄを物理量として画像処理部３４に出力する。

以下、推定部３３３が観測対象から未知のlogｎおよびlogｄを推定するための準備を説明する。まず、対照物体を用意する。本実施の形態では、人工的に作成した複数の基準ファントムを対照物体の例として説明する。基準ファントムは、散乱体として多数個の粒子を一様な媒体に混ぜ、固化させることによって形成される。散乱体は、予め、質量密度と音速と音響インピーダンスとが散乱体同士で等しくなるよう材質を同じに揃えられる。直径も等しくなるようふるいにかけられ、選別される。さらに、偏りがでないよう媒体に均一に混ぜられ、基準ファントム内において数密度も一様である。このように、基準ファントムにおける散乱体の材質、質量密度、音速、音響インピーダンス、直径、数密度は既知の値であり、基準ファントム内で一様である。一方、媒体の材質、質量密度、音速、音響インピーダンスもやはり既知の値であり、基準ファントム内で一様である。このように、基準ファントムにおいては、散乱体、媒体が均一、一様であるため、減衰率［ｄＢ／ｃｍ／ＭＨｚ］も一様である。そして、これらのパラメータを変えて、複数の基準ファントムを作成する。

次に、散乱体の数密度と散乱体の直径とが異なる組み合わせ（ｎ₁，ｄ₁）、（ｎ₂，ｄ₂）、（ｎ₃，ｄ₃）、・・・、（ｎ_N，ｄ_N）であることがわかっているＮ個の基準ファントムに超音波を送信して得られた超音波エコーに基づき特徴量ｃ₁，ｃ₂，ｃ₃，・・・，ｃ_N（ミッドバンドフィット）を得る。以下の説明では、Ｎ個のファントムのうち、ｉ番目の基準ファントムの散乱体の数密度をｎ_i、直径をｄ_i、減衰率ζ_iとし、このｉ番目の基準ファントムから得られた特徴量（ミッドバンドフィット）をｃ_iとする（１≦ｉ≦Ｎ）。この際、特徴量ｃ_iは、ｉ番目の基準ファントムの減衰率ζ_iを基に、ｉ番目の基準ファントムのスペクトルデータが算出された深度および式（５）を用いて減衰の影響が排除され、深度に依存しない。こうして、これら散乱体の数密度、散乱体の直径、特徴量からなるＮ個のデータセット（ｃ₁，ｎ₁，ｄ₁）、（ｃ₂，ｎ₂，ｄ₂）、（ｃ₃，ｎ₃，ｄ₃）、・・・、（ｃ_N，ｎ_N，ｄ_N）を得る。

次に、Ｎ個のデータセットから、散乱体の数密度と、散乱体の直径と、特徴量との関係を求める。図７は、ミッドバンドフィット（特徴量ｃ）と、散乱体の数密度の常用対数（logｎ）と、散乱体の直径の常用対数（logｄ）とを互いに直交する各軸にとった直交座標系を示している。本実施の形態１では、予め、１番目の基準ファントムに対応するプロットＰ₁（logｎ₁，logｄ₁，ｃ₁）、２番目の基準ファントムに対応するプロットＰ₂（logｎ₂，logｄ₂，ｃ₂）、３番目の基準ファントムに対応するプロットＰ₃（logｎ₃，logｄ₃，ｃ₃）、・・・、Ｎ番目の基準ファントムに対応するプロットＰ_N（logｎ_N，logｄ_N，ｃ_N）の回帰平面ＰＬを、最小二乗法を用いて求める。特徴量ｃに対するプロットＰ₁、Ｐ₂、Ｐ₃、・・・、Ｐ_Nの回帰平面ＰＬとは、プロットＰ₁、Ｐ₂、Ｐ₃、・・・、Ｐ_Nとの差分（誤差）、すなわち、線分Ｐ₁Ｐ₁＾、Ｐ₂Ｐ₂＾、Ｐ₃Ｐ₃＾、・・・、Ｐ_NＰ_N＾の距離の二乗和が最小となる平面である。ここで、Ｐ_i＾は、プロットＰ_iを通り、関心のある軸（ここでは、特徴量ｃの座標軸）に平行な直線とこの平面との交点である。このように、最小二乗法を用いて回帰平面ＰＬを求め、これに基づき分析することは重回帰分析に該当する。重回帰分析とは、独立変数が２種類以上の場合の回帰分析である。本実施の形態での重回帰分析の独立変数は散乱体の直径ｄと散乱体の数密度ｎにあたる。

以下、回帰平面ＰＬを具体的に求める方法を述べる。一般的に、平面の方程式は一次方程式で書くことができる。図７に示した特徴量ｃに対する回帰平面ＰＬは、三次元空間内にあるため、特徴量ｃを２変数（logｎ，logｄ）の一次結合で表現した方程式（式（７））として書くことができる。ここで、α、β、γは平面を決定する実数の定数である。特徴量ｃに対する回帰平面ＰＬを求めることと、α、β、γの値を求めることとは等価である。特徴量ａについても、特徴量ｃと同様にして、回帰平面ＰＬを考えることができ、特徴量ａを２変数（logｎ，logｄ）の一次結合で表現した方程式（式（８））として書くことができる。ここで、α´、β´、γ´は平面を決定する実数の定数である。特徴量ａに対する回帰平面ＰＬを求めることと、α´、β´、γ´の値を求めることとは等価である。

上式（７）、（８）から、下式（９）が得られる。

式（９）より下式（１０）が得られる。

ここで、定数α、β、γ、α´、β´、γ´の求め方について説明する。以下では、一例として定数α、β、γの求め方について説明する。Ｎ個の基準ファントムから特徴量ｃ_iをそれぞれ取得後、散乱体の数密度ｎ_i、直径ｄ_i、特徴量ｃ_iの全基準ファントムにわたる平均をそれぞれ求める。その後、下式（１１）、（１２）に示す行列Ｇおよび列ベクトルＹを定義する。回帰平面の方程式の係数αとβとは、下式（１１）、下式（１２）で定義される下式（１３）によって得ることができる。

ここで、左肩のｔは転置行列を示し、右肩の−１は逆行列を示している。

また、回帰平面ＰＬが、散乱体の数密度ｎ_i、直径ｄ_i、特徴量ｃ_iの各平均を座標とする点（Ｐ₁、Ｐ₂、・・・、Ｐ_Nの重心）を通過することから、式（７）より下式（１４）を得る。上式（１３）で得たα、βを下式（１４）へ代入して定数γを求めることができる。

上述したようにして、Ｎ個の基準ファントムに関するデータセット（ｃ₁，ｎ₁，ｄ₁）、（ｃ₂，ｎ₂，ｄ₂）、（ｃ₃，ｎ₃，ｄ₃）、・・・、（ｃ_N，ｎ_N，ｄ_N）を用いて、上式（７）の定数α、β、γを求めることができる。上式（８）の定数α´、β´、γ´についても同様に、Ｎ個の基準ファントムに関するデータセット（ａ₁，ｎ₁，ｄ₁）、（ａ₂，ｎ₂，ｄ₂）、（ａ₃，ｎ₃，ｄ₃）、・・・、（ａ_N，ｎ_N，ｄ_N）を用いて、定数α´、β´、γ´を求めることができる。

以上、定数α、β、γ、α´、β´、γ´が求められたことにより、推定部３３３が観測対象から未知のlogｎおよびlogｄを推定するための準備を説明できた。ここで、超音波内視鏡２を介して得られた観測対象の特徴量ｃ、特徴量ａ、および、観測対象内の散乱体の直径ｄの桁数−１であるlogｄ、散乱体の数密度ｎの桁数−１であるlogｎも、基準ファントムから得られた物理的な傾向に従うと考えられる。そのため、観測対象からのデータセット（ｎ，ｄ，ｃ）を図７にプロットすると、回帰平面ＰＬ上のどこかに推定されるはずである。よって、このデータセット（ｎ，ｄ，ｃ）は、回帰平面ＰＬの方程式（７）を満足する。同様に、観測対象からのデータセット（ｎ，ｄ，ａ）も式（８）を満足する。観測対象のｎ，ｄ，ｃ，ａは、式（７）と式（８）とを満足するので、式（１０）を満足する。よって、推定部３３３は、このように観測対象からの特徴量ｃ、特徴量ａを式（１０）へ代入し、logｎおよびlogｄを推定する。

なお、本実施の形態では、上述した基準ファントムからの超音波エコーに基づいたエコー信号の受信、特徴量ｃ、特徴量ａの算出、Ｎ個のデータセットの算出、定数α、β、γ、α´、β´、γ´の算出は、超音波診断装置３の外部で実施される。そして、取得した定数α、β、γ、α´、β´、γ´は、工場出荷前に、入力部３５を介して記憶部３７に内蔵された関係情報記憶部３７１に記憶される。本実施の形態で重要な点は、推定部３３３が、物理量が既知である対照物体からのエコー信号に基づく特徴量を利用して、重回帰分析の手法に基づき、logｎおよびlogｄを推定したことである。

図１に戻り、画像処理部３４は、エコー信号の振幅を輝度で表現するＢモード画像データを生成するＢモード画像データ生成部３４１と、Ｂモード画像データ生成部３４１が生成したＢモード画像データと、推定部３３３が推定した物理量（logｎおよび／またはlogｄ）を視覚情報と関連づけてＢモード画像とともに表示する物理量画像データを生成する物理量画像データ生成部３４２（物理量情報生成部）と、を有する。

Ｂモード画像データ生成部３４１は、信号処理部３２から受信したＢモード用受信データに対してゲイン処理、コントラスト処理等の公知の技術を用いた信号処理を行うとともに、表示装置４における画像の表示レンジに応じて定まるデータステップ幅に応じたデータの間引き等を行うことによってＢモード画像データを生成する。Ｂモード画像は、色空間としてＲＧＢ表色系を採用した場合の変数であるＲ（赤）、Ｇ（緑）、Ｂ（青）の値を一致させたグレースケール画像である。

Ｂモード画像データ生成部３４１は、信号処理部３２からのＢモード用受信データに走査範囲を空間的に正しく表現できるよう並べ直す座標変換を施した後、Ｂモード用受信データ間の補間処理を施すことによってＢモード用受信データ間の空隙を埋め、Ｂモード画像データを生成する。Ｂモード画像データ生成部３４１は、生成したＢモード画像データを物理量画像データ生成部３４２へ出力する。

物理量画像データ生成部３４２は、推定部３３３が推定した物理量に関連する視覚情報をＢモード画像データにおける画像の各画素に対して重畳することによって物理量画像データを生成する。物理量画像データ生成部３４２は、例えば図４に示す１つのサンプルデータ群Ｆ_j（ｊ＝１、２、・・・、Ｋ）のデータ量に対応する画素領域に対し、そのサンプルデータ群Ｆ_jから算出される周波数スペクトルの特徴量に対応する物理量に関連する視覚情報を割り当てる。物理量画像データ生成部３４２は、logｎおよびlogｄのいずれか一つに視覚情報としての色相を対応付けることによって物理量画像データを生成する。なお、物理量画像データ生成部３４２が、logｎおよび／またはlogｄの一方に色相を対応付けるとともに、他方に明暗を対応付けることによって物理量画像データを生成するようにしてもよい。物理量に関連する視覚情報としては、例えば色相、彩度、明度、輝度値、Ｒ（赤）、Ｇ（緑）、Ｂ（青）などの所定の表色系を構成する色空間の変数を挙げることができる。

ここで、物理量画像データ生成部３４２が生成する物理量画像データは、図３に示す走査領域Ｓにおいて、特定の深度幅および音線幅などで区切られる関心領域（Region of Interest：ＲＯＩ）に応じた領域の物理量画像が表示装置４に表示されるような物理量画像データを生成する。

制御部３６は、演算および制御機能を有するＣＰＵ等の汎用プロセッサ、またはＡＳＩＣもしくはＦＰＧＡ等の専用の集積回路等を用いて実現される。制御部３６は、記憶部３７が記憶、格納する情報を記憶部３７から読み出し、超音波診断装置３の作動方法に関連した各種演算処理を実行することによって超音波診断装置３を統括して制御する。なお、制御部３６を信号処理部３２および演算部３３と共通の汎用プロセッサまたは専用の集積回路等を用いて構成することも可能である。

記憶部３７は、減衰補正部３３２ｂが周波数スペクトルごとに算出した複数の特徴量や、画像処理部３４が生成した画像データを記憶する。また、記憶部３７は、推定部３３３が推定処理を行なう際の関係式（定数α、β、γ、α´、β´、γ´）を記憶する関係情報記憶部３７１を有する。

記憶部３７は、上記以外にも、例えば増幅処理に必要な情報（図２に示す増幅率と受信深度との関係）、対数変換処理に必要な情報（式（１）参照、例えばＡ、Ｖ_cの値）、周波数解析処理に必要な窓関数（Hamming、Hanning、Blackman等）の情報等を記憶する。

また、記憶部３７は、超音波診断装置３の作動方法を実行するための作動プログラムを含む各種プログラムを記憶する。作動プログラムは、ハードディスク、フラッシュメモリ、ＣＤ−ＲＯＭ、ＤＶＤ−ＲＯＭ、フレキシブルディスク等のコンピュータ読み取り可能な記録媒体に記録して広く流通させることも可能である。なお、上述した各種プログラムは、通信ネットワークを介してダウンロードすることによって取得することも可能である。ここでいう通信ネットワークは、例えば既存の公衆回線網、ＬＡＮ（Local Area Network）、ＷＡＮ（Wide Area Network）などによって実現されるものであり、有線、無線を問わない。

以上の構成を有する記憶部３７は、各種プログラム等が予めインストールされたＲＯＭ（Read Only Memory）、および各処理の演算パラメータやデータ等を記憶するＲＡＭ（Random Access Memory）やハードディスク等を用いて実現される。

図８は、以上の構成を有する超音波診断装置３が行う処理の概要を示すフローチャートである。まず、超音波診断装置３は、超音波内視鏡２から超音波振動子２１による観測対象の測定結果としてのエコー信号を受信する（ステップＳ１）。

超音波振動子２１からエコー信号を受信した信号増幅部３１１は、そのエコー信号の増幅を行う（ステップＳ２）。ここで、信号増幅部３１１は、例えば図２に示す増幅率と受信深度との関係に基づいてエコー信号の増幅（ＳＴＣ補正）を行う。次に、送受信部３１は、適当なサンプリング周波数（例えば５０ＭＨｚ）で増幅されたエコー信号と増幅されない元のエコー信号との双方をサンプリングして離散化し、それぞれからＲＦデータＡ、ＲＦデータＮを生成して、前者をＢモード画像データ生成部３４１へ、後者を周波数解析部３３１へ出力する。

続いて、Ｂモード画像データ生成部３４１は、送受信部３１から出力されたＲＦデータＡを用いてＢモード画像データを生成し、物理量画像データ生成部３４２へ出力する（ステップＳ３）。物理量画像データ生成部３４２はＢモード画像データには処理を施さず、そのまま、表示装置４へ出力する。Ｂモード画像データを受信した表示装置４は、そのＢモード画像データに対応するＢモード画像を表示する（ステップＳ４）。

この後、制御部３６は、術者等のユーザから入力部３５の図示しないボタンもしくはメニューを介して、物理量画像の「表示」もしくは「非表示」のどちらが選択されているのか確認する（ステップＳ５）。制御部３６は、「表示」の選択を確認した場合には演算部３３を構成する各部へ物理量画像作成開始命令を出力する（ステップＳ５：Ｙｅｓ）。「非表示」の選択を確認した場合は物理量画像作成開始命令を出さない（ステップＳ５：Ｎｏ）。演算部３３の各部は、物理量画像作成開始命令を受信すると、後述のステップＳ６以降の処理を実行する。なお、物理量画像作成開始命令の有無に関わらず、超音波診断装置３の送受信部３１、信号増幅部３１１、信号処理部３２、Ｂモード画像データ生成部３４１、物理量画像データ生成部３４２は上記ステップＳ１からＳ４までの処理を繰り返す。そのため、ユーザが入力部３５へ物理量画像の「非表示」を指示している間は、Ｂモード画像が超音波振動子２１による観測対象内の走査のたびに繰り返し表示装置４に表示される。

演算部３３の各部が物理量画像作成開始命令を受信した場合、まず、周波数解析部３３１は、ＲＦデータＮにＦＦＴ処理による周波数解析を行うことによって全てのサンプルデータ群に対するスペクトルデータを算出する（ステップＳ６：周波数解析ステップ）。図９は、ステップＳ６において周波数解析部３３１が実行する処理の概要を示すフローチャートである。以下、図９に示すフローチャートを参照して、周波数解析処理を詳細に説明する。

まず、周波数解析部３３１は、解析対象の音線を識別するカウンタｋをｋ₀とする（ステップＳ２１）。この初期値ｋ₀は、図３中、解析範囲の最右の音線の番号である。

続いて、周波数解析部３３１は、ＦＦＴ処理用に取得する一連のデータ群（サンプルデータ群）を代表するデータ位置（受信深度に相当）Ｚ^(k)の初期値Ｚ^(k) ₀を設定する（ステップＳ２２）。例えば、図４では、上述したように、音線ＳＲ_kの８番目のデータ位置を初期値Ｚ^(k) ₀として設定した場合を示している。この初期値Ｚ^(k) ₀は、音線ＳＲ_k上での解析範囲の最浅の受信深度である。

その後、周波数解析部３３１は、サンプルデータ群を取得し（ステップＳ２３）、取得したサンプルデータ群に対し、記憶部３７が記憶する窓関数を作用させる（ステップＳ２４）。このようにサンプルデータ群に対して窓関数を作用させることにより、サンプルデータ群が境界で不連続になることを回避し、アーチファクトが発生するのを防止することができる。

続いて、周波数解析部３３１は、データ位置Ｚ^(k)のサンプルデータ群が正常なデータ群であるか否かを判定する（ステップＳ２５）。図４を参照した際に説明したように、サンプルデータ群は、２のべき乗のデータ数を有している必要がある。以下、正常なサンプルデータ群のデータ数を２ⁿ（ｎは正の整数）とする。本実施の形態では、データ位置Ｚ^(k)が、できるだけＺ^(k)が属するサンプルデータ群の中心になるよう設定される。具体的には、サンプルデータ群のデータ数は２ⁿであるので、Ｚ^(k)はそのサンプルデータ群の中心に近い２ⁿ／２（＝２^n-1）番目の位置に設定される。この場合、サンプルデータ群が正常であるとは、データ位置Ｚ^(k)より浅い側に２^n-1−１（＝Ｎとする）個のデータがあり、データ位置Ｚ^(k)より深い側に２^n-1（＝Ｍとする）個のデータがあることを意味する。図４に示す場合、サンプルデータ群Ｆ₁、Ｆ₂、Ｆ₃、・・・、Ｆ_K-1はともに正常である。なお、図４ではｎ＝４（Ｎ＝７，Ｍ＝８）の場合を例示している。

ステップＳ２５における判定の結果、データ位置Ｚ^(k)のサンプルデータ群が正常である場合（ステップＳ２５：Ｙｅｓ）、周波数解析部３３１は、後述するステップＳ２７へ移行する。

ステップＳ２５における判定の結果、データ位置Ｚ^(k)のサンプルデータ群が正常でない場合（ステップＳ２５：Ｎｏ）、周波数解析部３３１は、不足分だけゼロデータを挿入することによって正常なサンプルデータ群を生成する（ステップＳ２６）。ステップＳ２５において正常でないと判定されたサンプルデータ群（例えば図４のサンプルデータ群Ｆ_K）は、ゼロデータを追加する前に窓関数が作用されている。このため、サンプルデータ群にゼロデータを挿入してもデータの不連続は生じない。ステップＳ２６の後、周波数解析部３３１は、後述するステップＳ２７へ移行する。

ステップＳ２７において、周波数解析部３３１は、サンプルデータ群を用いてＦＦＴ演算を行うことにより、振幅の周波数成分であるＶ（ｆ，Ｌ）を算出する。その後、周波数解析部３３１は、Ｖ(ｆ，Ｌ）に対数変換処理を施して、スペクトルデータＦ（ｆ，Ｌ）を得る（ステップＳ２７）。

続いて、周波数解析部３３１は、データ位置Ｚ^(k)をステップ幅Ｄで変化させる（ステップＳ２８）。ステップ幅Ｄについて、入力部３５を経由した術者の入力値を記憶部３７が予め記憶しているものとする。図４では、Ｄ＝１５の場合を例示している。ステップ幅Ｄは、できるだけ小さく、特に、Ｂモード画像データ生成部３４１がＢモード画像データを生成する際に利用するデータステップ幅と一致させることが望ましいが、周波数解析部３３１における演算量を削減したい場合には、ステップ幅Ｄとしてデータステップ幅より大きい値を設定してもよい。

その後、周波数解析部３３１は、データ位置Ｚ^(k)が音線ＳＲ_kにおける最大値Ｚ^(k) _maxより大きいか否かを判定する（ステップＳ２９）。この最大値Ｚ^(k) _maxは、音線ＳＲ_k上での解析範囲の最深の受信深度である。データ位置Ｚ^(k)が最大値Ｚ^(k) _maxより大きい場合（ステップＳ２９：Ｙｅｓ）、周波数解析部３３１はカウンタｋを１増加させる（ステップＳ３０）。これは、処理をとなりの音線へ移すことを意味する。一方、データ位置Ｚ^(k)が最大値Ｚ^(k) _max以下である場合（ステップＳ２９：Ｎｏ）、周波数解析部３３１はステップＳ２３へ戻る。

ステップＳ３０の後、周波数解析部３３１は、カウンタｋが最大値ｋ_maxより大きいか否かを判定する（ステップＳ３１）。カウンタｋがｋ_maxより大きい場合（ステップＳ３１：Ｙｅｓ）、周波数解析部３３１は一連の周波数解析処理を終了する。一方、カウンタｋがｋ_max以下である場合（ステップＳ３１：Ｎｏ）、周波数解析部３３１はステップＳ２２に戻る。この最大値ｋ_maxは、図３中、解析範囲の最左の音線の番号である。

このようにして、周波数解析部３３１は、解析対象領域内の（ｋ_max−ｋ₀＋１）本の音線の各々について深度別に複数回のＦＦＴ演算を行う。ＦＦＴ演算の結果は、受信深度および受信方向とともに記憶部３７に格納される。

なお、これら４種の値ｋ₀、ｋ_max、Ｚ^(k) ₀、Ｚ^(k) _maxについては、図３の全走査範囲を含むようなデフォルト値が記憶部３７にあらかじめ記憶されており、周波数解析部３３１は適宜これらの値を読み取って、図９の処理を行う。デフォルト値を読み取った場合、周波数解析部３３１は全走査範囲に対して周波数解析処理を行う。しかし、この４種の値ｋ₀、ｋ_max、Ｚ^(k) ₀、Ｚ^(k) _maxは、入力部３５を通じた術者等のユーザによる関心領域の指示入力によって変更可能である。変更されていた場合、周波数解析部３３１は指示入力された関心領域においてのみ周波数解析処理を行う。

以上説明したステップＳ６の周波数解析処理に続いて、単回帰分析部３３２は、周波数解析部３３１が取得した複数のスペクトルデータの補正前特徴量をそれぞれ算出し、各スペクトルデータの補正前特徴量に対して超音波の減衰の影響を排除する減衰補正を行って各スペクトルデータの特徴量を算出する（ステップＳ７〜Ｓ８）。

ステップＳ７において、近似部３３２ａは、周波数解析部３３１が生成した解析範囲内の位置に応じた複数のスペクトルデータをそれぞれ単回帰分析することにより、各スペクトルデータに対応する補正前特徴量を算出する（ステップＳ７）。具体的には、近似部３３２ａは、各スペクトルデータを単回帰分析することによって一次式で近似し、補正前特徴量として傾きａ₀、切片ｂ₀、ミッドバンドフィットｃ₀を算出する。例えば、図５に示す直線Ｌ₁₀は、近似部３３２ａが周波数帯域ＵのスペクトルデータＣ₁に対し単回帰分析によって近似した回帰直線である。

続いて、減衰補正部３３２ｂは、近似部３３２ａが各スペクトルデータに対して近似して得た補正前特徴量に対し、減衰率ζを用いて減衰補正を行うことにより、減衰補正後の特徴量を算出し、記憶部３７に格納する（ステップＳ８）。図６に示す直線Ｌ₁は、減衰補正部３３２ｂが減衰補正処理を行うことによって得られる直線の例である。

ステップＳ８において、減衰補正部３３２ｂは、上述した式（３）、（５）における受信深度ｚに、超音波信号の音線のデータ配列を用いて得られるデータ位置Ｚ＝ｖ_S／（２・ｆ_sp）・Ｄ・ｎ_S＋Ｚ₀を代入することによって算出する。ここで、ｆ_spはデータのサンプリング周波数、ｖ_sは音速、Ｄはデータステップ幅、ｎ_Sは処理対象のサンプルデータ群のデータ位置までの音線の１番目のデータからのデータステップ数、Ｚ₀は解析範囲の最浅の受信深度である。例えば、データのサンプリング周波数ｆ_spを５０ＭＨｚとし、音速ｖ_sを１５３０ｍ／ｓｅｃとし、図５に示すデータ配列を採用してデータステップ幅Ｄを１５とすると、ｚ＝０．２２９５ｎ_S＋Ｚ₀（ｍｍ）となる。

その後、推定部３３３は、単回帰分析部３３２が減衰補正を行った特徴量ａおよび特徴量ｃと、関係情報記憶部３７１に記憶されている定数α、β、γ、α´、β´、γ´とを用いて、数密度ｎおよび直径ｄの対数をとったlogｎおよび／またはlogｄを推定する（ステップＳ９：推定ステップ）。推定部３３３は、推定したlogｎおよび／またはlogｄを物理量として物理量画像データ生成部３４２に出力する。

物理量画像データ生成部３４２は、Ｂモード画像データ生成部３４１が生成したＢモード画像データにおける各画素に対し、ステップＳ９で推定された物理量に関連づけた視覚情報（例えば色相）を重畳することによって物理量画像データを生成する（ステップＳ１０：物理量情報生成ステップ）。

この後、表示装置４は、制御部３６の制御のもと、物理量画像データ生成部３４２が生成した物理量画像データに対応する物理量画像を表示する（ステップＳ１１）。図１０は、表示装置４における特徴量画像の表示例を模式的に示す図である。同図に示す特徴量画像２０１は、Ｂモード画像に物理量に関する視覚情報が重畳された画像を表示する重畳画像表示部２０２と、観測対象の識別情報などを表示する情報表示部２０３とを有する。なお、情報表示部２０３に、物理量または特徴量の情報、近似式の情報、ゲインやコントラスト等の情報等をさらに表示するようにしてもよい。また、物理量画像に対応するＢモード画像を物理量画像と並べて表示してもよい。

以上説明してきた一連の処理（ステップＳ１〜Ｓ１１）において、ステップＳ２〜Ｓ４の処理とステップＳ５〜Ｓ１０の処理とを並行して行うようにしてもよい。

以上説明してきた一連の処理において、観測対象の物理量として散乱体の直径と散乱体の数密度を、変数変換として対数変換を例に説明した。しかし、これはこの組み合わせによらず、他の物理量でもよく、他の変数変換でもよい。本来望んでいた物理量を、特徴量がそれに対し線形に変化する、もしくは、ゆるやかに変化する変数に変換できれば、変換後の変数を線形の重回帰分析に帰着、もしくは、近似して推定し、本来望んでいた物理量をある程度推定することができる。本実施の形態１では、観測対象の散乱体の直径と散乱体の数密度のオーダー（桁数−１）を推定することができた。

以上説明した本発明の実施の形態１によれば、単回帰分析部３３２により算出された特徴量と、基準ファントムを用いて算出された関係式（定数α、β、γ、α´、β´、γ´）とを用いて、推定部３３３が数密度ｎおよび直径ｄの対数をとったlogｎおよび／またはlogｄを推定するようにしたので、散乱体の数密度および／または大きさ（直径）を直接推定することを可能にした。これにより、熟練を要さずとも、特徴量として算出された値が、どういう病理に対応した組織性状であるかという病理解釈を容易、かつ確実に行うことができる。

また、本発明の実施の形態１によれば、基準ファントムからの減衰補正された特徴量を基に定数α、β、γ、α´、β´、γ´を求めるようにしたので、観測対象と基準ファントムとの間の減衰率の相違によらず、より正確に上記推定をすることができる。

なお、入力部３５を通じた術者等のユーザによる関心領域の指示入力によって、音線幅を決定するｋ₀、ｋ_max、および、深度幅を決定するＺ^(k) ₀、Ｚ^(k) _maxを、周波数解析部３３１が変更することで、指示入力された特定の深度幅および音線幅で区切られる関心領域に対してのみスペクトルデータを算出できるよう構成した。そのため、算出に関わる演算量を下げ、フレームレートを上げることができる。なお、ここでは、関心領域を深度幅と音線幅とで扇形に区切ったが、この例に限らず長方形、楕円であってもよく、他の形状でもよい。この場合、単回帰分析部３３２が、設定された関心領域内とその関心領域外とで個別に最適な減衰率を設定するようにしてもよい。

（実施の形態１の変形例）
続いて、本発明の実施の形態１の変形例について説明する。上述した実施の形態１では、推定部３３３が数密度ｎおよび直径ｄの対数をとったlogｎおよび／またはlogｄを推定し、物理量画像データ生成部３４２が、logｎおよび／またはlogｄに対応する視覚情報をＢモード画像に重畳するものとして説明したが、本変形例では、演算部３３が、推定部３３３が推定したlogｎおよび／またはlogｄを、さらに、ｎおよび／またはｄに変換する。図１１は、本発明の実施の形態１の変形例に係る超音波診断装置を備えた超音波診断システム１ａの構成を示すブロック図である。

本変形例に係る超音波診断システム１ａは、上述した実施の形態１に係る超音波診断システム１の構成に対して、超音波診断装置３の演算部３３が、変数変換部３３４をさらに備える。変数変換部３３４は、推定部３３３が推定したlogｎおよび／またはlogｄをｎおよび／またはｄに変換する。具体的には、logｎおよび／またはlogｄの推定値を１０の指数に代入してｎおよび／またはｄとする。変数変換部３３４は、変換後の数密度ｎおよび／または直径ｄを、物理量画像データ生成部３４２に出力する。

物理量画像データ生成部３４２は、変数変換部３３４が変換した物理量（数密度ｎおよび／または直径ｄ）に関連する視覚情報をＢモード画像データにおける画像の各画素に対して重畳することによって物理量画像データを生成する。物理量画像データ生成部３４２は、例えば図４に示す１つのサンプルデータ群Ｆ_j（ｊ＝１、２、・・・、Ｋ）のデータ量に対応する画素領域に対し、そのサンプルデータ群Ｆ_jから算出される周波数スペクトルの特徴量に対応する物理量に関連する視覚情報を割り当てる。

表示装置４は、制御部３６の制御のもと、物理量画像データ生成部３４２が生成した物理量画像データに対応する物理量画像を表示する。図１２は、表示装置４における特徴量画像の表示例を模式的に示す図である。同図に示す特徴量画像２０１は、物理量に関する視覚情報が重畳された画像をＢモード画像に表示する重畳画像表示部２０２と、観測対象の識別情報などを表示する情報表示部２０３とを有する。本変形例では、情報表示部２０３に、物理量の情報として、変数変換部３３４により変換された散乱体の数密度ｎおよび／または直径ｄを表示する。

本変形例によれば、推定部３３３が推定したlogｎおよび／またはlogｄを変数変換して散乱体の数密度ｎおよび直径ｄを物理量として表示するようにしたので、散乱体の数密度および／または大きさ（直径）に直接関連した視覚情報を用いて推定することを可能にした。これにより、熟練を要さずとも、特徴量として算出された値が、どういう病理に対応した組織性状であるかという病理解釈を容易、かつ確実に行うことができる。

また、本変形例によれば、観測対象や基準ファントムの散乱体の数密度ｎ、散乱体の直径ｄが直接には特徴量ｃ、特徴量ａと線形の関係ではなくとも、特徴量ｃ、特徴量ａとほぼ線形の関係になるlogｎ、logａを用いることができる。そして、logｎ、logａを変数変換することで本来望んでいた物理量である散乱体の数密度ｎ、散乱体の直径ｄについて、そのオーダー（桁数−１）だけではなく、値自身を線形の重回帰分析に帰着して推定することができる。

（実施の形態２）
続いて、本発明の実施の形態２について説明する。本実施の形態２は、上述した実施の形態１に係る超音波診断装置を備えた超音波診断システムの構成と同様の構成を備える。図１は共通である。上述した実施の形態１では、単回帰分析部３３２が単回帰分析により特徴量ａ、ｂ、ｃを算出するものとして説明したが、本実施の形態では単回帰分析部３３２が観測対象の特徴量として、これら３特徴量に加えて減衰率ζを算出する。減衰率ζは、例えば上述したスペクトルデータＦ（ｆ，Ｌ）を用いて算出することができる。具体的には、回帰分析によりスペクトルデータＦ（ｆ，Ｌ）の周波数ｆに対する回帰直線の傾きである特徴量ａを求め、次にその傾きの往復距離Ｌに対するさらなる傾きを求め、さらに−１を乗ずることで求められる。また、上述した実施の形態１では、推定部３３３が数密度ｎおよび直径ｄの対数をとったlogｎおよび／またはlogｄを推定するものとして説明したが、本実施の形態２では、推定部３３３が、数密度ｎおよび直径ｄの対数をとったlogｎおよび／またはlogｄに加えて、散乱強度（以下、散乱強度をｒとする）を推定する。ここでいう散乱強度とは、散乱体と媒体の振幅反射率、エネルギー反射率、および／または、それらの関数である。ここで、散乱強度は、以下のように定義される振幅反射率またはエネルギー反射率である。Ｚは、音響インピーダンスを示す。
振幅反射率 ＝ ｜Ｚ_散乱体−Ｚ_媒体｜／｜Ｚ_散乱体＋Ｚ_媒体｜
エネルギー反射率 ＝ ｜Ｚ_散乱体−Ｚ_媒体｜²／｜Ｚ_散乱体＋Ｚ_媒体｜²

本実施の形態２において、推定部３３３は、単回帰分析部３３２が減衰補正を行った特徴量ａおよび特徴量ｃおよび特徴量としての減衰率ζと、記憶部３７に記憶されている定数α、β、γ、δ、α´、β´、γ´、δ´、α´´、β´´、γ´´、δ´´とを用いて、数密度ｎおよび直径ｄの対数をとったlogｎおよびlogｄ、ならびに散乱強度ｒを推定する。具体的に、特徴量ａ（傾き）、特徴量ｃ（ミッドバンドフィット）および減衰率ζが、logｎ、logｄおよび散乱強度ｒと線形関係を有するものと仮定し、推定部３３３は、後述する４次元空間内における回帰平面（つまり３次元立体）の方程式（１５）〜（１７）から導いた下式（１９）をもとに、logｎ、logｄおよび散乱強度ｒを推定する。推定部３３３は、推定したlogｎ、logｄおよび／または散乱強度ｒを物理量として画像処理部３４に出力する。

上式（１５）〜（１７）から、下式（１８）が得られる。

式（１８）より下式（１９）が得られる。

ここで、定数α、β、γ、δ、α´、β´、γ´、δ´、α´´、β´´、γ´´、δ´´の求め方について説明する。以下では、一例として定数α、β、γ、δの求め方について説明する。定数α、β、γ、δは、散乱体の大きさ（直径）、数密度および散乱強度を調整した材料を均一に混ぜて分布させて固化することによって形成され、散乱体の大きさ（直径）および数密度が予めわかっているとともに、減衰率［ｄＢ／ｃｍ／ＭＨｚ］も均一であることが予めわかっている基準ファントムに超音波を送信し、得られた超音波エコーに基づいて算出された特徴量ｃ（ミッドバンドフィット）を算出してそれぞれ得ることができる。本実施の形態２では、上述した実施の形態１と同様に、散乱体の直径および数密度および散乱強度の組み合わせが異なるＮ個の基準ファントムを用いて、特徴量ｃをそれぞれ算出する。以下の説明では、Ｎ個のファントムのうち、ｉ番目の基準ファントムの散乱体の数密度をｎ_i、直径をｄ_i、散乱強度をｒ_i、減衰率をζ_iとし、このｉ番目の基準ファントムにより得られた特徴量（ミッドバンドフィット）をｃ_iとする（１≦ｉ≦Ｎ）。この際、特徴量ｃ_iは、ｉ番目の基準ファントムの減衰率ζ_iを基に、ｉ番目の基準ファントムのスペクトルデータが算出された深度および式（５）を用いて減衰の影響が排除され、深度に依存しない。

Ｎ個の基準ファントムから特徴量ｃ_iをそれぞれ取得後、散乱体の数密度ｎ_i、直径ｄ_i、散乱強度ｒ_iおよび特徴量ｃ_iの全基準ファントムにわたる平均をそれぞれ求める。その後、下式（２０）、（２１）に示す行列Ｇおよび列ベクトルＹを定義する。４次元平面内の回帰平面の方程式の定数αとβとγとは、式（２０）、（２１）で定義される下式（２２）によって得ることができる。

ここで、左肩のｔは転置行列を示し、右肩の−１は逆行列を示している。

また、回帰平面が、散乱体の数密度ｎ_i、直径ｄ_i、散乱強度ｒ_iおよび特徴量ｃ_iの各平均を座標とする点（重心）を通過することから、式（１５）より下式（２３）を得る。上式（２２）で得た定数α、β、γを下式（２３）へ代入して定数δを求めることができる。

上述したようにして、Ｎ個の基準ファントムを用いて、上式（１５）の定数α、β、γ、δを求めることができる。上式（１６）の定数α´、β´、γ´、δ´についても同様に、Ｎ個の基準ファントムを用いて、特徴量ａ_iを算出して、定数α´、β´、γ´、δ´を求めることができる。また、上式（１７）の定数α´´、β´´、γ´´、δ´´についても同様に、Ｎ個の基準ファントムを用いて、改めて測定した減衰率ζ_iもしくは既知の減衰率ζ_iを用いて、定数α´´、β´´、γ´´、δ´´を求めることができる。この取得した定数α、β、γ、δ、α´、β´、γ´、δ´、α´´、β´´、γ´´、δ´´は、記憶部３７に記憶される。

この後、超音波内視鏡２を介して得られた観測対象の特徴量ｃ、特徴量ａ、減衰率ζ、および、観測対象内の散乱体の直径ｄの桁数−１であるlogｄ、散乱体の数密度ｎの桁数−１であるlogｎ、散乱体の散乱強度ｒも、基準ファントムから得られた物理的な傾向に従うと考えられる。そのため、データセット（ｎ，ｄ，ｒ，ｃ）は、回帰平面の方程式（１５）を満足する。同様に、観測対象からのデータセット（ｎ，ｄ，ｒ，ａ）は式（１６）を満足し、観測対象からのデータセット（ｎ，ｄ，ｒ，ζ）は式（１７）を満足する。観測対象のｎ，ｄ，ｒ，ｃ，ａ，ζは、式（１５）と式（１６）と式（１７）とを満足するので、式（１９）を満足する。よって、推定部３３３は、このように観測対象からの特徴量ｃ、特徴量ａ、特徴量としての観測対象の減衰率ζを式（１９）へ代入し、logｎおよびlogｄおよび散乱強度ｒを推定する。

なお、本実施の形態２では、上述した基準ファントムからの超音波エコーに基づいたエコー信号の受信、特徴量ｃ、特徴量ａ、減衰率ζの算出、Ｎ個のデータセットの算出、定数α、β、γ、δ、α´、β´、γ´、δ´、α´´、β´´、γ´´、δ´´の算出は、超音波診断装置３の外部で実施される。そして、取得した定数α、β、γ、δ、α´、β´、γ´、δ´、α´´、β´´、γ´´、δ´´は、工場出荷前に、入力部３５を介して記憶部３７に内蔵された関係情報記憶部３７１に記憶される。

物理量画像データ生成部３４２は、推定部３３３が変換した物理量（logｎ、logｄおよび／またはｒ）に関連する視覚情報をＢモード画像データにおける画像の各画素に対して重畳することによって物理量画像データを生成する。物理量画像データ生成部３４２は、例えば図４に示す１つのサンプルデータ群Ｆ_j（ｊ＝１、２、・・・、Ｋ）のデータ量に対応する画素領域に対し、そのサンプルデータ群Ｆ_jから算出される周波数スペクトルの特徴量に対応する物理量に関連する視覚情報を割り当てる。

以上説明した本発明の実施の形態２によれば、単回帰分析部３３２により算出された特徴量と、減衰率と、基準ファントムを用いて算出された関係式（定数α、β、γ、δ、α´、β´、γ´、δ´、α´´、β´´、γ´´、δ´´）とを用いて、推定部３３３が数密度ｎおよび直径ｄの対数をとったlogｎ、logｄおよび／または散乱強度ｒを推定するようにしたので、散乱体の数密度、大きさ（直径）および／または散乱強度を直接推定することを可能にした。これにより、熟練を要さずとも、特徴量として算出された値が、どういう病理に対応した組織性状であるかという病理解釈を容易、かつ確実に行うことができる。

また、本発明の実施の形態２によれば、基準ファントムからの減衰補正された特徴量を基に定数α、β、γ、δ、α´、β´、γ´、δ´、α´´、β´´、γ´´、δ´´を求めるようにしたので、観測対象と基準ファントムとの間の減衰率の相違によらず、より正確に散乱体の数密度ｎの桁数−１であるlogｎ、および、散乱体の直径ｄの桁数−１であるlogｄおよび／または散乱強度ｒを推定することができる。

なお、上述した実施の形態２では、物理量として減衰率ζを用いるものとして説明したが、これに限らず、特徴量ａ、ｃの分散や、スペクトルデータにより重み付けした平均周波数、音速を用いるものであってもよい。ここで、音速は、複数の素子をアレイ状に設けた超音波振動子２１を電子的に走査させる構成を有し、送受信部３１が複数の素子に対応したビーム合成用の多チャンネル回路を有する場合であって、超音波振動子２１からのエコー信号のフォーカスの条件が最も良い場合の各素子に対する受信電圧の遅延時間から推定される。また、平均周波数は、以下の式（２４）により与えられる（１≦ｑ≦Ｎ_f：Ｎ_fは１より大きい整数）。

ここで、Σ_qｆ_q・Ｆ（ｆ_q，Ｌ）はスペクトル強度の重み付け加算を示している。

以上説明してきた一連の処理において、観測対象の物理量として散乱体の散乱強度のほか、散乱体の直径と散乱体の数密度を、変数変換として対数変換を例に説明した。しかし、これはこの組み合わせによらず、他の物理量でもよく、他の変数変換でもよい。本来望んでいた物理量を、特徴量がそれに対し線形に変化する、もしくは、ゆるやかに変化する変数に変換できれば、変換後の変数を線形の重回帰分析に帰着、もしくは、近似して推定し、本来望んでいた物理量をある程度推定することができる。本実施の形態２では、観測対象の散乱体の直径と散乱体の数密度のオーダー（桁数−１）を推定することができた。

（実施の形態３）
続いて、本発明の実施の形態３について説明する。上述した実施の形態１では、基準ファントムからの超音波エコーに基づいた定数α、β、γ、α´、β´、γ´が超音波診断装置ではなくその外部で算出され、関係情報記憶部３７１に予め記憶されているものとして説明したが、本実施の形態３では、超音波診断システム１ｂが、定数α、β、γ、α´、β´、γ´を算出可能な構成をなす。図１３は、本発明の実施の形態３に係る超音波診断装置を備えた超音波診断システムの構成を示すブロック図である。

本実施の形態３に係る超音波診断システム１ｂは、上述した実施の形態１に係る超音波診断システム１の構成に対して、超音波診断装置３が、変数変換部３８および重回帰分析部３９をさらに備える。

次に、変数変換部３８は、入力部３５または記憶部３７を介してＮ個の基準ファントムの散乱体の数密度ｎ_iおよび直径ｄ_iおよび減衰率ζ_iを取得する（１≦ｉ≦Ｎ）。変数変換部３８は、取得した数密度ｎ_iおよび直径ｄ_iを対数に変換してlogｎ_iおよびlogｄ_iを得る。変数変換部３８は、logｎ_iおよびlogｄ_iおよび減衰率ζ_iを記憶部３７に出力する。これらの処理をＮ個全ての基準ファントムについて実施する。

次に、超音波振動子２１がＮ個の基準ファントムのうちｉ番目の基準ファントムへ超音波を送信する。その後、送受信部３１、周波数解析部３３１、単回帰分析部３３２は、実施の形態１で観測対象からのエコー信号に施した処理と同様の処理をｉ番目の基準ファントムからのエコー信号に施す。このようにして、超音波診断システム１ｂは、ｉ番目の基準ファントムからの超音波エコーに基づいた特徴量ｃ_i、特徴量ａ_iを算出する（１≦ｉ≦Ｎ）。この際、減衰補正部３３２ｂは、ｉ番目の基準ファントムの減衰率ζ_iを記憶部３７から取得し、ｉ番目の基準ファントムのスペクトルデータが算出された深度および式（３）と、式（５）とを用い、減衰を補正することによって特徴量ｃ_i、特徴量ａ_iから減衰の影響を排除しておく。基準ファントムにおいて減衰は一様であるから、これらの特徴量ｃ_i、特徴量ａ_iは、深度に依存しない。そして、単回帰分析部３３２は、ｉ番目の基準ファントムの散乱体の特徴量ｃ_i、特徴量ａ_iを記憶部３７へ出力する。超音波診断システム１ｂは、これらの処理をＮ個全ての基準ファントムについて実施する。

重回帰分析部３９は、記憶部３７からＮ個の基準ファントムのlogｎ_i、logｄ_i、特徴量ｃ_i、特徴量ａ_iを全て取得する（１≦ｉ≦Ｎ）。そして、重回帰分析部３９は、Ｎ個の基準ファントムの散乱体の数密度ｎ_i、直径ｄ_i、特徴量ｃ_i、特徴量ａ_iの全基準ファントムにわたる平均をそれぞれ求める。その後、重回帰分析部３９は、実施の形態１と同様の演算方法で、上式（１１）〜（１４）を用いて上式（７）、（８）の定数α、β、γ、α´、β´、γ´を求める。重回帰分析部３９は、求めた定数α、β、γ、α´、β´、γ´または関係式を関係情報記憶部３７１に出力する。関係情報記憶部３７１は、受信した定数α、β、γ、α´、β´、γ´または関係式を記憶する。以上、定数α、β、γ、α´、β´、γ´または関係式を関係情報記憶部３７１に記憶するまでの処理は、観測対象への超音波の送信の前、望ましくは、工場出荷時に実施しておく。その後、観測対象への超音波の送信を行う。

入力部３５が観測対象の走査を指示する指示入力を術者から受け付けると、観測対象への超音波の送信が開始され、送受信部３１は超音波振動子２１からエコー信号を受信する。その後、送受信部３１、周波数解析部３３１、単回帰分析部３３２は、観測対象からのエコー信号に、実施の形態１と同様の処理を施す。こうして、単回帰解析部３３２は、観測対象からの超音波エコーを基に、減衰補正部３３２ｂにより適切に減衰を補正した観測対象の特徴量ｃ、特徴量ａを算出し、推定部３３３へ出力する。推定部３３３や物理量画像データ生成部３４２などの処理は、上述した実施の形態１と同様である。

以上説明した本発明の実施の形態３によれば、上述した実施の形態１の効果を得ることができるとともに、超音波診断装置３において、定数α、β、γ、α´、β´、γ´を求めるようにしたので、新たな超音波振動子を有する超音波内視鏡が追加されたり、基準ファントムが追加されたりした場合など、装置単位でアップデートを行うことができる。

また、本発明の実施の形態３によれば、超音波診断装置３において、定数α、β、γ、α´、β´、γ´を求めるようにした。そのため、超音波振動子ごとに送受信感度等の特性がばらつく場合でも超音波振動子の機体ごとに定数α、β、γ、α´、β´、γ´を算出することが容易になり、各超音波振動子に対応して、ばらつきによらず、より正確かつ容易に、散乱体の数密度ｎやその対数logｎ、および散乱体の直径ｄやその対数logｄを推定することができる。

また、超音波振動子の特性が経時的に変化しても、定期的に定数α、β、γ、α´、β´、γ´を算出することが容易になり、経時的変化によらず、より正確に、上記推定をすることができる。

また、本発明の実施の形態３によれば、基準ファントムからの減衰補正された特徴量を基に定数α、β、γ、α´、β´、γ´を求めるようしたので、観測対象と基準ファントムとの間の減衰率の相違によらず、より正確にlogｎ、logｄを推定することができる。

また、本発明の実施の形態３によれば、基準ファントムの散乱体の数密度ｎ_i、散乱体の直径ｄ_iが直接には特徴量ｃ_i、ａ_iには線形の関係がなくとも、特徴量ｃ_i、ａ_iとほぼ線形の関係になるlogｎ_i、logｄ_iを用いることができる。そして、本来望んでいた物理量である観測対象の散乱体の数密度ｎ、散乱体の直径ｄの情報の画像化を線形の重回帰分析に帰着して実現することができる。

なお、本実施の形態３に係る変数変換部３８および重回帰分析部３９を、上述した実施の形態１の変形例や、実施の形態２に適用することも可能である。

（実施の形態４）
続いて、本発明の実施の形態４について説明する。上述した実施の形態１では、関係情報記憶部３７１が、定数α、β、γ、α´、β´、γ´または関係式を予め記憶しているものとして説明したが、本実施の形態４では、関係情報記憶部３７１が、特徴量の値に応じて物理量を出力可能なルックアップテーブル（Lookup table：ＬＵＴ）を記憶している。図１４は、本発明の実施の形態４に係る超音波診断装置を備えた超音波診断システムの構成を示すブロック図である。

本実施の形態４に係る超音波診断システム１ｃは、上述した実施の形態１に係る超音波診断システム１の構成に対して、関係情報記憶部３７１が、ＬＵＴ３７１ａを記憶する。図１５は、本発明の実施の形態４に係る超音波診断装置が記憶するルックアップテーブルを説明するための図である。図１５では、一例として、特徴量ａ、ｃの入力により物理量として散乱体の数密度ｎの対数であるlogｎを出力するためのルックアップテーブルを示している。

ＬＵＴ３７１ａは、例えばlogｎを出力するためのルックアップテーブルの場合、上述した重回帰分析によって得られた定数α、β、γ、α´、β´、γ´と、所定の間隔で抽出される複数の仮の特徴量ａ、ｃとをもとに、上式（１０）よりlogｎおよびlogｄをそれぞれ算出し、縦軸を仮の特徴量ｃ、横軸を仮の特徴量ａとして、各セルにlogｎの値を代入することで生成される。logｄについても同様にして、定数α、β、γ、α´、β´、γ´と仮の特徴量ａ、ｃとをもとに、logｎを出力するためのルックアップテーブルを生成することができる。ＬＵＴ３７１ａは、超音波診断装置３で作成するものであってもよいし、外部の演算装置によって生成されたものを取得してもよいし、ネットワークを介して取得してもよい。なお、ルックアップテーブルを生成する際、Ｎ個の基準ファントムの特徴量ｃ_iと特徴量ａ_iと減衰率ζ_iとを用いる（１≦ｉ≦Ｎ）。これら特徴量は減衰率、基準ファントムのスペクトルデータが算出された深度、および式（５）、式（３）を用いて減衰の影響が排除され、深度に依存しないようにした上でルックアップテーブルを生成する。

関係情報記憶部３７１は、ＬＵＴ３７１ａとして、logｎを出力するためのルックアップテーブルと、logｄを出力するためのルックアップテーブルとを記憶する。推定部３３３は、単回帰分析部３３２から特徴量ａ、ｃを受信すると、例えば特徴量ａの小数点以下第１桁、特徴量ｃの小数点以下第２桁を丸める。その後、推定部３３３は、logｎを出力する場合は、関係情報記憶部３７１からlogｎを出力するためのルックアップテーブルを読出して参照し、logｎを推定する。また、推定部３３３は、logｄを出力する場合は、関係情報記憶部３７１からlogｄを出力するためのルックアップテーブルを読出して参照し、logｄを推定する。

以上説明した本発明の実施の形態４によれば、上述した実施の形態１の効果を得ることができるとともに、物理量を算出するために要する処理時間を短縮することができる。

また、本発明の実施の形態４によれば、ルックアップテーブルを生成する際、基準ファントムの特徴量から減衰の影響が排除され、深度に依存しないようにした上でルックアップテーブルを生成したため、深度別にルックアップテーブルを用意する必要がなく、データ量は小さく、処理も簡便である。

なお、上述した実施の形態１の変形例や、実施の形態２において、本実施の形態４に係るＬＵＴ３７１ａを記憶するようにすることも可能である。

ここまで、本発明を実施するための形態を説明してきたが、本発明は、上述した実施の形態によってのみ限定されるべきものではない。例えば、超音波診断装置において、各機能を有する回路同士をバスで接続することによって構成してもよいし、一部の機能が他の機能の回路構造に内蔵されるように構成してもよい。

また、上述した実施の形態１〜４では、散乱体の大きさが直径であるものとして説明したが、半径であってもよいし、体積であってもよい。

また、上述した実施の形態１〜４では、対照物体として、材質、質量密度、音速、音響インピーダンスが既知である媒体に、材質、質量密度、音速、音響インピーダンス、直径、数密度がやはり既知である散乱体を一様に混入させた基準ファントムを例に挙げて説明した。しかし、散乱体の直径、散乱体の散乱強度、散乱体の数密度等の物理量が既知で、かつ、分布が一様な対象であれば基準ファントムをこれに代えることができる。例えば、物理量が既知ないし正確に測定できれば動物の肝臓等、特定組織を用いてもよい。

また、上述した実施の形態１〜４では、物理量として、散乱体の数密度ｎ、散乱体の直径ｄ、散乱体の散乱強度ｒを例に挙げて説明した。さらに、特徴量として、スペクトルの傾きａ、スペクトルの切片ｂ、スペクトルのミッドバンドフィットｃ、減衰率ζを例に挙げて説明した。しかし、この中の一部であってもよいし、全部であってもよい。また、このほかの物理量や特徴量でもよい。例えば物理量として散乱体の形状を示す量、例えば異型性を示すフラクタル次元を用いてもよい。また、例えば物理量として散乱体の直径の分布の分散を用いてもよい。特徴量として音速や質量密度を用いてもよい。

また、上述した実施の形態１〜４において、推定部３３３が推定した物理量が複数存在する場合、表示装置４には、複数の物理量を同時に表示するものであってもよいし、入力部３５への指示入力により順次切り替えて表示してもよいし、異なるタイミング（フレーム）で表示するものであってもよい。

また、本実施の形態１〜４では、超音波プローブとしてライトガイド等の光学系を有する超音波内視鏡２を用いて説明したが、超音波内視鏡２に限らず、撮像光学系および撮像素子を有しない超音波プローブであってもよい。さらに、超音波プローブとして、光学系のない細径の超音波ミニチュアプローブを適用してもよい。超音波ミニチュアプローブは、通常、胆道、胆管、膵管、気管、気管支、尿道、尿管へ挿入され、その周囲臓器（膵臓、肺、前立腺、膀胱、リンパ節等）を観察する際に用いられる。

また、超音波プローブとして、観測対象の体表から超音波を照射する体外式超音波プローブを適用してもよい。体外式超音波プローブは、通常、腹部臓器（肝臓、胆嚢、膀胱）、乳房（特に乳腺）、甲状腺を観察する際に体表に直接接触させて用いられる。

また、超音波振動子は、リニア振動子でもラジアル振動子でもコンベックス振動子でも構わない。超音波振動子がリニア振動子である場合、その走査領域は矩形（長方形、正方形）をなし、超音波振動子がラジアル振動子やコンベックス振動子である場合、その走査領域は扇形や円環状をなす。また、超音波内視鏡は、超音波振動子をメカ的に走査させるものであってもよいし、超音波振動子として複数の素子をアレイ状に設け、送受信にかかわる素子を電子的に切り替えたり、各素子の送受信に遅延をかけたりすることで、電子的に走査させるものであってもよい。

このように、本発明は、請求の範囲に記載した技術的思想を逸脱しない範囲内において、様々な実施の形態を含みうるものである。

以上のように、本発明にかかる超音波診断装置、超音波診断装置の作動方法および超音波診断装置の作動プログラムは、特徴量に基づく組織性状の鑑別を容易かつ正確に行うのに有用である。

１、１ａ、１ｂ、１ｃ 超音波診断システム
２ 超音波内視鏡
３ 超音波診断装置
４ 表示装置
２１ 超音波振動子
３１ 送受信部
３２ 信号処理部
３３ 演算部
３４ 画像処理部
３５ 入力部
３６ 制御部
３７ 記憶部
３８、３３４ 変数変換部
３９ 重回帰分析部
２０１ 特徴量画像
２０２ 重畳画像表示部
２０３ 情報表示部
３３１ 周波数解析部
３３２ 単回帰分析部
３３２ａ 近似部
３３２ｂ 減衰補正部
３３３ 推定部
３４１ Ｂモード画像データ生成部
３４２ 物理量画像データ生成部
３７１ 関係情報記憶部

Claims (15)

1. 観測対象へ超音波を送信し、該観測対象で後方散乱された超音波を受信する超音波振動子を備えた超音波プローブが取得した超音波信号に基づいて超音波画像を生成する超音波診断装置であって、
   前記観測対象から受信した前記超音波信号に基づき、前記超音波信号に関連する特徴量を算出する特徴量算出部と、
   含有する散乱体の物理量が既知の対象の前記物理量、および該既知の対象から得られた特徴量に基づき導出された関係と、前記特徴量算出部が算出した前記観測対象の前記特徴量とを用いて、前記観測対象が含む散乱体の物理量を推定する推定部と、
   前記推定部により推定された推定結果を含み、表示部に表示させるための情報を生成する物理量情報生成部と、
   を備えたことを特徴とする超音波診断装置。
2. 前記物理量は、前記対象が含有する散乱体の数密度、前記散乱体の大きさ、および前記散乱体の散乱強度の少なくとも１つを含む、
   ことを特徴とする請求項１に記載の超音波診断装置。
3. 前記推定部は、前記特徴量算出部が算出した前記特徴量を、前記関係としての関係式に代入することにより前記観測対象が含む散乱体の物理量を推定する
   ことを特徴とする請求項１に記載の超音波診断装置。
4. 前記関係式、該関係式の係数、該関係式の定数項、および前記関係を記述したテーブルのうちの少なくとも一つを記憶する関係情報記憶部をさらに備え、
   前記推定部は、前記関係情報記憶部を参照して前記観測対象が含む散乱体の物理量を推定する
   ことを特徴とする請求項３に記載の超音波診断装置。
5. 前記関係式は、前記物理量および前記特徴量のうち少なくとも一部を重回帰分析することにより導出される
   ことを特徴とする請求項３に記載の超音波診断装置。
6. 前記既知の散乱体の物理量のうち、前記散乱体の数密度、および前記散乱体の大きさを含む場合、
   前記散乱体の数密度、および前記散乱体を非線形変換して前記重回帰分析を行う
   ことを特徴とする請求項５に記載の超音波診断装置。
7. 前記特徴量は、前記超音波信号をもとに算出される周波数特徴量を含む
   ことを特徴とする請求項１に記載の超音波診断装置。
8. 前記特徴量は、前記超音波信号をもとに算出される減衰率を含む
   ことを特徴とする請求項１に記載の超音波診断装置。
9. 前記特徴量は、前記超音波信号をもとに算出される音速を含む
   ことを特徴とする請求項１に記載の超音波診断装置。
10. 前記物理量情報生成部は、前記推定部が推定した前記物理量に応じて視覚情報を付与した画像データを生成する
    ことを特徴とする請求項１に記載の超音波診断装置。
11. 前記推定部が互いに異なる複数の物理量を推定した場合、
    前記物理量情報生成部は、前記複数の物理量を同時、順次または異なるタイミングで前記表示部に表示される前記情報を生成する
    ことを特徴とする請求項１に記載の超音波診断装置。
12. 前記非線形変換が施された前記散乱体の数密度、および／または前記散乱体の大きさに対してさらに非線形変換を施す変数変換部
    をさらに備えたことを特徴とする請求項６に記載の超音波診断装置。
13. 前記観測対象から受信した前記超音波信号に基づき算出した特徴量において減衰の影響を補正する減衰補正部をさらに備え、
    前記推定部が、前記含有する散乱体の物理量が前記既知の対象の前記物理量と、前記既知の対象から得られた特徴量において減衰の影響を補正した値とに基づき導出された関係と、前記減衰補正部が補正した前記観測対象の前記特徴量の値とを用いて、前記観測対象が含む散乱体の物理量を推定する
    ことを特徴とする請求項１に記載の超音波診断装置。
14. 観測対象へ超音波を送信し、該観測対象で後方散乱された超音波を受信する超音波振動子を備えた超音波プローブが取得した超音波信号に基づいて超音波画像を生成する超音波診断装置の作動方法であって、
    特徴量算出部が、前記観測対象から受信した前記超音波信号に基づき、前記超音波信号に関連する特徴量を算出する特徴量算出ステップと、
    推定部が、含有する散乱体の物理量が既知の対象の前記物理量、および該既知の対象から得られた特徴量に基づき導出された関係と、前記特徴量算出部が算出した前記観測対象の前記特徴量とを用いて、前記観測対象が含む散乱体の物理量を推定する推定ステップと、
    を含むことを特徴とする超音波診断装置の作動方法。
15. 観測対象へ超音波を送信し、該観測対象で後方散乱された超音波を受信する超音波振動子を備えた超音波プローブが取得した超音波信号に基づいて超音波画像を生成する超音波診断装置の作動プログラムであって、
    特徴量算出部が、前記観測対象から受信した前記超音波信号に基づき、前記超音波信号に関連する特徴量を算出する特徴量算出手順と、
    推定部が、含有する散乱体の物理量が既知の対象の前記物理量、および該既知の対象から得られた特徴量に基づき導出された関係と、前記特徴量算出部が算出した前記観測対象の前記特徴量とを用いて、前記観測対象が含む散乱体の物理量を推定する推定手順と、
    を前記超音波診断装置に実行させることを特徴とする超音波診断装置の作動プログラム。

Images