アテロームプラーク破裂のリスク評価のための画像処理方法
【課題】アテロームプラークの弾性勾配の正確な弾性画像、特にプラークが複数の互いに隣接する細胞の壊死体を有し、該壊死体にカルシウム含有物が含まれる可能性がある場合に、複数の構成成分における弾性を示す弾性画像を提供する。
【解決手段】本発明は‐血管組織が血圧機能により圧迫されることにより生じる内部変形を示す「エラストグラム」と呼ばれる受像工程と、事前に領域分割された複数の区域を有する画像を取得するために、前記エラストグラムを事前に領域分割する工程と、前記エラストグラムの少なくとも１の領域における弾性を示す弾性画像を演算する工程であって、前記領域は前記事前分割された画像における複数区域から選択された区域に対応し、前記弾性画像を用いて使用者はアテロームプラークの破裂をリスク評価することが可能となる前記弾性画像の演算工程とをそなえる。
【解決手段】本発明は‐血管組織が血圧機能により圧迫されることにより生じる内部変形を示す「エラストグラム」と呼ばれる受像工程と、事前に領域分割された複数の区域を有する画像を取得するために、前記エラストグラムを事前に領域分割する工程と、前記エラストグラムの少なくとも１の領域における弾性を示す弾性画像を演算する工程であって、前記領域は前記事前分割された画像における複数区域から選択された区域に対応し、前記弾性画像を用いて使用者はアテロームプラークの破裂をリスク評価することが可能となる前記弾性画像の演算工程とをそなえる。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、技術分野としては医療用画像処理に関し、さらに詳しくは血管組織のアテロームプラーク（血管内の粥状硬化部（脂質コア））破裂のリスク評価のための画像処理方法に関する。
【背景技術】
【０００２】
年々、心臓血管疾患による死者は癌による死者よりも多く発生している。狭心症、心筋梗塞、突然死及び脳卒中はいずれもアテローム性動脈硬化の影響によるものである。
【０００３】
組織学的研究によると、脆弱なアテロームプラークは壊死した大きな細胞外体１０として定義される。該細胞外体１０には、マクロファージ浸潤を伴う薄膜（例として膜厚１００μｍ以下）の線維性被膜２０が付随形成される。
【０００４】
心筋を灌漑する冠動脈３０は、アテロームプラークの突発的破裂により生じた血栓に狭窄され、したがって前記心筋における酸素が供給されない領域は損傷を受ける。この作用は突然死、不安定狭心症あるいは心筋梗塞の原因となる。
【０００５】
上記のような脆弱なプラークは、血管内超音波法（ＩＶＵＳ：Ｉｎｔｒａｖａｓｃｕｌａｒ Ｕｌｔｒａｓｏｕｎｄ）、光コヒーレンス・トモグラフィー（ＯＣＴ：Ｏｐｔｉｃａｌ Ｃｏｈｅｒｅｎｃｅ Ｔｏｍｏｇｒａｐｈｙ）、ＣＴスキャン及び核磁気共鳴画像法（ＭＲＩ：Ｍａｇｎｅｔｉｃ Ｒｅｓｏｎａｎｃｅ Ｉｍａｇｉｎｇ）を含む多様な技術により、臨床的に検出可能であることが、数件の報告によって明らかにされている。
【０００６】
しかしながら、形態及びプラーク組織の解析に基づく破裂予測は、極めて不正確なばかりでなく、十分な数のリスク指標が得られないという点においても好ましくない。
【０００７】
冠動脈プラークの破裂予測に関する目下の課題解決には、
‐線維層厚の正確な定量化及び
‐壊死体の形態解析
が必要となるだけでなく、特に前記線維層厚及び壊死体の力学特性に関する正確な知識、さらには動脈壁及びプラーク構成成分の弾性に関する正確な知識が求められる。
【０００８】
上記データ項目はプラーク破裂に関する好ましい生体力学的リスク指標となり得る基準であって、該基準により線維層に付加されるピーク応力を正確かつ効率的に評価することが可能となる。
【０００９】
しかしながら、生体内におけるプラークは不均質であり、さらにプラークは成長の過程で力学特性の変化を伴うため、プラークの力学的特性をｉｎ ｖｉｖｏ（生体内）で認識することは困難である。
【００１０】
さらに詳しくは、弾性（あるいは剛性）の画像化、すなわち血管組織の弾性を解析して画像化することは、応力分布を確実に算出するための必要前提条件とみなされる。
【００１１】
前記弾性画像の演算は、これまで様々なアプローチにより試みられてきた。
【００１２】
１心周期におけるプラーク負荷を入力データとして用い、前記血管弾性の画像を判断する研究が散見されている。この判断方法に関し、直接のアプローチあるいは相互作用手順が提案されている。
【００１３】
前記相互作用手順は、実測された負荷と、数値モデルを用いて演算された負荷との間の誤差を最小限に抑えることを目的とした最適化アルゴリズムを用いるものである。この文脈において、プラーク弾性の再現性向上は、前記最適化アルゴリズムの性能に依存している。
【００１４】
例えば、プラーク構成成分における弾性モジュライを抽出するための最適化アルゴリズム及び領域分割アルゴリズムが、いくつかの研究団体によって開発されている。しかしながらこれらの方法は、アテロームプラークが互いに隣接する複数の細胞の壊死体を含むという頻発事例において、弾性画像を取得するに十分な有効性あるいは性能レベルに達していない。
【００１５】
さらに加えて、これらの方法は以下に述べる部分においては十分に機能しない。第１に、該方法ではカルシウム含有物を抽出することが出来ないため、プラーク構成成分、特に剛体要素の特定が困難である。第２に、線維性被膜の厚さを正確に定量化することが困難であるため、プラークの脆弱性を正確に分析することに支障が生じる。
【００１６】
したがって、アテロームプラークの脆弱性の度合いを正確に分析するために、施術者は以下の項目を正確に把握している必要がある。
‐線維性被膜の厚さ
‐プラークが含有する様々な構成成分の大きさ
‐上記すべての成分における力学特性
【発明の概要】
【発明が解決しようとする課題】
【００１７】
本発明の目的は、アテロームプラークの形態上の詳細だけでなく、弾性勾配の正確な弾性画像、特にプラークが複数の互いに隣接する細胞の壊死体を有し、該壊死体にカルシウム含有物が含まれる可能性がある場合に、複数の構成成分における弾性を示す弾性画像を提供することにより、施術者にアテロームプラークの脆弱性の度合いを正確に評価させることである。
【００１８】
本発明の方法を用いれば、例えば、線維性被膜の厚さ及び壊死体の大きさ（厚さ）のような、アテロームプラークの破裂に関するリスク評価のための他の正確な基準を施術者に提供することも可能である。
【課題を解決するための手段】
【００１９】
上述の目的において、アテロームプラーク破裂のリスク評価のための画像処理方法を提案する。該方法は、
‐血管組織が血圧機能により圧迫されることにより生じる内部変形を示す「エラストグラム（ｅｌａｓｔｏｇｒａｍ）」と呼ばれる受像工程と、
‐事前に領域分割された（ｐｒｅ−ｓｅｇｍｅｎｔｅｄ）複数の区域を有する画像を取得するために、前記エラストグラムを事前に領域分割する（ｐｒｅ−ｓｅｇｍｅｎｔａｔｉｏｎ）工程と、
‐前記エラストグラムの少なくとも１の領域における弾性を示す弾性画像を演算する工程であって、前記領域は前記事前分割された画像における複数区域から選択された区域に対応し、前記弾性画像を用いて使用者はアテロームプラークの破裂をリスク評価することが可能となる前記弾性画像の演算工程と
をそなえた方法であって、
‐第１に、前記血管組織が非圧縮性の弾性固体であるとみなして事前領域分割を行い、ラグランジュ乗数の勾配及び空間依存性ベクトルの勾配の合計と等しくなるように、前記血管組織のヤング率（弾性係数）の勾配を決定し、
‐第２に、前記事前領域分割は、前記複数の区域（半径方向に分割した区域）における夫々の区域の輪郭を前記空間依存性ベクトルにより決定する工程を含む。
【００２０】
本発明の方法の好ましい特徴を以下に述べるが、該特徴は本発明を限定するものではない。
‐前記複数の区域における夫々の区域の輪郭は、前記空間依存性ベクトルの最大値と一致する。
‐前記空間依存性ベクトルＨは、以下の式で表すことができる。
Ｈ＝‐［ε］−１ ｄｉｖ［ε］
（‐［ε］は変形マトリックス、
‐［ε］−１は前記変形マトリックスの逆数、
‐ｄｉｖ［ε］は前記変形マトリックスとの差分を示す。）
‐前記弾性画像の演算工程は、前記血管組織の分水界分割（ｗａｔｅｒｓｈｅｄ ｓｅｇｍｅｎｔａｔｉｏｎ）のためのサブ工程を含み、該サブ工程は事前分割された画像における前記複数の区域の中から特定の区域が連続的に選択される工程である。
‐前記区域は、事前分割された画像においてＨ密度が最も高い区域から、該密度が最も低い区域へと連続的に選択される。
（参考：分水界分割（ウオータジェット処理）とは領域分割手法の一つで画像の輝度勾配を山、谷と見立てて、谷に流れる川を、山で囲まれた領域で分割する手法で、例えば団子状に２つの円がつながった部分をウオータジェット処理で２つの円に分離することが出来る。）
‐前記弾性画像の前記演算工程は、最適化のためのサブ工程を含む。該サブ工程は前記分割工程において夫々の新たなｎ＋１番目の区域が選択された後に行われ、前記最適化工程において、エラストグラムにおける実測された歪値と、前記選択されたｎ＋１番目の区域における演算された歪値との歪誤差（ｓｔｒａｉｎ ｅｒｒｏｒ）と呼ばれる誤差を最小値に抑えることにより、前記選択されたｎ＋１個の区域における弾性を評価する。
‐前記弾性画像の演算工程は、前記弾性画像の演算工程のための停止基準をそなえた前記弾性画像の比較工程を含む。
‐前記停止基準は、ｎ個の区域及びｎ＋１個の区域の弾性を特徴づける２つの連続した工程間において、好ましくは該２つの工程で演算された２つの歪誤差の差分が１０‐７より大きいことである。
‐前記停止基準は、好ましくは選択された区域の数が１０より多いことである。
【００２１】
本発明はさらに、アテロームプラークの破裂リスクを評価するための画像処理システムに関し、該システムは、上述した方法の実行手段を有する。
【００２２】
本発明はさらに、コンピューターで読み込み可能な媒体に記録されたプログラムコード指令を含むコンピュータープログラム製品として適用されるのが好ましく、該プログラム製品は、上述の方法の実行指令を含む。
【００２３】
以下、参照のための図面を示し、本発明の理解に供する。該図面は本発明を説明するためのものであって、本発明を限定するものではない。
【図面の簡単な説明】
【００２４】
【図１】図１はアテロームプラークの概要図である。
【図２】図２はアテロームプラークの破裂リスクを評価するための画像処理システムの１実施形態を示す。
【図３】図３はアテロームプラークの破裂リスクを評価するための画像処理方法の１実施形態を示す。
【発明を実施するための形態】
【００２５】
以下、本発明の方法を詳細に説明する。
【００２６】
図３に示すように、アテロームプラークの破裂リスクを評価するための画像処理方法は、
‐血管組織が血圧機能により圧迫されて生じる、血管組織の内部変形を示す「エラストグラム」と呼ばれる受像工程４０と、
‐複数の区域を有する事前に分割された領域画像を取得するために、前記エラストグラムにおいて行われる事前分割工程５０と、
‐前記エラストグラムの少なくとも１の半径方向に分割した区域における弾性を示す弾性画像を演算する工程６０であって、前記半径方向に分割した区域は前記事前分割された画像における複数区域から選択された区域に対応し、前記弾性画像を用いて使用者がアテロームプラークの破裂をリスク評価することが可能となる前記演算工程６０と
をそなえた方法であって、
該方法は以下に詳述するような特徴を有する。
‐第１に、前記血管組織が非圧縮性の弾性固体であるとみなして事前分割を行い、ラグランジュ乗数の勾配及び空間依存性ベクトルＨの勾配の合計と等しくなるように、ヤング率（縦弾性係数）の勾配を決定する。
‐第２に、前記事前分割は、前記空間依存性ベクトルにより前記事前分割した複数の区域における夫々の区域の輪郭を決定する工程を含む。
【００２７】
前記ラグランジュ乗数の勾配とは、前記媒体の非圧縮状態に基づくものであると理解される。
【００２８】
前記エラストグラムは、当業者に周知のあらゆる方法によって、例えば以下の方法により取得することが出来る。該方法とは、超音波探触子を患者の動脈の中で分析対象とする部位に挿入し、心拍によって動脈組織が圧縮／拡張されることにより、動脈を撮像した一連の超音波画像を取得するものである。該超音波画像における時系列の動態を把握することで、変形の分布を示す前記エラストグラムあるいは画像を評価することができる。
【００２９】
前記事前分割工程５０により、輪郭により画定された複数の区域（半径方向に分割した区域）を含む画像が取得可能となる。該複数の区域における夫々の区域の輪郭は、空間Ｈに依存するベクトルの最大値に対応する。
【００３０】
続いて、前記分割工程にて取得した画像に対し、前記弾性画像６０の前記演算工程が実行される。該弾性画像の演算工程は、いくつかのサブ工程を含む。
【００３１】
第１に、前記分割工程におけるサブ工程６１は、分水界線ごとに実行される。該サブ工程において、前記区域は事前分割された複数区域の中から連続的に選択される。前記区域は、事前分割された画像においてＨ密度が最も高い区域から、該密度が最も低い区域へと１つづつ順番に選択される。より詳しくは、画像における輪郭値の中で最も高い輪郭値を有する区域が最優先的に選択され、後述する複数のサブ工程は、該選択された区域に対して適用される。続いて該事前分割された画像の中から、未処理の区域の中で最も高い輪郭値を有する他の区域が選択され、その後も同様の選択が繰り返される。２つの区域が等しい輪郭値を有していた場合、どちらを優先させるかは別の基準に基づいて判断される。すなわち、例えば領域の大きさを基準として用いるか、あるいは当業者に周知の基準を用いてもよい。
【００３２】
次に、最適化のためのサブ工程６２が実行される。夫々選択された新たな区域において弾性が評価され、該新たな区域の弾性は、エラストグラムにおける実測された変形と、選択された区域における演算された変形との間に生じる変形誤差と呼ばれる誤差を最小限に抑える。
【００３３】
さらに詳しく説明するために、以下ではすでに選択済の区域をｎ個の区域とし、該ｎ個の区域における弾性がすでに評価されていることを前提とする。
【００３４】
新たなｎ＋１番目の区域は、前記サブ工程６１において選択される。当初のｎ個の区域における弾性はすでに演算されており、該ｎ個の区域の弾性演算は、選択された新たなｎ＋１番目の区域における弾性を決定することを目的として行われるものである。
【００３５】
前記選択された新たなｎ＋１番目の区域に割り当てる弾性値を、例として１とする。前記ｎ＋１個の区域における変形範囲は前記ｎ＋１番目の区域に割り当てられた弾性値（=１）を用いて演算される。
【００３６】
次に、エラストグラムに基づく実測された変形範囲と、前記選択されたｎ＋１個の区域における演算された変形範囲との間の変形誤差が算出される。この際、該選択されたｎ＋１番目の区域に割り当てられた弾性値は１とする。
【００３７】
続いて、前記ｎ＋１番目の区域に割り当てられる弾性値を２とし、前記変形範囲とそれに伴う変形誤差の演算が繰り返され、その後も同様の手順が繰り返される。
【００３８】
前記選択されたｎ＋１番目の区域に対して、前記変形誤差を最小化する弾性値が選択される。
【００３９】
当然ながら、前記選択された新たなｎ＋１番目の区域夫々の弾性値に対する評価は、実験的に弾性値を割り当てることにより実行されるわけではない。すなわち該評価は、数値解析における二分法ような相互作用的解析手法及び回帰分析等の検索プロセスによって行われる。
【００４０】
前記新たなｎ＋１番目の区域の弾性を評価した後、ｎ＋２番目の区域のさらなる選択が必要か否かを決定するテスト６３が実行される。
【００４１】
前記さらなる区域の選択要否を決定する停止基準としては、以下の２つの基準が選択可能である。
【００４２】
第１の停止基準は、選択された区域の数量関数であり、該選択された区域の数が１０を超えると、プロセスが停止し、前記事前分割された画像における他の区域はそれ以上選択されない。
【００４３】
第２の停止基準は、連続した２つの区域選択間で生じる変形誤差における変数関数である。夫々の新たなｎ＋１番目の区域の弾性を評価した後、ｎ＋１個の区域における変形誤差と、当初のｎ個の区域における変形誤差とを比較する。２つの連続したサブ工程（該サブ工程とは、すなわち区域を選択する工程、及び該選択された区域の弾性を評価する工程である）間において算出された前記２つの変形誤差の差分が１０−７以下である場合、プロセスが停止し、前記事前分割された画像における他の区域はそれ以上選択されない。
【００４４】
選択された複数の区域の弾性画像は、上記の方法により取得される。
【００４５】
前記選択された複数の区域のうちいくつかの区域は、前記基準とは別の基準に基づいて１つに統合される場合がある。該別の基準とは、領域の重複であり、例えば形状あるいは弾性値において重複が生じた場合に領域が統合される。
【００４６】
これにより、エラストグラムの中から選択された少なくとも１の関心領域における弾性画像が提供され、該夫々の範囲は前記事前分割された画像の複数の区域から選択された区域に対応する。該弾性画像は使用者に以下の情報を提供するため、使用者がアテロームプラークの破裂リスクを評価することが可能となる。
‐線維性被膜の厚さ（該厚さは例えば、使用者がカーソルを置くことにより容易に計測可能である）
‐プラークにおける選択された夫々の領域の大きさ（当業者に周知の方法により夫々の領域について算出可能である）
‐前記すべての領域に関する力学特性（本発明の方法によって該領域夫々の弾性が評価可能である）
【００４７】
図３に示すように、前記方法は、アテロームプラークの破裂リスクを評価するための画像処理システムにおいて実行可能である。該システムはデータ入力手段１と、処理手段２と、表示手段３とを含む。前記データ入力手段１により、使用者はデータを入力することが出来る。該データ入力手段１は、例えばコンピュータキーボード及び／あるいはマウスを備える。前記処理手段２は放射線画像の処理に用いられる。前記表示手段３により、特に前記処理手段によって処理された画像が表示可能である。
【００４８】
本発明の方法及びシステムに関する理論的側面における詳細を以下記載する。
【００４９】
（一般的原理）
本理論的研究は、入力データを事前に分割する工程と、分水界分割とアテロームプラークの構成成分の形態及びヤング率を抽出するための最適化手順とを組み合わせた新しいアプローチとによって、アテロームプラークの複雑なモジュログラムを決定するためのものである。
【００５０】
連続体力学の理論に基づく前記アプローチは、血管内超音波法（ＩＶＵＳ）により生体内撮像された患者の冠動脈病変に対して用いられ、好ましい結果を生じている。
【００５１】
アテロームプラークの脆弱性評価に影響する可能性のある様々な因子を検証することにより、本発明の方法におけるロバスト性及び性能が確認されている。
【００５２】
図３において、パラメトリック有限要素モデル（ＰＦＥＭ：ｐａｒａｍｅｔｒｉｃ ｆｉｎｉｔｅ ｅｌｅｍｅｎｔ ｍｏｄｅｌ）に基づく、本発明の方法が含む複数の異なる工程を示す。
【００５３】
第１の工程において、Ｓｕｍｉ他によって提案された理論モデルを拡張することにより導き出された弾性勾配のアプローチ（生命分子化学分野において細胞の弾性勾配の解析手法は知られている。）を用いて、アテロームプラークの不均質性を証明する基準が得られる。
【００５４】
さらに詳しくは、前記基準を用いると、半径方向の変形領域により、アテロームプラークの不均質性が証明される。この数学的プレコンディショニングデータはＳｕｍｉ修正変換（ＭＳＴ：Ｍｏｄｉｆｉｅｄ Ｓｕｍｉ’ｓ Ｔｒａｎｓｆｏｒｍ）と呼ばれる。
【００５５】
前記数学的プレコンディショニングデータにより、ヤング率の空間導関数を決定することが出来る。
【００５６】
第２の工程において、動的分水界分割（ＤＷＳ：Ｄｙｎａｍｉｃ ＷａｔｅｒＳｈｅｄ Ｓｅｇｍｅｎｔａｔｉｏｎ）は、弾性画像を取得するためのプレコンディショニング工程において取得したデータに適用される。
【００５７】
（連続体力学の理論に基づく輪郭検出工程）
（Ｓｋｏｒｏｖｏｄａ変換（Ｓｋｏｒｏｖｏｄａ ｔｒａｎｓｆｏｒｍ））
Ｓｋｏｒｏｖｏｄａ他は、ヤング率が夫々Ｅｉｎｔ及びＥｅｘｔである前記含有物と周囲媒体との間の境界面における局所的制約付きの連続方程式を用いて、含有物の輪郭をプラーク内変形領域から抽出するための数学的基準を提案している。
【００５８】
前記各輪郭点は、以下の割合を評価することによって検出される。
任意方向に沿った、且つ微小領域内における、
‐内部半径方向の変形εｉｊｉｎｔにおける構成成分の、
‐外部半径方向の変形εｉｊｅｘｔにおける構成成分に対する割合。
【００５９】
したがって、前記最大比の位置により、含有物の輪郭を確認することが可能となる。（血管の構成）素材が非圧縮性であり、応力テンソルεｒｒ−の動径成分のみを考慮に入れることを前提として得られた前記基準を示す、簡略化された数学的表現は以下のようになる。
【数１】
【００６０】
前記基準は本発明の方法によってもたらされる改善の度合いを評価するために用いられる。
【００６１】
（Ｓｕｍｉ修正変換（Ｓｕｍｉ’ｓ ｍｏｄｉｆｉｅｄ ｔｒａｎｓｆｏｒｍ））
空間内で力学特性が絶えず変化すると仮定し、さらに平面応力負荷下における等方向性の非圧縮性媒体に対して局所均衡方程式を立てることにより、Ｓｕｍｉはヤング率の勾配ベクトルと応力テンソル成分との間に関係性を見出した。この基準は、連続体（ｃｏｎｔｉｎｕｕｍ ｍｅｄｉｕｍ）内に埋没した円形状不均質成分のための剛性画像を再構成するために用いられ、好ましい結果が得られている。
【００６２】
しかしながら、平面応力条件は血管変形のモデリングとは関連を持たないことは明らかである。
【００６３】
したがって、このアプローチを、平面変形条件下において負荷される非圧縮性弾性媒体にまで拡張することが必要となる。アテロームプラークが非圧縮性であることを考慮に入れると、ヤング率の勾配表現は以下の式によって推定される。
線形弾性を仮定すると、不均一性媒体はフックの法則に基づき以下のように記述される。
【数２】
（‐［σ］及び［ε］は応力及び変形テンソル、
‐［Ｉ］は恒等行列、
‐Ｅは位置ベクトルｘの任意生産であるヤング率、
‐ｐは以下の運動学的応力に示される素材の非圧縮性により生じるラグランジュ乗数を示す。）
【数３】
【００６４】
重力及び慣性力を無視することにより、応力テンソル［σ］は以下に示す局所均衡方程式を満たす。
【数４】
【００６５】
前記数式２及び数式３を数式４に代入すると、ヤング率の弾性勾配に関する以下の方程式が得られる。
【数５】
【００６６】
ラグランジュ勾配を実験的に計測することは不可能なため、前記事前分割工程を考慮すると、上記方程式の第２項がより着目される。該第２項は含有物の輪郭位置を決定するための基準として用いられる。
【００６７】
したがって、空間Ｈに依存するベクトルは上記方程式の第２項に対応するように導入され、下記のように定義される。
【数６】
【００６８】
上記方程式のＨは空間依存性ベクトルを示す。含有物の輪郭画像を再構成することに成功すると、本発明のアプローチが帰納的に妥当化される。
【００６９】
極座標の２次元（２Ｄ）システムにおいて、空間Ｈに依存するベクトルの構成要素は以下の数式に示される。
【数７】
【００７０】
上記の数式において、εηθ及びεｒθは夫々応力テンソルの接線成分及びせん断成分である。
【００７１】
局所ヤング率の相対的変動に依存する基準を定義するためには、ベクトルＨにおける射影ｄｗは、原点位置ベクトルｄｘに導入され、以下の数式に示される。
【数８】
【００７２】
前記数式５、６及び８を用いると、ｄｗは有利には以下のように表現される。
【数９】
【００７３】
上記数式により、ｄｗは局所ヤング率の相対的変動に依存していることが分かる。ラグランジュ乗数の勾配を構成する成分が小さいほど、ｄＥ／Ｅによるｄｗの近似値算出は正確になる。
【００７４】
実際に現在の超音波技術を用いて応力テンソルを推定する際には、半径方向応力が最も簡単な指標となる。
【００７５】
それゆえ上記基準は、圧縮応力を無視して単純化される。したがって、基準ｄｗは以下の数式表現により与えられる。
【数１０】
【００７６】
半径方向応力構成成分の変動関数は、ソベルマスクを用いて応力場を畳み込むことにより算出される。したがって、ｄｗの振幅は、前記（数１０において）定義された規則的ポーラーメッシュの関数として算出され、それにより前記ＭＳＴ画像が取得可能となる。夫々の領域の区域は、ｄｗの振幅が最大値を示す線に対応している。この数学的プレコンディショニング工程は、Ｓｕｍｉ修正変換（ＭＳＴ：Ｍｏｄｉｆｉｅｄ Ｓｕｍｉ’ｓ Ｔｒａｎｓｆｏｒｍ）と呼ばれ、アテロームプラークの不均一性を証明するために用いられてきた。
【００７７】
（ヤング率弾性勾配の再構成）
図３は、本発明の方法が含む連続した諸工程を図式化したものである。該工程において、パラメトリック有限要素モデルを用いて半径方向応力場の弾性画像が特定される。
【００７８】
前記プレコンディショニング工程で得られたＭＳＴ領域に始まり、動的分水界分割（ＤＷＳ：Ｄｙｎａｍｉｃ ＷａｔｅｒＳｈｅｄ Ｓｅｇｍｅｎｔａｔｉｏｎ）を最適化手順と組み合わせて適用することにより、分析対象となるプラークの弾性マップ（モジュログラム）が得られる。
【００７９】
（例えば、プラークにおける一定数の含有物のために）分水界分割工程が反復されるが、該夫々の工程において、夫々の前記含有物は均一の剛性を有すると仮定される。
【００８０】
前記同一のヤング率弾性勾配は、勾配の最適化手順を用いて特定される。該最適化手順は、コンピュータ算出されたεｒｒｍｅａｎｓと、プラークの半径方向応力におけるεｒｒｃｏｍｐとの二乗平均平方根（ＲＭＳ：Ｒｏｏｔ ｍｅａｎ ｓｑｕａｒｅ）誤差を最小限に抑えるためのものであり、すなわち以下の数式で示される。
【数１１】
【００８１】
上記数式において、Ｎはプラークの網目構造における結節の総数を示し、Ｎｉは前記ポーラーメッシュ（構造を有する）プラークにおける結節の序数を示す。
【００８２】
上記最適化手順においては応力が負荷される。すなわち、未知のヤング率弾性勾配は０．１ｋｐａから１０４ｋｐａの間に分類される。該分類工程においては、許容手順において１０−８以下の閾値が許容された場合、あるいは領域の総数が３０である場合に、前記一連のヤング率に対する解は妥当であるとみなされる。
【００８３】
アテロームプラーク内の含有物数を最大１０個まで増加させることにより、より細かな領域分割が可能となる。
【００８４】
アテロームプラークの形態学における研究から、前記上限が現実的な値であることが知られている。
【００８５】
さらに、アルゴリズム的演算を考慮すると、前記上限が前記最適化工程によっては到達し得ない値であることは、帰納的に証明済である。
【００８６】
以上で述べた方法及びシステムに関し、本明細書の教示から実質的に逸脱することなく数多くの変更を行なうことができる。しかしながら、該変更は本願における特許請求の範囲に含まれるものである。
【技術分野】
【０００１】
本発明は、技術分野としては医療用画像処理に関し、さらに詳しくは血管組織のアテロームプラーク（血管内の粥状硬化部（脂質コア））破裂のリスク評価のための画像処理方法に関する。
【背景技術】
【０００２】
年々、心臓血管疾患による死者は癌による死者よりも多く発生している。狭心症、心筋梗塞、突然死及び脳卒中はいずれもアテローム性動脈硬化の影響によるものである。
【０００３】
組織学的研究によると、脆弱なアテロームプラークは壊死した大きな細胞外体１０として定義される。該細胞外体１０には、マクロファージ浸潤を伴う薄膜（例として膜厚１００μｍ以下）の線維性被膜２０が付随形成される。
【０００４】
心筋を灌漑する冠動脈３０は、アテロームプラークの突発的破裂により生じた血栓に狭窄され、したがって前記心筋における酸素が供給されない領域は損傷を受ける。この作用は突然死、不安定狭心症あるいは心筋梗塞の原因となる。
【０００５】
上記のような脆弱なプラークは、血管内超音波法（ＩＶＵＳ：Ｉｎｔｒａｖａｓｃｕｌａｒ Ｕｌｔｒａｓｏｕｎｄ）、光コヒーレンス・トモグラフィー（ＯＣＴ：Ｏｐｔｉｃａｌ Ｃｏｈｅｒｅｎｃｅ Ｔｏｍｏｇｒａｐｈｙ）、ＣＴスキャン及び核磁気共鳴画像法（ＭＲＩ：Ｍａｇｎｅｔｉｃ Ｒｅｓｏｎａｎｃｅ Ｉｍａｇｉｎｇ）を含む多様な技術により、臨床的に検出可能であることが、数件の報告によって明らかにされている。
【０００６】
しかしながら、形態及びプラーク組織の解析に基づく破裂予測は、極めて不正確なばかりでなく、十分な数のリスク指標が得られないという点においても好ましくない。
【０００７】
冠動脈プラークの破裂予測に関する目下の課題解決には、
‐線維層厚の正確な定量化及び
‐壊死体の形態解析
が必要となるだけでなく、特に前記線維層厚及び壊死体の力学特性に関する正確な知識、さらには動脈壁及びプラーク構成成分の弾性に関する正確な知識が求められる。
【０００８】
上記データ項目はプラーク破裂に関する好ましい生体力学的リスク指標となり得る基準であって、該基準により線維層に付加されるピーク応力を正確かつ効率的に評価することが可能となる。
【０００９】
しかしながら、生体内におけるプラークは不均質であり、さらにプラークは成長の過程で力学特性の変化を伴うため、プラークの力学的特性をｉｎ ｖｉｖｏ（生体内）で認識することは困難である。
【００１０】
さらに詳しくは、弾性（あるいは剛性）の画像化、すなわち血管組織の弾性を解析して画像化することは、応力分布を確実に算出するための必要前提条件とみなされる。
【００１１】
前記弾性画像の演算は、これまで様々なアプローチにより試みられてきた。
【００１２】
１心周期におけるプラーク負荷を入力データとして用い、前記血管弾性の画像を判断する研究が散見されている。この判断方法に関し、直接のアプローチあるいは相互作用手順が提案されている。
【００１３】
前記相互作用手順は、実測された負荷と、数値モデルを用いて演算された負荷との間の誤差を最小限に抑えることを目的とした最適化アルゴリズムを用いるものである。この文脈において、プラーク弾性の再現性向上は、前記最適化アルゴリズムの性能に依存している。
【００１４】
例えば、プラーク構成成分における弾性モジュライを抽出するための最適化アルゴリズム及び領域分割アルゴリズムが、いくつかの研究団体によって開発されている。しかしながらこれらの方法は、アテロームプラークが互いに隣接する複数の細胞の壊死体を含むという頻発事例において、弾性画像を取得するに十分な有効性あるいは性能レベルに達していない。
【００１５】
さらに加えて、これらの方法は以下に述べる部分においては十分に機能しない。第１に、該方法ではカルシウム含有物を抽出することが出来ないため、プラーク構成成分、特に剛体要素の特定が困難である。第２に、線維性被膜の厚さを正確に定量化することが困難であるため、プラークの脆弱性を正確に分析することに支障が生じる。
【００１６】
したがって、アテロームプラークの脆弱性の度合いを正確に分析するために、施術者は以下の項目を正確に把握している必要がある。
‐線維性被膜の厚さ
‐プラークが含有する様々な構成成分の大きさ
‐上記すべての成分における力学特性
【発明の概要】
【発明が解決しようとする課題】
【００１７】
本発明の目的は、アテロームプラークの形態上の詳細だけでなく、弾性勾配の正確な弾性画像、特にプラークが複数の互いに隣接する細胞の壊死体を有し、該壊死体にカルシウム含有物が含まれる可能性がある場合に、複数の構成成分における弾性を示す弾性画像を提供することにより、施術者にアテロームプラークの脆弱性の度合いを正確に評価させることである。
【００１８】
本発明の方法を用いれば、例えば、線維性被膜の厚さ及び壊死体の大きさ（厚さ）のような、アテロームプラークの破裂に関するリスク評価のための他の正確な基準を施術者に提供することも可能である。
【課題を解決するための手段】
【００１９】
上述の目的において、アテロームプラーク破裂のリスク評価のための画像処理方法を提案する。該方法は、
‐血管組織が血圧機能により圧迫されることにより生じる内部変形を示す「エラストグラム（ｅｌａｓｔｏｇｒａｍ）」と呼ばれる受像工程と、
‐事前に領域分割された（ｐｒｅ−ｓｅｇｍｅｎｔｅｄ）複数の区域を有する画像を取得するために、前記エラストグラムを事前に領域分割する（ｐｒｅ−ｓｅｇｍｅｎｔａｔｉｏｎ）工程と、
‐前記エラストグラムの少なくとも１の領域における弾性を示す弾性画像を演算する工程であって、前記領域は前記事前分割された画像における複数区域から選択された区域に対応し、前記弾性画像を用いて使用者はアテロームプラークの破裂をリスク評価することが可能となる前記弾性画像の演算工程と
をそなえた方法であって、
‐第１に、前記血管組織が非圧縮性の弾性固体であるとみなして事前領域分割を行い、ラグランジュ乗数の勾配及び空間依存性ベクトルの勾配の合計と等しくなるように、前記血管組織のヤング率（弾性係数）の勾配を決定し、
‐第２に、前記事前領域分割は、前記複数の区域（半径方向に分割した区域）における夫々の区域の輪郭を前記空間依存性ベクトルにより決定する工程を含む。
【００２０】
本発明の方法の好ましい特徴を以下に述べるが、該特徴は本発明を限定するものではない。
‐前記複数の区域における夫々の区域の輪郭は、前記空間依存性ベクトルの最大値と一致する。
‐前記空間依存性ベクトルＨは、以下の式で表すことができる。
Ｈ＝‐［ε］−１ ｄｉｖ［ε］
（‐［ε］は変形マトリックス、
‐［ε］−１は前記変形マトリックスの逆数、
‐ｄｉｖ［ε］は前記変形マトリックスとの差分を示す。）
‐前記弾性画像の演算工程は、前記血管組織の分水界分割（ｗａｔｅｒｓｈｅｄ ｓｅｇｍｅｎｔａｔｉｏｎ）のためのサブ工程を含み、該サブ工程は事前分割された画像における前記複数の区域の中から特定の区域が連続的に選択される工程である。
‐前記区域は、事前分割された画像においてＨ密度が最も高い区域から、該密度が最も低い区域へと連続的に選択される。
（参考：分水界分割（ウオータジェット処理）とは領域分割手法の一つで画像の輝度勾配を山、谷と見立てて、谷に流れる川を、山で囲まれた領域で分割する手法で、例えば団子状に２つの円がつながった部分をウオータジェット処理で２つの円に分離することが出来る。）
‐前記弾性画像の前記演算工程は、最適化のためのサブ工程を含む。該サブ工程は前記分割工程において夫々の新たなｎ＋１番目の区域が選択された後に行われ、前記最適化工程において、エラストグラムにおける実測された歪値と、前記選択されたｎ＋１番目の区域における演算された歪値との歪誤差（ｓｔｒａｉｎ ｅｒｒｏｒ）と呼ばれる誤差を最小値に抑えることにより、前記選択されたｎ＋１個の区域における弾性を評価する。
‐前記弾性画像の演算工程は、前記弾性画像の演算工程のための停止基準をそなえた前記弾性画像の比較工程を含む。
‐前記停止基準は、ｎ個の区域及びｎ＋１個の区域の弾性を特徴づける２つの連続した工程間において、好ましくは該２つの工程で演算された２つの歪誤差の差分が１０‐７より大きいことである。
‐前記停止基準は、好ましくは選択された区域の数が１０より多いことである。
【００２１】
本発明はさらに、アテロームプラークの破裂リスクを評価するための画像処理システムに関し、該システムは、上述した方法の実行手段を有する。
【００２２】
本発明はさらに、コンピューターで読み込み可能な媒体に記録されたプログラムコード指令を含むコンピュータープログラム製品として適用されるのが好ましく、該プログラム製品は、上述の方法の実行指令を含む。
【００２３】
以下、参照のための図面を示し、本発明の理解に供する。該図面は本発明を説明するためのものであって、本発明を限定するものではない。
【図面の簡単な説明】
【００２４】
【図１】図１はアテロームプラークの概要図である。
【図２】図２はアテロームプラークの破裂リスクを評価するための画像処理システムの１実施形態を示す。
【図３】図３はアテロームプラークの破裂リスクを評価するための画像処理方法の１実施形態を示す。
【発明を実施するための形態】
【００２５】
以下、本発明の方法を詳細に説明する。
【００２６】
図３に示すように、アテロームプラークの破裂リスクを評価するための画像処理方法は、
‐血管組織が血圧機能により圧迫されて生じる、血管組織の内部変形を示す「エラストグラム」と呼ばれる受像工程４０と、
‐複数の区域を有する事前に分割された領域画像を取得するために、前記エラストグラムにおいて行われる事前分割工程５０と、
‐前記エラストグラムの少なくとも１の半径方向に分割した区域における弾性を示す弾性画像を演算する工程６０であって、前記半径方向に分割した区域は前記事前分割された画像における複数区域から選択された区域に対応し、前記弾性画像を用いて使用者がアテロームプラークの破裂をリスク評価することが可能となる前記演算工程６０と
をそなえた方法であって、
該方法は以下に詳述するような特徴を有する。
‐第１に、前記血管組織が非圧縮性の弾性固体であるとみなして事前分割を行い、ラグランジュ乗数の勾配及び空間依存性ベクトルＨの勾配の合計と等しくなるように、ヤング率（縦弾性係数）の勾配を決定する。
‐第２に、前記事前分割は、前記空間依存性ベクトルにより前記事前分割した複数の区域における夫々の区域の輪郭を決定する工程を含む。
【００２７】
前記ラグランジュ乗数の勾配とは、前記媒体の非圧縮状態に基づくものであると理解される。
【００２８】
前記エラストグラムは、当業者に周知のあらゆる方法によって、例えば以下の方法により取得することが出来る。該方法とは、超音波探触子を患者の動脈の中で分析対象とする部位に挿入し、心拍によって動脈組織が圧縮／拡張されることにより、動脈を撮像した一連の超音波画像を取得するものである。該超音波画像における時系列の動態を把握することで、変形の分布を示す前記エラストグラムあるいは画像を評価することができる。
【００２９】
前記事前分割工程５０により、輪郭により画定された複数の区域（半径方向に分割した区域）を含む画像が取得可能となる。該複数の区域における夫々の区域の輪郭は、空間Ｈに依存するベクトルの最大値に対応する。
【００３０】
続いて、前記分割工程にて取得した画像に対し、前記弾性画像６０の前記演算工程が実行される。該弾性画像の演算工程は、いくつかのサブ工程を含む。
【００３１】
第１に、前記分割工程におけるサブ工程６１は、分水界線ごとに実行される。該サブ工程において、前記区域は事前分割された複数区域の中から連続的に選択される。前記区域は、事前分割された画像においてＨ密度が最も高い区域から、該密度が最も低い区域へと１つづつ順番に選択される。より詳しくは、画像における輪郭値の中で最も高い輪郭値を有する区域が最優先的に選択され、後述する複数のサブ工程は、該選択された区域に対して適用される。続いて該事前分割された画像の中から、未処理の区域の中で最も高い輪郭値を有する他の区域が選択され、その後も同様の選択が繰り返される。２つの区域が等しい輪郭値を有していた場合、どちらを優先させるかは別の基準に基づいて判断される。すなわち、例えば領域の大きさを基準として用いるか、あるいは当業者に周知の基準を用いてもよい。
【００３２】
次に、最適化のためのサブ工程６２が実行される。夫々選択された新たな区域において弾性が評価され、該新たな区域の弾性は、エラストグラムにおける実測された変形と、選択された区域における演算された変形との間に生じる変形誤差と呼ばれる誤差を最小限に抑える。
【００３３】
さらに詳しく説明するために、以下ではすでに選択済の区域をｎ個の区域とし、該ｎ個の区域における弾性がすでに評価されていることを前提とする。
【００３４】
新たなｎ＋１番目の区域は、前記サブ工程６１において選択される。当初のｎ個の区域における弾性はすでに演算されており、該ｎ個の区域の弾性演算は、選択された新たなｎ＋１番目の区域における弾性を決定することを目的として行われるものである。
【００３５】
前記選択された新たなｎ＋１番目の区域に割り当てる弾性値を、例として１とする。前記ｎ＋１個の区域における変形範囲は前記ｎ＋１番目の区域に割り当てられた弾性値（=１）を用いて演算される。
【００３６】
次に、エラストグラムに基づく実測された変形範囲と、前記選択されたｎ＋１個の区域における演算された変形範囲との間の変形誤差が算出される。この際、該選択されたｎ＋１番目の区域に割り当てられた弾性値は１とする。
【００３７】
続いて、前記ｎ＋１番目の区域に割り当てられる弾性値を２とし、前記変形範囲とそれに伴う変形誤差の演算が繰り返され、その後も同様の手順が繰り返される。
【００３８】
前記選択されたｎ＋１番目の区域に対して、前記変形誤差を最小化する弾性値が選択される。
【００３９】
当然ながら、前記選択された新たなｎ＋１番目の区域夫々の弾性値に対する評価は、実験的に弾性値を割り当てることにより実行されるわけではない。すなわち該評価は、数値解析における二分法ような相互作用的解析手法及び回帰分析等の検索プロセスによって行われる。
【００４０】
前記新たなｎ＋１番目の区域の弾性を評価した後、ｎ＋２番目の区域のさらなる選択が必要か否かを決定するテスト６３が実行される。
【００４１】
前記さらなる区域の選択要否を決定する停止基準としては、以下の２つの基準が選択可能である。
【００４２】
第１の停止基準は、選択された区域の数量関数であり、該選択された区域の数が１０を超えると、プロセスが停止し、前記事前分割された画像における他の区域はそれ以上選択されない。
【００４３】
第２の停止基準は、連続した２つの区域選択間で生じる変形誤差における変数関数である。夫々の新たなｎ＋１番目の区域の弾性を評価した後、ｎ＋１個の区域における変形誤差と、当初のｎ個の区域における変形誤差とを比較する。２つの連続したサブ工程（該サブ工程とは、すなわち区域を選択する工程、及び該選択された区域の弾性を評価する工程である）間において算出された前記２つの変形誤差の差分が１０−７以下である場合、プロセスが停止し、前記事前分割された画像における他の区域はそれ以上選択されない。
【００４４】
選択された複数の区域の弾性画像は、上記の方法により取得される。
【００４５】
前記選択された複数の区域のうちいくつかの区域は、前記基準とは別の基準に基づいて１つに統合される場合がある。該別の基準とは、領域の重複であり、例えば形状あるいは弾性値において重複が生じた場合に領域が統合される。
【００４６】
これにより、エラストグラムの中から選択された少なくとも１の関心領域における弾性画像が提供され、該夫々の範囲は前記事前分割された画像の複数の区域から選択された区域に対応する。該弾性画像は使用者に以下の情報を提供するため、使用者がアテロームプラークの破裂リスクを評価することが可能となる。
‐線維性被膜の厚さ（該厚さは例えば、使用者がカーソルを置くことにより容易に計測可能である）
‐プラークにおける選択された夫々の領域の大きさ（当業者に周知の方法により夫々の領域について算出可能である）
‐前記すべての領域に関する力学特性（本発明の方法によって該領域夫々の弾性が評価可能である）
【００４７】
図３に示すように、前記方法は、アテロームプラークの破裂リスクを評価するための画像処理システムにおいて実行可能である。該システムはデータ入力手段１と、処理手段２と、表示手段３とを含む。前記データ入力手段１により、使用者はデータを入力することが出来る。該データ入力手段１は、例えばコンピュータキーボード及び／あるいはマウスを備える。前記処理手段２は放射線画像の処理に用いられる。前記表示手段３により、特に前記処理手段によって処理された画像が表示可能である。
【００４８】
本発明の方法及びシステムに関する理論的側面における詳細を以下記載する。
【００４９】
（一般的原理）
本理論的研究は、入力データを事前に分割する工程と、分水界分割とアテロームプラークの構成成分の形態及びヤング率を抽出するための最適化手順とを組み合わせた新しいアプローチとによって、アテロームプラークの複雑なモジュログラムを決定するためのものである。
【００５０】
連続体力学の理論に基づく前記アプローチは、血管内超音波法（ＩＶＵＳ）により生体内撮像された患者の冠動脈病変に対して用いられ、好ましい結果を生じている。
【００５１】
アテロームプラークの脆弱性評価に影響する可能性のある様々な因子を検証することにより、本発明の方法におけるロバスト性及び性能が確認されている。
【００５２】
図３において、パラメトリック有限要素モデル（ＰＦＥＭ：ｐａｒａｍｅｔｒｉｃ ｆｉｎｉｔｅ ｅｌｅｍｅｎｔ ｍｏｄｅｌ）に基づく、本発明の方法が含む複数の異なる工程を示す。
【００５３】
第１の工程において、Ｓｕｍｉ他によって提案された理論モデルを拡張することにより導き出された弾性勾配のアプローチ（生命分子化学分野において細胞の弾性勾配の解析手法は知られている。）を用いて、アテロームプラークの不均質性を証明する基準が得られる。
【００５４】
さらに詳しくは、前記基準を用いると、半径方向の変形領域により、アテロームプラークの不均質性が証明される。この数学的プレコンディショニングデータはＳｕｍｉ修正変換（ＭＳＴ：Ｍｏｄｉｆｉｅｄ Ｓｕｍｉ’ｓ Ｔｒａｎｓｆｏｒｍ）と呼ばれる。
【００５５】
前記数学的プレコンディショニングデータにより、ヤング率の空間導関数を決定することが出来る。
【００５６】
第２の工程において、動的分水界分割（ＤＷＳ：Ｄｙｎａｍｉｃ ＷａｔｅｒＳｈｅｄ Ｓｅｇｍｅｎｔａｔｉｏｎ）は、弾性画像を取得するためのプレコンディショニング工程において取得したデータに適用される。
【００５７】
（連続体力学の理論に基づく輪郭検出工程）
（Ｓｋｏｒｏｖｏｄａ変換（Ｓｋｏｒｏｖｏｄａ ｔｒａｎｓｆｏｒｍ））
Ｓｋｏｒｏｖｏｄａ他は、ヤング率が夫々Ｅｉｎｔ及びＥｅｘｔである前記含有物と周囲媒体との間の境界面における局所的制約付きの連続方程式を用いて、含有物の輪郭をプラーク内変形領域から抽出するための数学的基準を提案している。
【００５８】
前記各輪郭点は、以下の割合を評価することによって検出される。
任意方向に沿った、且つ微小領域内における、
‐内部半径方向の変形εｉｊｉｎｔにおける構成成分の、
‐外部半径方向の変形εｉｊｅｘｔにおける構成成分に対する割合。
【００５９】
したがって、前記最大比の位置により、含有物の輪郭を確認することが可能となる。（血管の構成）素材が非圧縮性であり、応力テンソルεｒｒ−の動径成分のみを考慮に入れることを前提として得られた前記基準を示す、簡略化された数学的表現は以下のようになる。
【数１】
【００６０】
前記基準は本発明の方法によってもたらされる改善の度合いを評価するために用いられる。
【００６１】
（Ｓｕｍｉ修正変換（Ｓｕｍｉ’ｓ ｍｏｄｉｆｉｅｄ ｔｒａｎｓｆｏｒｍ））
空間内で力学特性が絶えず変化すると仮定し、さらに平面応力負荷下における等方向性の非圧縮性媒体に対して局所均衡方程式を立てることにより、Ｓｕｍｉはヤング率の勾配ベクトルと応力テンソル成分との間に関係性を見出した。この基準は、連続体（ｃｏｎｔｉｎｕｕｍ ｍｅｄｉｕｍ）内に埋没した円形状不均質成分のための剛性画像を再構成するために用いられ、好ましい結果が得られている。
【００６２】
しかしながら、平面応力条件は血管変形のモデリングとは関連を持たないことは明らかである。
【００６３】
したがって、このアプローチを、平面変形条件下において負荷される非圧縮性弾性媒体にまで拡張することが必要となる。アテロームプラークが非圧縮性であることを考慮に入れると、ヤング率の勾配表現は以下の式によって推定される。
線形弾性を仮定すると、不均一性媒体はフックの法則に基づき以下のように記述される。
【数２】
（‐［σ］及び［ε］は応力及び変形テンソル、
‐［Ｉ］は恒等行列、
‐Ｅは位置ベクトルｘの任意生産であるヤング率、
‐ｐは以下の運動学的応力に示される素材の非圧縮性により生じるラグランジュ乗数を示す。）
【数３】
【００６４】
重力及び慣性力を無視することにより、応力テンソル［σ］は以下に示す局所均衡方程式を満たす。
【数４】
【００６５】
前記数式２及び数式３を数式４に代入すると、ヤング率の弾性勾配に関する以下の方程式が得られる。
【数５】
【００６６】
ラグランジュ勾配を実験的に計測することは不可能なため、前記事前分割工程を考慮すると、上記方程式の第２項がより着目される。該第２項は含有物の輪郭位置を決定するための基準として用いられる。
【００６７】
したがって、空間Ｈに依存するベクトルは上記方程式の第２項に対応するように導入され、下記のように定義される。
【数６】
【００６８】
上記方程式のＨは空間依存性ベクトルを示す。含有物の輪郭画像を再構成することに成功すると、本発明のアプローチが帰納的に妥当化される。
【００６９】
極座標の２次元（２Ｄ）システムにおいて、空間Ｈに依存するベクトルの構成要素は以下の数式に示される。
【数７】
【００７０】
上記の数式において、εηθ及びεｒθは夫々応力テンソルの接線成分及びせん断成分である。
【００７１】
局所ヤング率の相対的変動に依存する基準を定義するためには、ベクトルＨにおける射影ｄｗは、原点位置ベクトルｄｘに導入され、以下の数式に示される。
【数８】
【００７２】
前記数式５、６及び８を用いると、ｄｗは有利には以下のように表現される。
【数９】
【００７３】
上記数式により、ｄｗは局所ヤング率の相対的変動に依存していることが分かる。ラグランジュ乗数の勾配を構成する成分が小さいほど、ｄＥ／Ｅによるｄｗの近似値算出は正確になる。
【００７４】
実際に現在の超音波技術を用いて応力テンソルを推定する際には、半径方向応力が最も簡単な指標となる。
【００７５】
それゆえ上記基準は、圧縮応力を無視して単純化される。したがって、基準ｄｗは以下の数式表現により与えられる。
【数１０】
【００７６】
半径方向応力構成成分の変動関数は、ソベルマスクを用いて応力場を畳み込むことにより算出される。したがって、ｄｗの振幅は、前記（数１０において）定義された規則的ポーラーメッシュの関数として算出され、それにより前記ＭＳＴ画像が取得可能となる。夫々の領域の区域は、ｄｗの振幅が最大値を示す線に対応している。この数学的プレコンディショニング工程は、Ｓｕｍｉ修正変換（ＭＳＴ：Ｍｏｄｉｆｉｅｄ Ｓｕｍｉ’ｓ Ｔｒａｎｓｆｏｒｍ）と呼ばれ、アテロームプラークの不均一性を証明するために用いられてきた。
【００７７】
（ヤング率弾性勾配の再構成）
図３は、本発明の方法が含む連続した諸工程を図式化したものである。該工程において、パラメトリック有限要素モデルを用いて半径方向応力場の弾性画像が特定される。
【００７８】
前記プレコンディショニング工程で得られたＭＳＴ領域に始まり、動的分水界分割（ＤＷＳ：Ｄｙｎａｍｉｃ ＷａｔｅｒＳｈｅｄ Ｓｅｇｍｅｎｔａｔｉｏｎ）を最適化手順と組み合わせて適用することにより、分析対象となるプラークの弾性マップ（モジュログラム）が得られる。
【００７９】
（例えば、プラークにおける一定数の含有物のために）分水界分割工程が反復されるが、該夫々の工程において、夫々の前記含有物は均一の剛性を有すると仮定される。
【００８０】
前記同一のヤング率弾性勾配は、勾配の最適化手順を用いて特定される。該最適化手順は、コンピュータ算出されたεｒｒｍｅａｎｓと、プラークの半径方向応力におけるεｒｒｃｏｍｐとの二乗平均平方根（ＲＭＳ：Ｒｏｏｔ ｍｅａｎ ｓｑｕａｒｅ）誤差を最小限に抑えるためのものであり、すなわち以下の数式で示される。
【数１１】
【００８１】
上記数式において、Ｎはプラークの網目構造における結節の総数を示し、Ｎｉは前記ポーラーメッシュ（構造を有する）プラークにおける結節の序数を示す。
【００８２】
上記最適化手順においては応力が負荷される。すなわち、未知のヤング率弾性勾配は０．１ｋｐａから１０４ｋｐａの間に分類される。該分類工程においては、許容手順において１０−８以下の閾値が許容された場合、あるいは領域の総数が３０である場合に、前記一連のヤング率に対する解は妥当であるとみなされる。
【００８３】
アテロームプラーク内の含有物数を最大１０個まで増加させることにより、より細かな領域分割が可能となる。
【００８４】
アテロームプラークの形態学における研究から、前記上限が現実的な値であることが知られている。
【００８５】
さらに、アルゴリズム的演算を考慮すると、前記上限が前記最適化工程によっては到達し得ない値であることは、帰納的に証明済である。
【００８６】
以上で述べた方法及びシステムに関し、本明細書の教示から実質的に逸脱することなく数多くの変更を行なうことができる。しかしながら、該変更は本願における特許請求の範囲に含まれるものである。
【特許請求の範囲】
【請求項１】
アテロームプラーク破裂のリスク評価のための画像処理方法であって、
該方法は、
‐血管組織が血圧機能により圧迫されることにより生じる内部変形を示す「エラストグラム」と呼ばれる受像工程（４０）と、
‐複数の区域を有する、事前に分割された画像を取得するために、前記エラストグラムを事前に分割する工程（５０）と、
‐前記エラストグラムの少なくとも１の領域における弾性を示す弾性画像を演算する工程（６０）であって、前記領域は前記事前分割された画像における複数区域から選択された区域に対応し、前記弾性画像を用いて使用者はアテロームプラークの破裂をリスク評価することが可能となる前記弾性画像の演算工程（６０）と
をそなえた方法であって、
‐第１に、前記血管組織が非圧縮性の弾性固体であるとみなして事前分割を行い、ラグランジュ乗数の勾配及び空間依存性ベクトルの勾配の合計と等しくなるように、ヤング率の弾性勾配を決定し、
‐第２に、前記事前分割は、前記複数の区域における夫々の区域の輪郭を前記空間依存性ベクトルにより決定する工程を含むアテロームプラーク破裂のリスク評価のための画像処理方法。
【請求項２】
前記複数の区域における夫々の区域の輪郭は、前記空間依存性ベクトルの最大値と一致することを特徴とする請求項１に記載の方法。
【請求項３】
前記空間依存性ベクトルＨが、以下の式で示されることを特徴とする請求項１または２に記載の方法。
Ｈ＝‐［ε］−１ ｄｉｖ［ε］
（‐［ε］は変形マトリックス、
‐［ε］−１は前記変形マトリックスの逆数、
‐ｄｉｖ［ε］は前記変形マトリックスとの差分を示す。）
【請求項４】
前記弾性画像の演算工程（６０）は、分水界分割（６１）のためのサブ工程を含み、該サブ工程は事前分割された画像における前記複数の区域の中から特定の区域が連続的に選択される工程であることを特徴とする請求項１乃至３のいずれかに記載の方法。
【請求項５】
前記区域は、事前分割された画像においてＨ密度が最も高い区域から、該密度が最も低い区域へと連続的に選択されることを特徴とする請求項４に記載の方法。
【請求項６】
前記弾性画像の前記演算工程（６０）は、最適化のためのサブ工程（６２）をさらに含み、該サブ工程は前記分割工程において夫々の新たなｎ＋１番目の区域が選択された後に行われることを特徴とする請求項４または５に記載の方法。
【請求項７】
前記最適化工程（６２）において、エラストグラムにおける実測された歪値と、前記選択されたｎ＋１個の区域における演算された歪値との歪誤差（ｓｔｒａｉｎ ｅｒｒｏｒ）と呼ばれる誤差を最小化することにより、前記分割工程において選択された前記新たなｎ＋１番目の区域の弾性を評価することを特徴とする請求項５に記載の方法。
【請求項８】
前記弾性画像の前記演算工程（６０）は、該演算工程を停止するための停止基準をそなえた前記弾性画像の比較工程（６３）を含むことを特徴とする請求項４乃至７のいずれかに記載の方法。
【請求項９】
前記停止基準は、ｎ個の区域及びｎ＋１個の区域における夫々の弾性を特徴づける２つの連続した工程に間において、演算された２つの歪誤差の差分が１０‐７を超えることであることを特徴とする請求項８に記載の方法。
【請求項１０】
前記停止基準は、前記選択された区域の数が１０を超えることであることを特徴とする請求項８に記載の方法。
【請求項１１】
アテロームプラーク破裂のリスク評価のための画像処理システムであって、該システムは、請求項１乃至１０のいずれかに記載の方法を実行するための手段（１、２、３）を有することを特徴とするアテロームプラーク破裂のリスク評価のための画像処理システム。
【請求項１２】
コンピューターで読み込み可能な媒体に記録されたプログラムコード指令を含むコンピュータープログラム製品であって、該コンピュータープログラム製品は、請求項１乃至１０いずれかに記載の方法の実行指令を含むことを特徴とするコンピュータープログラム製品。
【請求項１】
アテロームプラーク破裂のリスク評価のための画像処理方法であって、
該方法は、
‐血管組織が血圧機能により圧迫されることにより生じる内部変形を示す「エラストグラム」と呼ばれる受像工程（４０）と、
‐複数の区域を有する、事前に分割された画像を取得するために、前記エラストグラムを事前に分割する工程（５０）と、
‐前記エラストグラムの少なくとも１の領域における弾性を示す弾性画像を演算する工程（６０）であって、前記領域は前記事前分割された画像における複数区域から選択された区域に対応し、前記弾性画像を用いて使用者はアテロームプラークの破裂をリスク評価することが可能となる前記弾性画像の演算工程（６０）と
をそなえた方法であって、
‐第１に、前記血管組織が非圧縮性の弾性固体であるとみなして事前分割を行い、ラグランジュ乗数の勾配及び空間依存性ベクトルの勾配の合計と等しくなるように、ヤング率の弾性勾配を決定し、
‐第２に、前記事前分割は、前記複数の区域における夫々の区域の輪郭を前記空間依存性ベクトルにより決定する工程を含むアテロームプラーク破裂のリスク評価のための画像処理方法。
【請求項２】
前記複数の区域における夫々の区域の輪郭は、前記空間依存性ベクトルの最大値と一致することを特徴とする請求項１に記載の方法。
【請求項３】
前記空間依存性ベクトルＨが、以下の式で示されることを特徴とする請求項１または２に記載の方法。
Ｈ＝‐［ε］−１ ｄｉｖ［ε］
（‐［ε］は変形マトリックス、
‐［ε］−１は前記変形マトリックスの逆数、
‐ｄｉｖ［ε］は前記変形マトリックスとの差分を示す。）
【請求項４】
前記弾性画像の演算工程（６０）は、分水界分割（６１）のためのサブ工程を含み、該サブ工程は事前分割された画像における前記複数の区域の中から特定の区域が連続的に選択される工程であることを特徴とする請求項１乃至３のいずれかに記載の方法。
【請求項５】
前記区域は、事前分割された画像においてＨ密度が最も高い区域から、該密度が最も低い区域へと連続的に選択されることを特徴とする請求項４に記載の方法。
【請求項６】
前記弾性画像の前記演算工程（６０）は、最適化のためのサブ工程（６２）をさらに含み、該サブ工程は前記分割工程において夫々の新たなｎ＋１番目の区域が選択された後に行われることを特徴とする請求項４または５に記載の方法。
【請求項７】
前記最適化工程（６２）において、エラストグラムにおける実測された歪値と、前記選択されたｎ＋１個の区域における演算された歪値との歪誤差（ｓｔｒａｉｎ ｅｒｒｏｒ）と呼ばれる誤差を最小化することにより、前記分割工程において選択された前記新たなｎ＋１番目の区域の弾性を評価することを特徴とする請求項５に記載の方法。
【請求項８】
前記弾性画像の前記演算工程（６０）は、該演算工程を停止するための停止基準をそなえた前記弾性画像の比較工程（６３）を含むことを特徴とする請求項４乃至７のいずれかに記載の方法。
【請求項９】
前記停止基準は、ｎ個の区域及びｎ＋１個の区域における夫々の弾性を特徴づける２つの連続した工程に間において、演算された２つの歪誤差の差分が１０‐７を超えることであることを特徴とする請求項８に記載の方法。
【請求項１０】
前記停止基準は、前記選択された区域の数が１０を超えることであることを特徴とする請求項８に記載の方法。
【請求項１１】
アテロームプラーク破裂のリスク評価のための画像処理システムであって、該システムは、請求項１乃至１０のいずれかに記載の方法を実行するための手段（１、２、３）を有することを特徴とするアテロームプラーク破裂のリスク評価のための画像処理システム。
【請求項１２】
コンピューターで読み込み可能な媒体に記録されたプログラムコード指令を含むコンピュータープログラム製品であって、該コンピュータープログラム製品は、請求項１乃至１０いずれかに記載の方法の実行指令を含むことを特徴とするコンピュータープログラム製品。
【図１】
【図２】
【図３】
【図２】
【図３】
【公表番号】特表２０１２−５０９１２２（Ｐ２０１２−５０９１２２Ａ）
【公表日】平成２４年４月１９日（２０１２．４．１９）
【国際特許分類】
【出願番号】特願２０１１−５３６８８４（Ｐ２０１１−５３６８８４）
【出願日】平成２１年１１月２３日（２００９．１１．２３）
【国際出願番号】ＰＣＴ／ＥＰ２００９／０６５６１７
【国際公開番号】ＷＯ２０１０／０６６５６１
【国際公開日】平成２２年６月１７日（２０１０．６．１７）
【出願人】（５０９３５２５５４）ユニヴェルシテ ジョセフ フーリエ−グレノーブル アン (2)
【氏名又は名称原語表記】ＵＮＩＶＥＲＳＩＴＥ ＪＯＳＥＰＨ ＦＯＵＲＩＥＲ−ＧＲＥＮＯＢＬＥ １
【出願人】（５００５３１１４１）セントレ・ナショナル・デ・ラ・レシェルシェ・サイエンティフィーク (84)
【Ｆターム（参考）】
【公表日】平成２４年４月１９日（２０１２．４．１９）
【国際特許分類】
【出願日】平成２１年１１月２３日（２００９．１１．２３）
【国際出願番号】ＰＣＴ／ＥＰ２００９／０６５６１７
【国際公開番号】ＷＯ２０１０／０６６５６１
【国際公開日】平成２２年６月１７日（２０１０．６．１７）
【出願人】（５０９３５２５５４）ユニヴェルシテ ジョセフ フーリエ−グレノーブル アン (2)
【氏名又は名称原語表記】ＵＮＩＶＥＲＳＩＴＥ ＪＯＳＥＰＨ ＦＯＵＲＩＥＲ−ＧＲＥＮＯＢＬＥ １
【出願人】（５００５３１１４１）セントレ・ナショナル・デ・ラ・レシェルシェ・サイエンティフィーク (84)
【Ｆターム（参考）】
[ Back to top ]