(57)【特許請求の範囲】
【請求項１】
連続したフレームからなる超音波Ｂモード画像列を取得して、所定のフレームに描画される頸動脈短軸断面の頸動脈壁を推定する頸動脈壁推定手段と、
推定された頸動脈壁及びその周辺組織を含む領域の画像をテンプレート画像とし、このテンプレート画像の頸動脈壁と前記周辺組織を含む領域の画像を、変形させて、変形後のテンプレート画像と、前記取得した各フレームの画像との誤差を最小にすることにより頸動脈壁の径の大きさを推定して頸動脈壁の時間変化を取得する心拍変動追跡手段を備えることを特徴とする画像処理装置。
【請求項２】
前記頸動脈壁推定手段は、
頸動脈の位置と形状を決める頸動脈パラメータが含まれて頸動脈であるらしさを表わす評価関数を最大にすることにより、頸動脈壁を推定することを特徴とする請求項１に記載の画像処理装置。
【請求項３】
前記心拍変動追跡手段は、前記フレームが、フレーム毎に頸動脈の中心が変動している場合には、頸動脈の中心座標のずれも推定することを特徴とする請求項１又は請求項２に記載の画像処理装置。
【請求項４】
コンピュータを、
連続したフレームからなる超音波Ｂモード画像列を取得して、所定のフレームに描画される頸動脈短軸断面の頸動脈壁を推定する頸動脈壁推定手段と、
推定された頸動脈壁及びその周辺組織を含む領域の画像をテンプレート画像とし、このテンプレート画像の頸動脈壁と前記周辺組織を含む領域の画像を、変形させて、変形後のテンプレート画像と、前記取得した各フレームの画像との誤差を最小にすることにより頸動脈壁の径の大きさを推定して頸動脈壁の時間変化を取得する心拍変動追跡手段として機能させることを特徴とする画像処理プログラム。
【請求項５】
前記頸動脈壁推定手段が、
頸動脈の位置と形状を決める頸動脈パラメータが含まれて頸動脈であるらしさを表わす評価関数を最大にすることにより、頸動脈壁を推定することを特徴とする請求項４に記載の画像処理プログラム。
【請求項６】
前記心拍変動追跡手段が、
前記フレームが、フレーム毎に頸動脈の中心が変動している場合には、頸動脈の中心座標のずれも推定することを特徴とする請求項４又は請求項５に記載の画像処理プログラム。
【請求項７】
請求項４乃至請求項６のいずれか１項に記載の画像処理プログラムを記憶した記憶媒体。
【請求項８】
請求項１乃至請求項３のうちいずれか１項に記載の画像処理装置に接続され、探触子の超音波送受信により得られたエコー信号を基に生体内の超音波Ｂモード画像列を生成し、前記画像処理装置に前記超音波Ｂモード画像列を出力することを特徴とする超音波診断装置。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、画像処理装置、画像処理プログラム、記憶媒体及び超音波診断装置に係り、特に、頸動脈の超音波画像の画像処理装置、画像処理プログラム、記憶媒体及び超音波診断装置に関する。
【背景技術】
【０００２】
動脈硬化症は、生活習慣の不摂生が危険因子となる生活習慣病の一つであり、自覚症状がほとんどないまま若年層から進行し、放っておくと脳梗塞、心筋梗塞、狭心症などに発展する可能性がある。動脈硬化症は治癒が困難であるため、早期に動脈硬化の兆候をとらえて患者に生活習慣の改善を促し、病気の進行を遅らせることが重要である。
【０００３】
動脈硬化症の検査の中で、超音波エコー検査は、非侵襲に血管を直接観察できることから、現在の臨床診断において必要不可欠な検査として定着しつつある。動脈硬化のための超音波エコー検査では、血管内部に付着したプラークと呼ばれる脂肪分の有無や内膜中膜複合体厚(IMT: intima-media complex thickness)を評価する。しかし、動脈硬化初期段階ではプラークもなく、IMTも薄くて計測が困難であるため、こうした形態学的観察では、初期の兆候をとらえることは難しい。
【０００４】
これに対し、動脈硬化は動脈が硬化し弾力性を失うことであるから、心拍動に伴う動脈の動きの良さを評価すればよいという考え方がある。プラーク、IMTなどの形態学的計測とは異なり、動脈の動きの定量的計測を行う場合、計算機の援用が不可欠である。
【０００５】
特許文献１では頸動脈径の心拍変動を測定することにより、動脈硬化の指標を得ることが提案されている。さらに、本出願人は特願２００７−２３９１９８号で、超音波Ｂモード長軸断面像を用いて頸動脈径の心拍変動を推定する方法を提案している。非特許文献１では、Ｂモード画像上で頸動脈壁の１点を手動で指定したのち、ＲＦ信号を用いて指定した頸動脈壁の動きを１０μｍの精度で追跡する手法を提案し実用化している。非特許文献１では、ＲＦ信号を用いているため、一部の超音波診断装置でしか利用できないが、現在のところ最も高い精度を有していると推察される。
【０００６】
なお、非特許文献２〜５は、出願時の技術水準が記載されている。
【特許文献１】特許第３８８２０８４号
【非特許文献１】原田烈光,「エコートラッキング法による動脈硬化指標の計測」，超音波検査技術，Vol.28(3):pp.29-35,2003.
【非特許文献２】ジーエッチ・グランルンド、「フーリエ  プレプロセシング  フォー  ハンド  プリント  キャラクタ  レコグニッシヨン(Fourier preprocessing for hand character recognition)」、IEEE Trans. on Computers, Vol. C-21, pp.195-201,1972.
【非特許文献３】クリス・エル・デ・コルテ他、「キャラクタリゼーション  オブ  プラーク  コンポーネンツ  ウイズ  イントラバスキュラ  ウルトラサウンド  エラストグラフィ  イン  ヒューマン  フェモラル  アンド  コロナリイ  アーテェリーズ  イン  ヴィトロ(Characterization of Plaque Components With Intravascular Ultrasound Elastgraphy in Human Femoral and Coronary Arteries In Vitro)」、
【非特許文献４】ラジ・エー・バルデュシング，クリス・エル・デ・コルテ他，「コンパリソン  オブ  ファイナイト  エレメンツ  モデル  エラストグラムス  アイヴィユーエス  エラストグラムス  アクワイヤド  フロム  ファントム  アンド  アーテェリーズ(Comparison of Finite Elements Model Elastograms and IVUS Elastograms acquied from Phantoms and Arteries)」、IEEE Symposium, pp.1873-1875,2002
【非特許文献５】林香苗、「日本人並びに日本産医学実験動物の解剖及び生理学係数」、第２版（訂正復刻版）、文京書院、1965,pp9
【発明の開示】
【発明が解決しようとする課題】
【０００７】
本出願人が提案した超音波Ｂモード長軸断面像を用いて頸動脈径の心拍変動を推定する方法は、頸動脈は長軸（ｙ軸）座標に関わらず長軸においては同質であり、等しい心拍変動を示すという前提条件で行う。多くの場合は、こうした条件で行うことができるが、局所的に動脈硬化が進行するケースもあり、このようなケースでは、前記前提条件を満たさないため、前述の方法を適用することはできない。なお、本明細書では、頸動脈が延びる方向を長軸とし、長軸に直交する平面（横断面）に含まれる軸を短軸という。
【０００８】
ところで、局所的に頸動脈に動脈硬化が生ずる場合は、超音波Ｂモード画像である頸動脈の短軸断面像を使用して頸動脈の心拍変動を推定することが考えられるが、超音波Ｂモード画像の１画素のサイズは、おおよそ０．０７ｍｍ〜０．１１ｍｍであり、十分に高い解像度とはいえない。そのため、１画素よりも細かい精度で心拍変動を捉える技術が不可欠であると考えられていた。
【０００９】
又、非特許文献１のようにＲＦ信号を用いる場合には、高いサンプリングレートで行って、位相差追跡法を用いれば高い推定精度を得ることができるが、Ｂモード画像には適用できない。さらに、非特許文献１では、頸動脈の動きの追跡では、精度の観点からＲＦ信号を用いることが有利であるものの、ＲＦ信号を出力可能な超音波診断装置は一部の機種に限られていることから、現時点では一般的とはいえず、汎用性に欠ける。
【００１０】
特願２００７−２３９１９８号では、ｙ軸に延びる頸動脈画像の空間的形状を利用して推定精度を改善している。しかし、短軸（ｘ−ｚ）断面では、下記の理由により形状情報を利用することが難しいため、これまでの方法では、十分な精度で頸動脈の心拍変動を追跡できない問題がある。１つ目は、頸動脈短軸断面形状はおおよそ円形であるが、個人差が大きいため、統一した形状では表現しにくい。２つ目は短軸断面では、超音波画像中の頸動脈が占める面積は、長軸断面の場合に比して小さいため、利用できる空間情報が少ない。
【００１１】
本発明の目的は、局所的に頸動脈に動脈硬化が生ずる場合にも、超音波Ｂモード画像を使用して十分な精度で頸動脈の心拍変動を追跡できる画像処理装置を提供することにある。
【００１２】
又、本発明の他の目的は、局所的に頸動脈に動脈硬化が生ずる場合にも、超音波Ｂモード画像を使用して十分な精度で頸動脈の心拍変動を追跡できる画像処理プログラムを提供することにある。
【００１３】
又、本発明の他の目的は、局所的に頸動脈に動脈硬化が生ずる場合にも、超音波Ｂモード画像を使用して十分な精度で頸動脈の心拍変動を追跡できる画像処理プログラムを記憶した記憶媒体を提供することにある。
【００１４】
さらに、本発明の他の目的は、上記画像処理装置を備えた超音波診断装置を提供することにある。
【課題を解決するための手段】
【００１５】
上記目的を達成するために、請求項１に記載の発明は、連続したフレームからなる超音波Ｂモード画像列を取得して、所定のフレームに描画される頸動脈短軸断面の頸動脈壁を推定する頸動脈壁推定手段と、推定された頸動脈壁及びその周辺組織を含む領域の画像をテンプレート画像とし、このテンプレート画像の頸動脈壁と前記周辺組織を含む領域の画像を、変形させて、変形後のテンプレート画像と、前記取得した各フレームの画像との誤差を最小にすることにより頸動脈壁の径の大きさを推定して頸動脈壁の時間変化を取得する心拍変動追跡手段を備えることを特徴とする画像処理装置を要旨とするものである。
【００１６】
請求項２の発明は、請求項１において、前記頸動脈壁推定手段は、頸動脈の位置と形状を決める頸動脈パラメータが含まれて頸動脈であるらしさを表わす評価関数を最大にすることにより、頸動脈壁を推定することを特徴とする。
【００１７】
請求項３の発明は、請求項１又は請求項２において、前記心拍変動追跡手段は、前記フレームが、フレーム毎に頸動脈の中心が変動している場合には、頸動脈の中心座標のずれも推定することを特徴とする。
【００１８】
請求項４の発明は、コンピュータを、連続したフレームからなる超音波Ｂモード画像列を取得して、所定のフレームに描画される頸動脈短軸断面の頸動脈壁を推定する頸動脈壁推定手段と、推定された頸動脈壁及びその周辺組織を含む領域の画像をテンプレート画像とし、このテンプレート画像の頸動脈壁と前記周辺組織を含む領域の画像を、変形させて、変形後のテンプレート画像と、前記取得した各フレームの画像との誤差を最小にすることにより頸動脈壁の径の大きさを推定して頸動脈壁の時間変化を取得する心拍変動追跡手段として機能させることを特徴とする画像処理プログラムを要旨とするものである。
【００１９】
請求項５の発明は、請求項４において、前記頸動脈壁推定手段が、頸動脈の位置と形状を決める頸動脈パラメータが含まれて頸動脈であるらしさを表わす評価関数を最大にすることにより、頸動脈壁を推定することを特徴とする。
【００２０】
請求項６の発明は、請求項４又は請求項５において、前記心拍変動追跡手段が、前記フレームが、フレーム毎に頸動脈の中心が変動している場合には、頸動脈の中心座標のずれも推定することを特徴とする。
【００２１】
請求項７の発明は、請求項４乃至請求項６のいずれか１項に記載の画像処理プログラムを記憶した記憶媒体を要旨とするものである。
請求項８の発明は、請求項１乃至請求項３のうちいずれか１項に記載の画像処理装置に接続され、探触子の超音波送受信により得られたエコー信号を基に生体内の超音波Ｂモード画像列を生成し、前記画像処理装置に前記超音波Ｂモード画像列を出力することを特徴とする超音波診断装置を要旨とするものである。
【発明の効果】
【００２２】
請求項１の発明によれば、局所的に頸動脈に動脈硬化が生ずる場合にも、超音波Ｂモード画像を使用して十分な精度で頸動脈の心拍変動を追跡できる効果を奏する。特に、超音波Ｂモード画像の１画素のサイズは、おおよそ０．０７ｍｍ〜０．１１ｍｍであり、十分に高い解像度とはいえないにも拘わらず、頸動脈とともにその周辺組織をも考慮することにより十分な精度で頸動脈の心拍変動を追跡できる意義はおおきいといえる。又、超音波Ｂモード画像を取得することができる一般的な超音波診断装置を使用することができるため、汎用性に富むことができる。
【００２３】
請求項２の発明によれば、頸動脈壁推定手段が、頸動脈の位置と形状を決める頸動脈パラメータが含まれて頸動脈であるらしさを表わす評価関数を最大にすることにより、頸動脈壁を容易に推定することができる。
【００２４】
請求項３の発明によれば、心拍変動追跡手段が、前記フレームが、フレーム毎に頸動脈の中心が変動している場合には、頸動脈の中心座標のずれも推定することにより、フレーム毎に中心が変動している場合においても、請求項２、又は請求項１の効果を容易に奏することができる。
【００２５】
請求項４の発明によれば、局所的に頸動脈に動脈硬化が生ずる場合にも、超音波Ｂモード画像を使用して十分な精度で頸動脈の心拍変動を追跡できる画像処理装置としてコンピュータを機能させることができる画像処理プログラムを提供できる。
【００２６】
請求項５の発明によれば、コンピュータを、頸動脈の位置と形状を決める頸動脈パラメータが含まれて頸動脈であるらしさを表わす評価関数を最大にすることにより、頸動脈壁を容易に推定することができる画像処理装置として機能させることができる画像処理プログラムを提供できる。
【００２７】
請求項６の発明によれば、コンピュータを、フレームが、フレーム毎に頸動脈の中心が変動している場合には、頸動脈の中心座標のずれも推定することにより、頸動脈の心拍変動を容易に追跡することができる画像処理装置として機能させることができる画像処理プログラムを提供できる。
【００２８】
請求項７の発明によれば、局所的に頸動脈に動脈硬化が生ずる場合にも、超音波Ｂモード画像を使用して十分な精度で頸動脈の心拍変動を追跡できる画像処理装置としてコンピュータを機能させることができる画像処理プログラムを記憶した記憶媒体を提供できる。
【００２９】
請求項８の発明によれば、請求項１乃至請求項３のうちいずれか１項に記載の効果を備えた超音波診断装置を提供できる。
【発明を実施するための最良の形態】
【００３０】
以下、本発明を具体化した一実施形態を図１〜１３を参照して説明する。
図１に示すように、本実施形態の画像処理装置１１は、頸動脈壁推定手段、心拍変動追跡手段としてのコンピュータ１２、ディスプレイ１３、プリンタ１４、キーボード１５、及び超音波診断装置１６を備えている。
【００３１】
前記超音波診断装置１６は、探触子２１、表示部２２及び図示しない超音波発生源が備えられている。超音波診断装置１６は、探触子２１の先端面から生体内にパルス波を送波するとともに頸動脈２３からの反射波（エコー信号）を受波する。この超音波診断装置１６は、前記探触子２１の先端面が長四角形状に形成されており、前記エコー信号に基づいて頸動脈２３のＢモードで取得された短軸断面形状（すなわち、頸動脈短軸断面形状）を含むエコー動画像（すなわち、Ｂモード画像）を生成する。このエコー動画像は、超音波Ｂモード画像列に相当する。超音波診断装置１６の表示部２２には、被検者の頸動脈２３が拡張収縮変形する際の動画像が表示される。なお、図５には、頸動脈短軸断面での座標系を示している。
【００３２】
超音波診断装置１６のビデオ出力から前記Ｂモード画像が出力され、図示しないＤＶコンバータにより所定の解像度（ｍ×ｎ画素）、所定数のフレーム／秒でビデオキャプチャされた後、ＩＥＥＥ１３９４信号として、コンピュータ１２に入力される。所定数のフレーム／秒は、本実施形態では３０フレーム／秒としているが、この値が限定されるものではなく、頸動脈の心拍動の変化を明確に観測できる値であればよい。
【００３３】
また、ディスプレイ１３で表示される超音波画像領域の画素数は、水平方向にＮｘ、垂直方向にＮｙである。本実施形態では、ディスプレイ１３で表示される超音波画像領域の画素数は、Ｎｘ＝５４０（pixel）,Ｎｙ＝４２０（pixel）としているが、この数値に限定されるものではない。
【００３４】
前記超音波診断装置１６により取得されたエコー動画像がコンピュータ１２に入力されたコンピュータ１２が備える記憶装置４４に格納される。前記記憶装置４４は、例えばハードディスクや半導体記憶装置からなり、前述した前記各種情報の読み出し及び書き込みが可能である。
【００３５】
図２に示すように、コンピュータ１２はＣＰＵ４１（中央処理装置）、ＲＯＭ４２、ＲＡＭ４３を備え、記憶媒体としてのＲＯＭ４２に格納された画像処理プログラムを実行する。ＲＡＭ４３は、前記プログラムを実行する際の作業用メモリである。
【００３６】
（実施形態の作用）
次に上記のように構成された画像処理装置１１において、ＣＰＵ４１が実行する画像処理プログラムの処理を説明するが、このプログラム処理の説明の前に「頸動脈、及び周辺組織の応力・歪み特性、二層円筒モデル、並びに心拍変動追跡法の概略」について説明する。
【００３７】
（頸動脈、及び周辺組織の応力・歪み特性）
頸動脈は心拍動により拡張・収縮を繰り返すが、それに伴い頸動脈の周辺にある脂肪などの周辺組織もまた押し縮められ、復元する。一般にこうした物質の応力と歪みの関係は材料力学の立場から示される。
【００３８】
図３に、超音波Ｂモード画像における頸動脈の短軸断面画像（以下、単に断面画像又は単に画像という）が示されている。図３において、ほぼ中央の円形のものが頸動脈５０であり、黒い部分（図においてはハッチング部分）が頸動脈内腔５０ａ、周辺の白い部分環状部分が頸動脈の血管壁（以下では頸動脈壁５０ｂという）であり、その周囲に脂肪などの周辺組織がある。
【００３９】
（二層円筒モデル）
ここで頸動脈と周辺組織を、図４に示すように一様な二層円筒（材料力学的には二層厚肉円筒）でモデル化する。頸動脈壁の内径（内半径）をｒ１とし、外半径（外径）をｒ２とする。周辺組織は、ｒ２から無限遠まで存在するものとする。頸動脈壁の弾性係数をＥ１とする。頸動脈壁の外側には、ｒ１，ｒ２に比べ十分大きな半径まで脂肪などの周辺組織があり、その周辺組織の弾性係数をＥ２とする。また、頸動脈壁、周辺組織のポアソン比をいずれも０．５とする。一般的には、この半径ｒ１，ｒ２は、内圧０の時を考えるが、生体の場合、内圧すなわち血圧が０の状態を考えられないため、ある適当な内圧（血圧）がかかっている場合を基準内圧とする。この基準内圧からさらに内圧（血圧）Ｐが頸動脈に加えられた際の半径ｒでの半径方向の変位ｄ（ｒ）は、二層円筒モデルの理論式（１）で表わすことができる。
【００４０】
【数１】
0004879872.tif 000002
すなわち、半径ｒにある組織は、内圧Ｐにより、半径ｒ＋ｄ（ｒ）に移動する。内圧Ｐが正の場合、変位ｄは正、つまり拡張し、内圧Ｐが負の場合、変位ｄは負、すなわち収縮を表わす。
【００４１】
式（１）にｒ＝ｒ１を代入すれば、頸動脈壁の内面（以下、頸動脈内壁という）での変位は、
【００４２】
【数２】
0004879872.tif 000003
となる。ここで、生理学的見地から、頸動脈の弾性係数Ｅ１と周辺組織Ｅ２の弾性係数は、それほど大きく違わないことが非特許文献３，４から分かっている。また、頸動脈壁の厚みW=ｒ２−ｒ１は、２０，３０歳代で約0.73mm、４０歳代で約0.79mm、５０歳代で約0.98mm、６０歳代で約1.00mm、７０歳代で約1.04mm程度であることが、非特許文献５で分かっており、頸動脈の直径（6-12mm）や周辺組織の厚みに比べて十分に薄いといえる。以上のことから、式（１）は、
【００４３】
【数３】
0004879872.tif 000004
と近似できる。
【００４４】
この式（３）は、心収縮、すなわち、血圧上昇により、頸動脈は拡張するが、周辺組織もまた拡張し、その拡張の度合い、つまり変位は、頸動脈中心からの距離ｒに反比例することを意味している。この一例として、Ｅ１＝200kPa、Ｅ２＝100kPa、ｒ１＝4.0mm、ｒ２＝4.8mm、Ｐ＝40mmHg（＝5.3kPa）の時の半径ｒと変位ｄ（r; P, r1, r2, E1, E2）の関係を図６に示す。この条件において、基準血圧からの内圧をＰ＝−10,0,20,40mmHgとした際の「半径−変位曲線」の一例を図７に示す。同図に示すように内径ｒと変位ｄとの関係は、内圧Ｐ（血圧）に関わりなく、反比例で表わされる。
【００４５】
（心拍変動追跡法の概略）
ここで、最も重要なことは式（３）は反比例を表わす式であるため、図６，７のグラフ中、例えば、頸動脈壁の内径ｒ１とその点での変位ｄ１が与えられれば、周辺組織も含めて任意の半径ｒでの変位は完全に決定できる。すなわち、１フレームの断面画像Img１（テンプレート画像）があったとき、その画像における頸動脈壁の内径ｒ１が分かっていれば、その画像を用いて式（３）のように変形させれば、その他の任意の血圧時における断面画像Img２（変形後のテンプレート画像）を模擬的に作成できる。その際、拡張、或いは収縮の度合いを制御するパラメータはｄ１のみである。
【００４６】
断面画像Img１からｄ１により模擬的に作成した画像（変形後のテンプレート画像）と、計測された断面画像Img２の差を最小にするようにパラメータｄ１を決定すれば、それは、元の断面画像Img１からの断面画像Img２における変位ｄ１の推定値となっている。すなわち、断面画像Img１では頸動脈壁の内径がｒ１であったが、計測された断面画像Img２では頸動脈壁の内径がｒ１＋ｄ１になったと解釈できる。
【００４７】
これを計測されたすべての断面画像について行えば、頸動脈壁の内径の時間変化を追跡できたことになる。
この心拍変動追跡法のポイントは、頸動脈壁の心拍変動を追跡するために、頸動脈の周辺組織の心拍動による変位を利用していることである。頸動脈壁の動きを追跡することは、線の追跡であるのに対し、こうした方法では、面の追跡となっており、追跡の精度を大幅に改善でき、さらに、利用できる面（すなわち、周辺組織）の面積が大きければ、それだけ精度改善ができる。
【００４８】
（画像処理プログラムの処理）
次に画像処理装置１１において、ＣＰＵ４１が実行する画像処理プログラムの処理を説明する。ＣＰＵ４１は、図１０のＳ１０〜Ｓ４０の処理を行うことにより、頸動脈の心拍変動追跡を行う。又、心拍変動追跡法の概念を図１２に示す。
【００４９】
Ｓ１０では、ＣＰＵ４１は、断面画像列（すなわち、動画像）を取得して計測する。ここで、取得された超音波Ｂモード画像列をＩ（ｘ，ｙ，ｔ），ｘ＝１，．．．，Ｎｘ，ｙ＝１，．．．，Ｎｙ，ｔ＝１，．．．，Ｎｔとする。ｘ，ｙ，ｔはそれぞれ水平方向の画素番号、垂直方向の画素番号、フレーム番号を表わす。たたじ、画像の左上をｘ＝１，ｙ＝１とする。また、Ｎｘ，Ｎｙ，Ｎｔは水平方向の画素総数、垂直方向の画素総数、フレーム総数を表わす。１枚の画像に着目する場合は、Ｉ（ｘ，ｙ）、ｘ＝１，．．．，Ｎｘ，ｙ＝１，．．．，Ｎｙと書く。
【００５０】
Ｓ２０では、断面画像列の中の適当な１つのフレーム（所定フレーム）の画像Img1について、頸動脈内壁の位置と形状ｒ１（θ）を推定する（図９参照）。なお、詳細は後述する。
【００５１】
Ｓ３０では、ＣＰＵ４１は、所定のフレームの画像を式（３）により心拍変動を与えた画像とその他のフレームの画像の誤差を最小にするように、変位ｄ１を推定する。そして、上記の処理を全てのフレームについて行い、変位の時間変化ｄ１（ｔ）を得る。なお、詳細は後述する。
【００５２】
Ｓ４０では、頸動脈形状、つまり、頸動脈壁の内径の時間変化を、ｒ１（θ）＋ｄ１（ｔ）としてディスプレイ１３やプリンタ１４の出力装置に出力する。
（一つのフレームの画像での頸動脈壁の内径ｒ１(θ)の推定；Ｓ２０の詳細）
前記Ｓ２０においては、ＣＰＵ４１は、画像列中の適当な一つのフレームの画像Ｉmg1を選択し、その画像中の頸動脈壁の内径ｒ１(θ)を推定する。この選択されるフレームは、限定されるものではなく、無作為に選択しても、最初のフレームでもよい。或いは、この選択されるフレームは、頸動脈が最も鮮明に映っているフレーム、頸動脈収縮期にあるフレーム、頸動脈拡張期にあるフレームとしてもよい。
【００５３】
頸動脈内壁は、真円、もしくは、多少歪んだ円である。そこで、頸動脈内壁の内径ｒ１（θ）、すなわち、頸動脈短軸断面での頸動脈内壁形状は、閉曲線を表現するためのフーリェ記述子でモデル化する。なお、フーリェ記述子については非特許文献２に記載されている。
【００５４】
フーリェ記述子で頸動脈内壁形状、すなわち頸動脈壁をモデル化するにあたり、画像の座標ｘ，ｙを複素数ｘ＋ｉｙで表わすものとする。ただし、ｉは虚数単位を表わす。フーリェ記述子ｆは媒介変数ｓを用いて、複素平面上の閉曲線を
【００５５】
【数４】
0004879872.tif 000005
として表現する。ｃｋはフーリェ記述子と呼ばれ、閉曲線の位置や形状を決めるパラメータである。また、整数Ｋを大きな値にするほど入り組んだ複雑な閉曲線を表現できるようになる。なお、頸動脈の場合、Ｋ＝５程度とすればよいが、この値に限定されるものではない。これらのフーリェ記述子により頸動脈の形状が規定されるから、これらのフーリェ記述子をまとめて頸動脈パラメータp＝{ck|k=0,±1,±2,....,±K}と呼ぶことにする。
【００５６】
頸動脈内腔は、閉曲線ｆ（ｓ）で囲まれる領域であり、これらの領域に含まれる画素（ｘ，ｙ）の集合をＴとする。また、頸動脈壁の厚みをＷとして、Ｔから距離Ｗ以下の距離にある画素（ｘ，ｙ）の集合を頸動脈壁を表わす画素（ｘ，ｙ）の集合と考え、これを集合Ｂとおく。集合Ｂは次式で与えられる。
【００５７】
【数５】
0004879872.tif 000006
画像において、頸動脈壁は強いエコーのために白く、頸動脈内腔はエコーが弱いため黒く描かれる。すなわち、頸動脈壁を表わすＢでは画素値が大きく、頸動脈壁内腔を表わすＴでは画素値は小さくなる。そこで、位置・形状を与える頸動脈パラメータｐに対し、それが頸動脈であるらしさを表わす以下の評価関数を導入する。
【００５８】
【数６】
0004879872.tif 000007
ただし、α、βは、いずれも正値の定数であり、|・|は、集合の要素数を表わす。右辺第１項は、頸動脈壁の想定する位置で画素値が大きい、つまり白いほど大きな値をとり、第２項は、頸動脈内腔と想定される位置で画素値が小さい、つまり黒いほど大きな値をとる。第３項において、ｇ（Ａ；Ａmin,Ａmax）は次式で表わされる関数である。
【００５９】
【数７】
0004879872.tif 000008
ここで、Ａは内腔部、すなわち、集合Ｔの面積を表わし、Ａmin，Ａmaxは、それぞれ頸動脈内腔面積の考え得る最小値と最大値である。式（６）の第３項は、その範囲を外れた場合のペナルティを表わす。Ａmin，Ａmaxは収集した画像サンプルに基づいて予め固定値として入力して設定される。例えば、Ａmin＝２０(mm２)，Ａmax＝１１０(mm２)の値等が設定されるが、この値に限定されるものではない。また、βはコンピュータにおいて、丸め誤差が生じない程度に十分な大きな正数を選択すればよい。式（６）において、αは頸動脈壁と頸動脈内腔でのエコーの強さの比、超音波診断装置１６の画像出力などの仕様に依存するが、予備試験によりα＝２．０程度を選択すればよい。
【００６０】
なお、式（６）は、重み
【００６１】
【数８】
0004879872.tif 000009
を定義すれば、
【００６２】
【数９】
0004879872.tif 000010
と書くこともできる。
【００６３】
この重み関数ｗ（ｘ，ｙ；ｐ）の例を図８に示す。
画像Ｉ（ｘ，ｙ）における頸動脈の位置と形状、すなわち頸動脈パラメータｐの値は、式（６）、若しくは式（９）で与えられるＥ（ｐ）を最大にするように決定する。こうしたパラメータの値は解析的には求められないため、最急降下法などの非線形最適化法を利用して求める。なお、画像中には、頸動脈に類似した円形の組織はないため、局所解に陥ることは少ない。
【００６４】
ＣＰＵ４１は、得られた頸動脈内壁の輪郭線、すなわち、閉曲線ｆ（ｓ），０≦ｓ≦１から頸動脈壁のθ毎の内径ｒ１（θ）は以下のように求める。まず、頸動脈の中心座標は、０次フーリェ記述子c0を用いて（Ｒｅ（c0），Ｉｍ（c0））で与えられる。ただし、Ｒｅ（・），Ｉｍ（・）はそれぞれ「・」の実部、虚部を表わす。従って、頸動脈の中心から見た輪郭線はｆ（ｓ）−c0，０≦ｓ≦１となるから、ある媒介変数ｓでの値での動径ｒ（ｓ）と偏角θ（ｓ）は、それぞれ
【００６５】
【数１０】
0004879872.tif 000011
【００６６】
【数１１】
0004879872.tif 000012
となる。そして、媒介変数ｓを外し、ｒ（ｓ）をθ（ｓ）の関数と見たものが頸動脈内壁の形状の極座標表示ｒ１（θ），０≦θ≦２πとなる。
【００６７】
ＣＰＵ４１により推定された頸動脈壁（頸動脈内壁）の輪郭線、すなわち、内径ｒ１（θ）の一例を図１１に示す。図１１において、Ｑが頸動脈の輪郭線を示している。
（他のフレームの画像での頸動脈の変位ｄ１の推定；Ｓ３０の詳細）
上記のように推定された頸動脈内壁とその周辺組織を適当に含む領域をＣＰＵ４１は、切り出す。すなわち、上記フレームの画像Img1において、推定された頸動脈の中心を原点として極座標表示されたものをテンプレート画像template(r,θ)とする。
【００６８】
又、ＣＰＵ４１は、任意のフレーム（ｔ番目のフレーム）の画像について、前記Ｓ２０で推定された頸動脈の中心を原点として極座標表示したものをＩ（ｒ，θ，ｔ）で表わす。
【００６９】
原理的には、任意のフレームでの頸動脈、及び周辺組織を含む画像Ｉ（ｒ，θ，ｔ）は、テンプレート画像を適当なｄ１により変形、すなわち、拡張、或いは収縮させて
【００７０】
【数１２】
0004879872.tif 000013
と近似できる。この近似の２乗誤差
【００７１】
【数１３】
0004879872.tif 000014
を最小にするように、すなわち、
【００７２】
【数１４】
0004879872.tif 000015
として、ＣＰＵ４１は、ｄ１を決定する。このｄ１が時刻ｔにおける変位となる（図９参照）。
【００７３】
ところで、ｒは２乗誤差を評価するために利用される領域の半径であり、ｒを大きくすると、すなわち、周辺組織をより多く含むように切り出すと、頸動脈形状の心拍変動の推定に利用できる領域が増え、推定精度も高くなることが期待できる。しかし、一方、むやみに周辺組織をより多く含むように切り出すと、周辺組織の非一様性などの影響もでてきて、逆に推定精度を悪化させることになる。この場合、頸動脈壁の内径ｒ１（θ）の２〜３倍程度の半径の周辺組織を含むように切り出すのが好ましい。
【００７４】
ここで、一例として、あるフレームの画像について、テンプレート画像との２乗誤差を算出した結果が図１２に示されている。ただし、ｒ１は、平均の半径
【００７５】
【数１５】
0004879872.tif 000016
を表わす。この図１２では、変化率ｄ１／ｒ１＝0.08で２乗誤差Ｊが最小になっており、この値が変位の推定値となる。図１２では、変化率と２乗誤差Ｊの関係を示している。
【００７６】
こうした処理をＣＰＵ４１は、テンプレート画像に使用したフレームを除いて、取得したすべてのフレームについて行うことにより、変位の時間変化ｄ１（ｔ）を得る。この結果、頸動脈壁の角θ方向の半径変化は、ｒ１（θ）＋ｄ１（ｔ）となり、ＣＰＵ４１により、頸動脈形状の心拍変動の追跡が行われたことになる。
【００７７】
（頸動脈の中心が変動している場合）
ＣＰＵ４１は、頸動脈の中心位置が変動している場合には、ｄ１だけでなく、中心座標のずれ(x0,y0)も合わせて推定する。
【００７８】
この場合、Ｉ（ｒ，θ，ｔ）の中心座標を(x0,y0)だけシフトさせた画像をＩ’（ｒ，θ，ｔ）とし、式（１３），式（１４）の代わりに
【００７９】
【数１６】
0004879872.tif 000017
【００８０】
【数１７】
0004879872.tif 000018
とする。
【００８１】
また、元のテンプレート画像template(r,θ)から心拍変動を与えた画像template(r＋ｄ（ｒ；ｒ１（θ），ｄ１）,θ)を作成する際、標本点（画素）の位置は、画像を伸縮させるため、元のテンプレート画像とは異なる。このままでは、計測画素Ｉ（ｒ，θ，ｔ）との誤差を評価できないため、標本点を一致させるリ・サンプリングを行う必要がある。上記の説明では、２次元線形補間を利用しているが、その他、モーフィングなどで利用される様々な補間（内挿）技術、例えば、２次元キュービック補間、２次元スプライン補間などの補間（内挿）法を利用できる。
【００８２】
（推定結果）
２０歳と６０歳の健常男性に対する頸動脈超音波Ｂモード画像列の動画像を取得して、それらの頸動脈形状の心拍変動（θ＝０度方向の半径変化ｒ１（θ）＋ｄ１（ｔ）），０≦ｔ≦４（秒）の推定結果を図１３（ａ）〜（ｌ）に示す。図１３（ａ）〜（ｆ）は２０歳健常男性、図１３（ｇ）〜（ｌ）は６０歳健常男性に対する頸動脈形状の心拍変動の推定結果である。図１３（ａ）〜（ｌ）において、縦軸は頸動脈の半径（mm）であり、横軸は時間（ｓ）である。図１３（ａ）〜（ｆ）及び図１３（ｇ）〜（ｌ）のそれぞれにおいて順に末端側（脳側）から心臓側へ５mm間隔毎に頸動脈短軸断面を移動させて頸動脈の心拍変動を推定した結果が示されている。なお、図１３（ｇ）〜（ｌ）の縦軸のスケールは、図１３（ａ）〜（ｌ）の縦軸のスケールに対して２倍ほど異なる。
【００８３】
２０歳男性に対する頸動脈の心拍変動は、大きくはっきりとしている。一方、６０歳男性に対する頸動脈の心拍変動は、小さく、周期性も乏しい。又、２０歳男性においてはどの断面でも同じような心拍変動を示しているが、図１３（ｈ）及び図１３（ｋ）に示すように６０歳男性での頸動脈の短軸断面では、その他の部位に比べ心拍変動がやや小さく、局所的に動脈硬化が進んでいることが予想される。
【００８４】
本実施形態によって発揮される効果について、以下に記載する。
（１）  本実施形態の画像処理装置１１は、コンピュータ１２が頸動脈壁推定手段として、連続したフレームからなる超音波Ｂモード画像列を取得して、所定のフレームに描画される頸動脈短軸断面の頸動脈壁を推定するようにしている。そして、コンピュータ１２が心拍変動追跡手段として推定された頸動脈壁及びその周辺組織を含む領域のテンプレート画像の頸動脈壁と前記周辺組織を含む領域の画像を変形させ、変形後のテンプレート画像と、取得した各フレームの画像との誤差を最小にすることにより頸動脈壁の径の大きさを推定して頸動脈壁の時間変化を取得するようにしている。
【００８５】
この結果、本実施形態では、局所的に頸動脈に動脈硬化が生ずる場合にも、超音波Ｂモード画像を使用して十分な精度で頸動脈の心拍変動を追跡できる。特に、超音波Ｂモード画像の１画素のサイズは、おおよそ０．０７ｍｍ〜０．１１ｍｍであり、十分に高い解像度とはいえないにも拘わらず、頸動脈とともにその周辺組織をも考慮することにより十分な精度で頸動脈の心拍変動を追跡できる。さらに、超音波Ｂモード画像を取得することができる一般的な超音波診断装置を使用することができるため、汎用性に富むことができる。
【００８６】
（２）  また、本実施形態の画像処理装置１１は、コンピュータ１２が、頸動脈壁推定手段として、頸動脈の位置と形状を決める頸動脈パラメータが含まれて頸動脈であるらしさを表わす評価関数を最大にすることにより、頸動脈壁を推定する。この結果、頸動脈壁を容易に推定することができる。
【００８７】
（３）  本実施形態の画像処理装置１１は、コンピュータ１２が、心拍変動追跡手段として、フレームが、フレーム毎に頸動脈の中心が変動している場合においても、頸動脈の中心座標のずれも推定するようにした。この結果、フレームが、フレーム毎に頸動脈の中心が変動している場合においても、頸動脈壁を容易に推定することができる。
【００８８】
（４）  本実施形態の画像処理プログラムは、コンピュータ１２を、頸動脈壁推定手段として、連続したフレームからなる超音波Ｂモード画像列を取得して、所定のフレームに描画される頸動脈短軸断面の頸動脈壁を推定させる。又、画像処理プログラムは、コンピュータ１２を、心拍変動追跡手段として、前記推定された頸動脈壁及びその周辺組織を含む領域の画像をテンプレート画像とし、同画像の頸動脈壁と前記周辺組織を含む領域の画像を変形させる。そして、画像処理プログラムは、コンピュータ１２を、心拍変動追跡手段として、変形後のテンプレート画像と、取得した各フレームの画像との誤差を最小にすることにより頸動脈壁の径の大きさを推定して頸動脈壁の時間変化を取得させる。この結果、局所的に頸動脈に動脈硬化が生ずる場合にも、超音波Ｂモード画像を使用して十分な精度で頸動脈の心拍変動を追跡できる画像処理装置としてコンピュータを機能させることができる画像処理プログラムを提供できる。
【００８９】
（５）  本実施形態の画像処理プログラムは、コンピュータ１２を、頸動脈壁推定手段として、頸動脈の位置と形状を決める頸動脈パラメータが含まれて頸動脈であるらしさを表わす評価関数を最大にすることにより、頸動脈壁を推定させる。この結果、頸動脈壁を容易に推定することができる画像処理プログラムを提供できる。
【００９０】
（６）  本実施形態の画像処理プログラムは、コンピュータ１２を、心拍変動追跡手段として、フレームが、フレーム毎に頸動脈の中心が変動している場合には、頸動脈の中心座標のずれも推定させる。この結果、頸動脈の心拍変動を容易に追跡することができる画像処理装置として機能させることができる画像処理プログラムを提供できる。
【００９１】
（７）  本実施形態のＲＯＭ４２は、前記画像処理プログラムを記憶した記憶媒体とした。この結果、局所的に頸動脈に動脈硬化が生ずる場合にも、超音波Ｂモード画像を使用して十分な精度で頸動脈の心拍変動を追跡できる画像処理プログラムを記憶した記憶媒体を提供できる。
【００９２】
（８）  本実施形態の超音波診断装置１６は、画像処理装置１１に接続されているため、上記（１）乃至（３）の効果を容易に実現できる超音波診断装置となる。
なお、前記実施形態を次のように変更して構成することもできる。
【００９３】
○  前記実施形態では、ＲＯＭ４２を画像処理プログラムを記憶する記憶媒体としたが、記憶装置４４を画像処理プログラムを記憶する記憶媒体にしてもよく、或いは、ＣＤや、ＤＶＤ、ＵＳＢメモリ等の記憶媒体に画像処理プログラムを格納して、ＣＤやＤＶＤ等の記憶媒体を読込みするドライバ装置等にて読込みして実行するようにしてもよい。
【００９４】
○  前記動脈硬化解析結果は、ディスプレイ１３またはプリンタ１４のいずれか一方のみから出力されているように構成してもよい。
【図面の簡単な説明】
【００９５】
【図１】本発明を具体化した一実施形態の画像処理装置１１の概略図。
【図２】画像処理装置１１の概略電気ブロック図。
【図３】超音波Ｂモード画像における頸動脈の短軸断面画像。
【図４】二層円筒モデルの説明図。
【図５】頸動脈短軸断面での座標系の説明図。
【図６】半径ｒと変位d(r; P,r1,r2,E1,E2)の関係を示すグラフ。
【図７】基準血圧との差を−10,0,20,40mmHgとした際の半径−変位曲線の説明図。
【図８】重み関数ｗの説明図。
【図９】心拍変動追跡法の概念図。
【図１０】コンピュータ１２が処理する画像処理プログラムのフローチャート。
【図１１】ＣＰＵ４１により推定された頸動脈壁（頸動脈内壁）の輪郭線Ｑの説明図。
【図１２】変化率ｄ１／ｒ１と変化率と２乗誤差Ｊの関係を示すグラフ。
【図１３】（ａ）〜（ｆ）は２０歳の健常男性に対する頸動脈に対する推定結果の説明図、（ｇ）〜（ｌ）は６０歳の健常男性に対する頸動脈に対する推定結果の説明図。
【符号の説明】
【００９６】
１１…画像処理装置、
１２…コンピュータ（頸動脈壁推定手段、心拍変動追跡手段）、
１３…ディスプレイ、１６…医用超音波装置、
２１…探触子、４２…ＲＯＭ（記憶媒体）。
【図１】
0004879872.tif 000019
【図２】
0004879872.tif 000020
【図３】
0004879872.tif 000021
【図４】
0004879872.tif 000022
【図５】
0004879872.tif 000023
【図６】
0004879872.tif 000024
【図７】
0004879872.tif 000025
【図８】
0004879872.tif 000026
【図９】
0004879872.tif 000027
【図１０】
0004879872.tif 000028
【図１１】
0004879872.tif 000029
【図１２】
0004879872.tif 000030
【図１３】
0004879872.tif 000031