JP4598824B2 - 信号解析法 - Google Patents

Description

本発明は信号を解析するための方法に関する。その解析は、持続時間制約を組み込む隠れマルコフモデリング技法を利用して、心電図（ＥＣＧ）又は他の生物医学信号のような信号を、その要素波形特徴に自動的にセグメンテーションすることを含むことができ、任意選択で、それと合わせて、セグメンテーションの精度及び頑強性を指示する、関連する信頼度指標を導出することを含む。

ＥＣＧ（頭字語ＥＫＧによっても知られている）は、心臓の電気的活動を測定する重要な非侵襲性の信号である。個々の心拍は多数の異なる心臓学的段階から構成され、それらの段階は、ＥＣＧ波形において１組の異なる特徴を生み出す。これらの特徴は、心臓の特定の領域における筋細胞の消極（放電）又は再分極（再充電）のいずれかを表す。図１は、人のＥＣＧ波形及び関連する特徴を示す。ＥＣＧ波形の標準的な特徴は、Ｐ波、ＱＲＳ群及びＴ波である。さらに、小さなＵ波が（Ｔ波の後に）時折存在する。

心周期はＰ波（その開始点及び終了点はＰ_ｏｎ及びＰ_ｏｆｆと呼ばれる）で開始し、それは心臓の心房消極(atrial depolarization)の周期に対応する。この後にＱＲＳ群が続いており、ＱＲＳ群は一般的に、ＥＣＧ波形の最もはっきりとした特徴であり、心室消極(ventricular depolarization)の周期に対応する。ＱＲＳ群の開始点及び終了点はＱ点及びＪ点と呼ばれる。ＱＲＳ群の後にＴ波が続いており、それは心室再分極の周期に対応する。Ｔ波の終了点はＴ_ｏｆｆと呼ばれ、心周期の終了を表す（Ｕ波が存在しないものと仮定する場合）。

ＥＣＧ信号を詳細に調べることによって、特徴的なＥＣＧ波形から多数の有益な測定値を導出することができる。その後、これらの測定値を用いて、患者の病状を評価することができ、心調律(cardiac rhythm)内にある任意の潜在的な異常を検出することができる。

特に重要な測定値は「ＱＴ間隔」であり、それはＱＲＳ群の開始からＴ波の終了までの時間、すなわちＴ_ｏｆｆ−Ｑと定義される。この時間間隔は、心室内の電気的活動（消極及び再分極の両方）の総持続時間に対応する。

ＱＴ間隔は、通常の間隔よりも長くなると、ＱＴ延長症候群（ＬＱＴＳ）の有効な指標となるので特に重要である。これは、顕在はしていないが、致命的な状態であり、その患者は、「トルサード・ド・ポワント(torsade de pointes)」として知られている非常に速く、異常な心調律（不整脈）に陥りやすくなる。このリズムが生じるとき、心臓は有効に拍動できなくなり、脳への血流が大幅に低下する。その結果、意識を突然失い、心臓死にいたる可能性もある。

別の重要な指標はＰＲ間隔であり、それは、Ｐ波の開始からＱＲＳ群の開始までの時間、すなわちＱ−Ｐ_ｏｎと定義される。これは、心房消極の開始から心室消極の開始までの時間に対応し、房室心ブロックのような状態を検出するために用いることができる。

ＱＴ間隔及びＰＲ間隔の測定及び評価に加えて、その要素特徴へのＥＣＧ波形の正確且つ頑強なセグメンテーションも、ＳＴ上昇のような状態を検出するために、さらには、連続した心拍内のＱＲＳ群のピーク（Ｒ点として知られている）の位置を特定することによって心拍数を計算するために重要である。

それゆえ、ＥＣＧ波形において特徴境界を正確に求めることに関する上記の重要性を踏まえて、この過程を自動化するために多大な努力がなされてきた。この分野における初期の仕事は、閾値処理及びテンプレートマッチングのような古典的な信号処理及びパターン認識手法に焦点を合わせていた。最近になって、隠れマルコフモデル（ＨＭＭ）に基づく統計的な手法が、ＥＣＧセグメンテーションの問題に適用されている。

隠れマルコフモデリングでは、根底を成す、対象となる状態シーケンスが存在するが、そのシーケンスは直に観測することはできない。一方では、さらに、この状態シーケンスに確率的に関連付けられる信号又は観測シーケンスも存在し、そこから、所与の状態シーケンスの最も確からしい値を推測することが望ましい。

隠れマルコフモデルの中心には、２つの確率関数があり、そのうちの一方は、特定の時刻ステップの直前の時刻ステップにおける状態の値を与えられるときの、その特定の時刻ステップにおける対象となる状態の確率に関連しており、もう一方は、その特定の状態値を与えられるときの信号又は観測値の確率に関連している。

隠れマルコフモデリング手法は、ＥＣＧ波形特徴の固有の統計的な特徴を利用するので、ＥＣＧ解析の場合に特に適している。さらに、そのモデルは、波形特徴の連続的な性質（すなわち、Ｐ波の後にＱＲＳ群が続き、その後にＴ波が続いている）も巧みに利用することができ、それは雑音のあるＥＣＧ信号のセグメンテーションを改善するのを助ける。

残念なことに、隠れマルコフモデルは、極めて信頼度の低いセグメンテーションを引き起こす傾向があるという重大な短所を抱えている。これらのセグメンテーションは、セグメント化された特徴が典型的には少数の時刻サンプルしか続かないという特徴を有する。この問題は、非常に短い持続時間を有する特徴に向かってモデルが偏ることから起こり、それゆえ、標準的なＨＭＭは、高い度合いの頑強性を必要とするＥＣＧセグメンテーションのような応用形態には適していない。

図２及び図４は、専門家であるＥＣＧ解析者（analyst）によって求められるような波形境界（垂直な実線）及び標準的な隠れマルコフモデルを用いて推測される波形境界（垂直な破線）とともに、２つのサンプルＥＣＧ波形を示す。いずれの事例でも、元の信号内にそれぞれ１つしか存在しないにもかかわらず、そのモデルは２つのＱＲＳ群、及び２つのＴ波を見いだしている。詳細には、モデルセグメンテーションは、非常に短いＱＲＳ群及びＴ波を含み、それは、Ｐ波の開始直後、及び本当のＱＲＳ群位置の直前に生じる。そのような「二重拍動」セグメンテーションは、モデルの性能及び信頼度に深刻な影響を及ぼす。

この背景に留意して、本発明の１つの目的は、ＥＣＧ信号のような信号のセグメンテーションの信頼度を高めることに関する。

本発明のさらに別の目的は、そのモデルから信頼度指標を導出することである。この信頼度指標を用いて、セグメンテーションの品質及び頑強性を評価することができ、さらに、自動化システムによる信頼度のあるセグメンテーションが不可能であるほど雑音が多すぎるか、又は著しく異常である、いかなるＥＣＧ波形をも検出することができる。

本発明の一態様によれば、隠れマルコフモデルを用いた生物医学信号のセグメンテーションのための、コンピュータによって実施される方法であって、当該モデルは複数の状態を含み、
当該方法は、
状態のうちの少なくとも１つについて、最小持続時間制約ｄ_ｍｉｎを指定するステップと、
指定された最小持続時間を有するモデル内の状態毎に、当該状態を１組のサブ状態で置き換えるステップであって、当該サブ状態の総数は最小持続時間ｄ_ｍｉｎの値に等しい、置き換えるステップと、
１組のサブ状態を連結するステップであって、それによって、左−右マルコフ連鎖を形成し、最初のｄ_ｍｉｎ−１個のサブ状態はそれぞれ、０の自己遷移確率と、その右隣の状態に遷移する１の遷移確率と、モデル内の任意の他の状態に遷移する０の遷移確率とを有する、連結するステップと、
モデルを、生物医学信号を表すデータに適用するステップであって、それによって、信号の状態へのセグメンテーションについて情報を得る、モデルを適用するステップと
を含む方法が提供される。

本発明の別の態様によれば、ビタビアルゴリズムに対する変更形態を用いて、観測のシーケンスを含む信号を、有限状態離散時間マルコフ過程の状態のシーケンスにセグメンテーションするための、コンピュータによって実施される方法であって、変更形態は、
ビタビアルゴリズムに組み込まれる有限状態離散時間マルコフ過程の少なくとも１つの状態について持続時間変数及び持続時間制約を定義することであって、持続時間制約は当該少なくとも１つの状態について最小持続時間を指定する、定義することと、
信号に変更されたビタビアルゴリズムを適用することであって、それによって、観測のシーケンスを説明する最も確からしい持続時間制約付き状態シーケンスを計算する、適用することと
を含む方法において、
各時刻ステップにおいて、有限状態離散時間マルコフ過程内の状態毎に、その時刻ステップまでの観測のシーケンスを説明するとともに状態で終わっている、最も確からしい状態シーケンスを計算する際に、
持続時間制約を有する状態毎に、その状態のための持続時間変数を用いて、その状態だけから成るとともに直前の時刻ステップにおいてその状態で終わっている、一連の先行状態のシーケンスの長さを絶えず注視し、
その状態のための持続時間変数が、その状態のための指定された持続時間制約以上である場合には、所与の時刻ステップにおける状態シーケンス計算において、その状態からマルコフ過程内の任意の他の所与の状態への遷移が考慮され、
その状態のための持続時間変数が、その状態のための指定された持続時間制約より小さい場合には、所与の時刻ステップにおける状態シーケンス計算において、その状態からマルコフ過程内の任意の他の状態への遷移が考慮されず、
かつ、方法において、
所与の時刻ステップまで１組の最も確からしい状態シーケンスを計算した後に、その特定の状態だけから成るとともに考慮されたばかりの時刻ステップにおいてその状態で終わっている、一連の先行状態シーケンスの長さを絶えず注視するために、持続時間制約を有する状態毎に持続時間変数を更新する、方法が提供される。

本発明のさらなる態様は、確率的なセグメンテーションアルゴリズムによってセグメンテーションされた信号を解析するための、コンピュータによって実施される方法であって、
複数のセグメンテーションされた信号特徴毎に信頼度指標を計算することと、
個々のセグメンテーションされた信号特徴長に対して信頼度指標をプロットすることと、
密度モデリング技法を適用することであって、それによって、高い信頼度特徴に関連付けられるデータ空間の適切な領域を求める、適用することと、
特定の信号特徴のための信頼度指標がこの領域の外側に該当するか否かを判定することと
を含む方法を提供する。

ここで、一例にすぎないが、添付の図面を参照しながら、本発明の実施形態を記述する。

本発明の以下に記述される実施形態は、ＥＣＧ信号の特徴のうちの１つ又は複数を認識するようにトレーニングされる隠れマルコフモデル（ＨＭＭ）の基礎を基にする。そのモデルは複数の状態から構成され、各状態はＥＣＧ信号の特定の領域を表す。ＥＣＧセグメンテーションのための隠れマルコフモデルの１つの形態のアーキテクチャ又は「トポロジ」の図式的な表示が図７ａに示される。このモデルは６つの固有の状態から構成され、それらの状態は、Ｐ波、Ｐ波の終了とＱＲＳ群の開始との間の基線の部分（「基線１」と呼ばれる）、ＱＲＳ群、Ｔ波、Ｕ波及びＴ波（又はもし存在する場合にはＵ波）の終了と後続の心拍のＰ波の開始との間の基線の部分を表す。ＥＣＧ波形のこれらの部分は図１に示される。ＥＣＧセグメンテーションのための隠れマルコフモデルアーキテクチャの多数の代替の形態が図７ｂ〜図７ｅに示される。それぞれの場合に、そのモデルは、「隠れ」状態シーケンス（白抜きのノードによって示される）から構成され、そのシーケンスは観測される信号（陰影付きのノードによって示される）に確率論的に関連付けられる。ＥＣＧセグメンテーションの場合、隠れ状態ｓ_ｔは、時刻ｔにおいて活動している特定の波形特徴に対応し、観測される信号サンプルＯ_ｔは、ＥＣＧ波形の関連する信号サンプルに対応する。

そのモデルは３つのパラメータ、すなわちｉ）初期状態分布π、ｉｉ）遷移行列Ａ及びｉｉｉ）ＨＭＭの状態ｋ毎の観測確率モデルｂ_ｋによって支配される。ＨＭＭをトレーニングすることは、特定の誤差関数を最小にするために、これらのパラメータを調整することを含む。

そのモデルパラメータのために適した１組の値を与えられるとき、ＨＭＭはＥＣＧ信号のための生成確率モデルの形をとる。より厳密には、そのモデルはＥＣＧ波形を生成するための確率的な手順を定義する。この手順は初期状態の選択で開始し、それは、初期状態分布πからのサンプリングによって達成される。その際、選択された特定の初期状態は、ＥＣＧ信号の第１の波形特徴（すなわち、時刻サンプルｔ＝１において活動している特徴）に対応する。その後、この波形特徴のための関連するＥＣＧ信号サンプルが、その特定の状態のための対応する観測確率モデルｂ_ｋからのサンプリングによって生成される。

その手順における次のステージは、次のモデル状態、すなわち時刻サンプルｔ＝２において活動している波形特徴を確率的に選択することである。これは、遷移行列からのサンプリングによって達成される。この行列は、そのモデル内の状態毎に、次の状態が生じる確率を定義する簡単な表の形をとる。一旦、このようにして次の状態が選択されたなら、そのモデルは、適切な観測確率モデルからのサンプリングによって、関連するＥＣＧ信号サンプルを生成する。この手順は、残りの時刻ステップについて連続的に繰り返される。

特定のＥＣＧ信号を与えられるとき、その信号を、モデルアーキテクチャによって定義される、対象となる種々の特徴又は領域にセグメント化すること（すなわち、各モデル状態の信号境界の位置を特定すること）は、上記の「生成手順」を逆にすることを含む。したがって、その目的は、所与のＥＣＧ信号を引き起こす隠れ状態の最も確からしいシーケンスを推測することである。その際、この状態シーケンスは、ＥＣＧ信号内の対象となる領域の境界を、それゆえ対応するセグメンテーションを定義する。

このセグメンテーション手順は、ビタビアルゴリズムを用いることによって達成することができる。より厳密には、L. R. Rabiner著「A Tutorial on Hidden Markov Models and Selected Application in Speech Recognition」（Proceedings of the IEEE, 77: 257-286, 1989）に記述されるように、ビタビアルゴリズムは、トレーニングされた隠れマルコフモデルとともに所与のＥＣＧ信号を入力として取り込み、最も確からしい隠れ状態シーケンスを返す。

そのモデルがＥＣＧ波形のセグメンテーションという作業を実行するのに成功するために、ＥＣＧの、頑強でかつ情報を与える符号化を用いる必要がある。特に有効な表現は、N. P. Hughes、L. Tarassenko及びS. J. Roberts著「Markov Models for Automated ECG Interval Analysis」（Advances in Neural Information Processing Systems 16, MIT Press, 2003）に記述されるような、ＥＣＧ信号のアンデシメーテッド(undecimaetd)ウェーブレット変換（ＵＷＴ）からの１組の係数である。

一旦、ＥＣＧのための表現が選択されたなら、ＨＭＭをトレーニングすることは、ＥＣＧデータの確率（又は尤度）、すなわちｐ（Ｏ｜λ）を最大にするために、そのモデルのパラメータを調整することを含む。ただし、ＯはＥＣＧ信号のデータセットを表し、λは１組のＨＭＭパラメータを表す。トレーニングは、教師あり(supervised)で、又は教師なし(unsupervised)で達成することができる。

教師あり学習の場合、トレーニングデータは、多数のＥＣＧ波形および対応する状態ラベルから成る。したがって、所与のＥＣＧ波形からの信号サンプルＯ_ｔ毎に、そのサンプルがどの特定の状態（すなわち、波形特徴）に属するかを指示するラベルｌ_ｔも入手できるようにする。その後、以下の最尤推定を用いて、初期状態分布及び遷移行列を計算することができる。

ただし、ｍ_ｉは、第１の信号サンプル（全てのＥＣＧ波形にわたる）が状態ｉに属する全回数であり、ｎ_ｉｊは全てのラベルシーケンスにおける状態ｉから状態ｊへの全遷移数であり、Ｎはデータセット内の個々のＥＣＧ信号の数である。

特定の状態に属する全ての信号サンプルを抽出し、その後、観測モデルをそれらのサンプルに直に当てはめることによって、その状態ｋ毎の観測確率モデルｂ_ｋを学習することができる。ガウス観測モデルの場合には、これは単に、与えられたサンプルの平均及び分散を計算することを含む。ガウス混合モデルの場合には、与えられたサンプルから観測モデルのパラメータを学習するために、ＥＭアルゴリズム（期待値最大化アルゴリズム）を用いる必要がある。

教師なし学習の場合、トレーニングデータは、ＥＣＧ波形だけから成る。その後、バウム−ウェルチアルゴリズム（それは、より一般的なＥＭアルゴリズムの特殊な場合である）を用いて、隠れマルコフモデルをトレーニングすることができる。

そのモデルが、信頼度のないセグメンテーションを引き起こすのを防ぐために、本発明の実施形態は、そのモデルに１組の持続時間制約を課す。これらの制約は、そのモデルが生理的に（関連する波形特徴の持続時間に関して）起こりそうもないセグメンテーションを推測する能力を制限する。

詳細には、持続時間制約は、モデル内の状態毎にただ１つの数ｄ_ｍｉｎの形をとり、それは、その特定の状態のための（「サンプル」内の）最小許容持続時間を表す。これらの持続時間制約は、多数の異なる方法で推定することができる。ラベル付きのＥＣＧデータが入手可能である場合には、波形特徴毎の持続時間制約は、データセット内に存在するその特徴の最小持続時間の適切な割合として簡単に推定することができる。たとえば、最小Ｔ波持続時間は、データセット内に存在する最小Ｔ波持続時間の８０％に設定することができる。

ラベル付きのデータが入手できない場合に、且つ使用できるものが生のＥＣＧデータだけである場合には、持続時間制約は、そのデータから学習されなければならない。これは、トレーニングされたモデルを用いてＥＣＧデータをデータセットにセグメンテーションし、その後、結果としてもたらされたセグメンテーションから（ラベル付きのデータの場合に先に説明されたのと同じようにして）持続時間制約を推定することによって達成することができる。この場合には、最初に信頼度のない「二重拍動」セグメンテーションをすべて除去する必要があり、それを除去しなければ、推定された持続時間制約の品質に悪影響を及ぼすであろう。

一旦、持続時間制約が計算されたなら、それらを用いて、隠れマルコフモデルのセグメンテーション性能を改善することができる。これは、持続時間制約をビタビアルゴリズムに組み込むことによって、又は持続時間制約をＨＭＭアーキテクチャに組み込むことによって達成することができる。ここで、これらの両方の方法を詳細に検討する。

ビタビアルゴリズムは、L. R. Rabiner著「A Tutorial on Hidden Markov Models and Selected Applications in Speech Recognition」（Proceedings of the IEEE, 77: 257-286, 1989）に詳述されるように、隠れマルコフモデルを実際にうまく利用できるようにする重要な推測手順である。ＥＣＧ解析との関連では、ビタビアルゴリズムによれば、所与のＥＣＧ波形におけるＰ波、ＱＲＳ群及びＴ波のための最も確からしい境界位置（すなわち、オンセット及びオフセット）を特定できるようになる。

持続時間制約をビタビアルゴリズムに組み込むためには、関連する波形特徴のそれぞれの持続時間が所与の持続時間制約に従うように、そのアルゴリズムによって返される状態シーケンスを制限する必要がある。

より形式的には、ＨＭＭλ（Ｋ個の状態を有する）及びＥＣＧ信号Ｏ＝Ｏ_１Ｏ_２．．．Ｏ_Ｔを与えられるとき、ｋ＝１．．．Ｋ毎に持続時間制約ｄ_ｍｉｎ（ｋ）を満たす最も確からしい状態シーケンスＳ＝ｓ_１ｓ_２．．．ｓ_Ｔを推測したい。各持続時間制約は、状態ｋの任意の「実行(run)」の持続時間が少なくともｄ_ｍｉｎ（ｋ）サンプルでなければならないという要件を具現する。それゆえ、その目的は、（個々の状態実行の持続時間に関して）無効であると予めわかっている状態シーケンスを推測しないようにすることを確実にすることである。

標準的なビタビアルゴリズムは、最初のｔ回の観測を説明するとともに状態ｉにおいて終了する、最も確からしい状態シーケンスの尤度を計算するための効率的な手順に基づく。

そのアルゴリズムの中心には、δ_ｔ（ｉ）を計算するための以下の漸化関係がある。

ただし、上記の関係は、状態シーケンスが一次マルコフであるという事実を利用している（すなわち、δ_ｔ（ｉ）を計算する際に、直前の状態しか考慮に入れる必要がない）。その際、任意の時点ｔにおいて、関連付けられた遷移確率が０でないなら、状態トレリスの中の最も確からしいパスは、その現在の状態から異なる状態に常に遷移することができる。したがって、持続時間制約をビタビアルゴリズムに組み込むために、このパスは、対象となる現在の状態が最小許容持続時間以上の時間にわたって占有されている場合にのみ、１つの状態から別の状態への任意の遷移が許されるように制限されなければならない。

本発明のこの第１の実施形態による持続時間制約付きビタビアルゴリズムは以下のように進む。漸化式の各時刻ステップｔにおいて、時刻ｔ−１において終了する、状態ｉの現在の「実行」の持続時間ｄ（ｉ）を絶えず注視する。より厳密には、状態毎の持続時間変数を用いて、その特定の状態だけから成るとともに直前の時刻ステップにおいてその状態で終わっている、一連の先行状態シーケンスの長さを絶えず注視する。この持続時間が、最小許容持続時間ｄ_ｍｉｎ（ｉ）以上である場合には、δ_ｔ（ｊ）を計算するときに、この状態から任意の他の状態ｊへの遷移を考慮する（標準的なビタビ漸化式の場合と同様）。しかしながら、現在の状態の占有が最小許容値より低い場合には、この状態から異なる状態への遷移は許されず、その状態は（その最小持続時間制約が満たされるまで）自己遷移しかできない。δ_ｔ（ｉ）がｔの全ての値について計算されると、標準的な後戻り(backtracking)法を用いて、最終的な状態シーケンスを生成することができる。

持続時間制約付きビタビアルゴリズムの手順全体が、図６の流れ図において概略的に示される。そのアルゴリズムの第１のステージは、ＨＭＭパラメータ、観測された信号、及び状態毎の最小持続時間制約を入力することである。これは、ボックス１０に示される。

ボックス２０において示される、そのアルゴリズムの次のステージは、対象となる変数の初期化に関する。より厳密には、状態ｉ毎に、δ_１（ｉ）＝π_ｉｂ_ｉ（Ｏ_１）、ψ_１（ｉ）＝０及びｄ（ｉ）＝１が設定される。標準的なビタビ手順と同じように、変数ψ_ｔ（ｉ）は、そのモデルが時刻ｔにおいて状態ｉにあった場合に、時刻ｔ−１における最も確からしい状態を格納するために用いられる。変数ｄは持続時間制約付きビタビ手順に固有であり、時刻サンプルｔ−１において終了する、各状態の現在の「実行」の持続時間を絶えず注視するために用いられる。より厳密には、状態ｉ毎に、持続時間変数ｄ（ｉ）を用いて、状態ｉだけから成るとともに時刻ステップｔ−１においてその状態ｉにおいて終了する、一連の先行状態シーケンスの長さが絶えず注視される。

持続時間制約付きビタビアルゴリズムの主要部がボックス３０及び４０において示される。これら２つの手順は、２〜Ｔの全ての時刻サンプルについて連続的に繰り返される（ただし、Ｔは信号の長さである）。ボックス３０は、所与の時刻ステップｔにおける状態ｉ毎のδ_ｔ（ｉ）及びψ_ｔ（ｉ）の値を計算する。標準的なビタビアルゴリズムでは、δ_ｔ（ｉ）及びψ_ｔ（ｉ）を計算する際にそれぞれ用いられる「ｍａｘ」演算子及び「ａｒｇｍａｘ」演算子が、全ての可能な直前の状態（１〜Ｋ）について評価される。しかしながら、持続時間制約付きビタビアルゴリズムでは、これらの演算子は、それらの最小持続時間制約ｄ_ｍｉｎを満たす状態についてのみ評価される。さらに、それらの演算子は、最も確からしい直前の状態を現在考慮中の状態とすることが常に許されるように、現在の状態ｉについても評価される。

これを仮定すると、それゆえ、ボックス３０の最初のステップは、その最小持続時間制約を満たす１組の状態ｊ（デフォルトで状態ｉを含む）を計算することである。すなわち、それらの状態ｋはｄ（ｋ）≧ｄ_ｍｉｎ（ｋ）を満たす。その後、δ_ｔ（ｉ）及びψ_ｔ（ｉ）の計算がこの組の状態ｊについてのみ評価される。

持続時間制約付きビタビアルゴリズムにおける次のステージは、状態毎の持続時間変数を更新することである。より厳密には、最も確からしい直前の状態ψ_ｔ（ｉ）も状態ｉである場合には、この状態のための持続時間変数が１だけインクリメントされる。すなわちｄ（ｉ）＝ｄ（ｉ）＋１に設定される。これは、「自己遷移」の場合に対応する。その場合、そのモデルは連続した時刻ステップにわたって同じ状態のままであり、それゆえ、時刻サンプルｔにおいて終了する状態ｉの現在の「実行」の持続時間が１だけ増加する。しかしながら、最も確からしい直前の状態ψ_ｔ（ｉ）が状態ｉでない場合には、この状態の持続時間変数は１にリセットされる。すなわちｄ（ｉ）＝１に設定される。これは、異なる状態から状態ｉへの遷移の場合に対応し、それゆえ、状態ｉの自己遷移の現在の「実行」を終了する。

ボックス３０及び４０が時刻サンプルＴまで計算されたなら、標準的なビタビ後戻り手順によって、以前に計算されたδ変数及びψ変数を用いて最も確からしい状態シーケンスを抽出することができる。この後戻り手順の最初のステップは、ボックス５０に示されるように、δ_Ｔの値を最大にする、すなわち以下の式を満たす状態を見つけることである。

この値が計算されたなら、ψ_Ｔ（ｓ_Ｔ）の値から、直前の時刻ステップ（すなわち、Ｔ−１）において最も確からしい状態を簡単に見つけることができる。その後、この帰納的な手順を用いて、全ての残りの時刻サンプルについて最も確からしい状態値を計算することができる。

本発明の第２の実施形態によれば、ＨＭＭのセグメンテーション性能を改善する代替の方法は、持続時間制約をモデルアーキテクチャに直に組み込むことである。このようにして、一旦、制約が組み込まれたなら、その後、標準的なビタビアルゴリズムを用いて、所与のＥＣＧ波形の最適な持続時間制約付きセグメンテーションを推測することができる。

持続時間制約をＨＭＭアーキテクチャに組み込むことは、そのモデルの各状態を１組のｄ_ｍｉｎ個の新たな状態で置き換えることを含む。これらの状態は、単純に左右に連結され、共通の観測密度を共有する。

より厳密には、そのモデルは、ｄ_ｍｉｎ（ｋ）（ただし、ｄ_ｍｉｎ（ｋ）>１）の最小持続時間を有する状態ｋ毎に、元の状態ｋの直前にあるｄ_ｍｉｎ（ｋ）−１個の付加的な状態で拡張される。付加的な状態はそれぞれ、０の自己遷移確率と、その右隣の状態に遷移する１の遷移確率と、そのモデル内の任意の他の状態に遷移する０の遷移確率とを有する。こうして、まとめて考えると、これらの状態は単純な左−右マルコフ連鎖を形成し、その連鎖内の各状態は、（その連鎖を用いる任意の実行中に）長くとも１つの時刻サンプルの間だけ占有される。

この連鎖の最も重要な特徴は、状態毎の観測確率モデルのパラメータが元の状態ｋのパラメータと同じであることである。このようにモデルパラメータを共有(share)することは一般的に、「共有化(tying)」と呼ばれる。こうして、元の状態ｋの代わりに用いるｄ_ｍｉｎ個の新たな状態に関連付けられる観測は、元の状態と同じ組のパラメータによって支配される。

持続時間制約をＨＭＭアーキテクチャに組み込むための手順全体が図８に図式的に示される。その図の左側には、元の隠れマルコフモデルアーキテクチャからの状態（図７に示される）のうちの１つがある。その図の右側には、この状態が持続時間制約を組み込むことによって変形される新たな構造がある。詳細には、ｄ_ｍｉｎ（ｋ）−１個の新たな状態（自己遷移ループを持たない）は、単純に左から右へと元の状態に連結され、これらの状態のための観測確率モデルのパラメータは元の状態の対応するパラメータに「共有化される」。この後者の手順は、その図の破線のボックスによって指示される。

上記の過程は、元の隠れマルコフモデル内の元の状態（すなわち、Ｐ、Ｂ１、ＱＲＳ、Ｔ等）毎に繰り返される。したがって、持続時間制約をモデルアーキテクチャに組み込む全体的な効果は、そのモデル内の元の単一の状態がそれぞれ、新たな「持続時間制約」構造によって置き換えられることである。図９は、結果として構成される持続時間制約ＨＭＭアーキテクチャの一部を示す。そのモデル内の各波形特徴（たとえば、Ｐ波、ＱＲＳ群等）が、ここでは複数の状態に関連付けられる。その際、この複数の状態は、ＥＣＧ波形特徴毎の最小持続時間要件をカプセル化（encapsulate）する。

持続時間制約が、先に説明されたようにしてＨＭＭに組み込まれたなら、その後、そのモデルは、標準的なビタビアルゴリズムを用いた新たなＥＣＧ波形のセグメンテーションに用いることができる。その際、そのモデルはもはや、その持続時間が許容される持続時間より短い状態シーケンスを生成することができないので、結果として生成されるセグメンテーションは、持続時間制約に従うことを保証される。

図３及び図５は、図２及び図４において先に検討されたのと同じＥＣＧ波形を示す。図３及び図５では、専門家であるＥＣＧ解析者によって求められるような波形境界は、垂直な実線によって示されており、本発明の一実施形態による、持続時間制約を利用する隠れマルコフモデルを用いて自動的に求められる状態境界は、垂直な破線として示される。そこでは、モデルと専門家である解析者とが概ね一致するだけでなく、図２及び図４に示されるような、従来モデルの信頼度のない「二重拍動」セグメンテーションが除去されており、結果的な解析は、はるかに高い頑強性を有する。それらの境界は、特定の状態の持続時間、たとえばＱＴ間隔を求める等の、さらに別の解析に、再び自動的にかけることができ、一連のパルスを含む信号にわたる、そのような間隔の変化の統計的な解析を実行することができる。それらの結果は、生の境界情報又はさらに別の解析の結果の形で、たとえば、コンピュータシステムのディスプレイ上に表示することによって出力される。

本発明の第３の実施形態によれば、ＥＣＧ波形の正確なセグメンテーションに加えて、トレーニングされたモデルからの「信頼度指標」（又はスコア）を生成することもできる。この信頼度指標は、結果として生成されたセグメンテーションにおいて、そのモデルがどの程度信用できるかを定量化し、それを用いて、セグメンテーションされたＥＣＧ波形の品質を評価することができる。それゆえ、この指標を用いるとき、雑音が多すぎるか、又は「正常でない」（波形形態に関して）ＥＣＧ波形を自動的に検出し、信頼度のあるセグメンテーション、すなわちＱＴ間隔の信頼度のある推定値を提供することができる。信頼度指標は、確率的に動作する任意のセグメンテーションアルゴリズムとともに用いることができる。

信頼度指標は、モデルを与えられるときの信号の対数尤度(log likelihood)、すなわちｌｏｇ ｐ（Ｏ_１Ｏ_２．．．Ｏ_Ｔ｜λ）と定義することができるか、又は（モデルを与えられるときの）信号及び最も確からしい状態シーケンスｓ_１ｓ_２．．．ｓ_Ｔの結合(joint)対数尤度、すなわちｌｏｇ ｐ（Ｏ_１Ｏ_２．．．Ｏ_Ｔ，ｓ_１ｓ_２．．．ｓ_Ｔ｜λ）と定義することができる。対数尤度は、生の尤度とは対照的に、結果として生成される値をコンピュータのダイナミックレンジ内に保持するために用いられる。

信頼度指標のいずれの定義も、隠れマルコフモデルの場合に、さらには隠れセミマルコフモデル（L. Thoraval、G. Carrault及びF. Mora著「Continuously variable durations HMMs for ECG segmentation」（Proceedings IEEE-EMBS, 1992, pp. 529-530）に記述される）、階乗(factorial)隠れマルコフモデル（Z. Ghahramani及びM. I. Jordan著「Factorial hidden Markov models」（Machine Learning, Vol. 29, Issue 2, Nov. 1997）に記述される）のような、これらのモデルのさらに高度な形式の場合に、効率的に計算することができる。最初の形式は、そのモデルを与えられるときの信号の対数尤度と定義され、以下の式によって与えられる、標準的なＨＭＭフォワード変数(forwards variable)から計算することができる。

結果として生成される値をコンピュータのダイナミックレンジ内に保持するために、各時刻ステップにおいてこの変数をスケーリングするのが一般的である。これは、フォワード変数を、以下の式によって与えられるスケーリング係数で割ることによって達成される。

その際、そのモデルを与えられるときの信号の対数尤度は、以下のように簡単に計算することができる。

第２の形式の信頼度指標は、（モデルを与えられるときの）信号と最も確からしい状態シーケンスとの結合対数尤度と定義され、ビタビアルゴリズムとの関連で先に説明されたように、デルタ変数δ_ｔ（ｉ）から計算することができる。同様に、以下の式によって与えられるスケーリング係数を用いて各時刻ステップにおいてこの変数をスケーリングするのが一般的である。

その際、（モデルを与えられるときの）信号と最も確からしい状態シーケンスとの対数尤度は、以下のように簡単に計算することができる。

上記のスケーリング手順を用いることへの代替形態は、L. R. Rabiner著「A Tutorial on Hidden Markov Models and Selected Applications in Speech Recognition」（Proceedings of the IEEE, 77: 257-286, 1989）に記述されるように、ｌｏｇ α_ｔ（ｉ）又はｌｏｇ δ_ｔ（ｉ）の値を直に計算することである。この場合、信頼度指標は、帰納中に直に計算され、さらに別の処理は不要である。

多くの状況において、所与のＥＣＧ信号は、対象となる多数のＥＣＧ特徴（多数のＥＣＧ拍動又は多数のＱＴ間隔等）を含むことになり、それゆえ、信号内の対象となる個々の特徴毎に個別の信頼度指標を計算することが望ましい（信号全体のための集団の指標ではない）。時刻サンプルｔ１で開始し、時刻サンプルｔ２において終了する特定の特徴について、その特徴のための信頼度指標は、時刻ｔ２までの信号部分に対する信頼度指標から、時刻ｔ１まで（ただし、ｔ１を含まない）の信号部分に対する信頼度指標を減算すること、すなわち以下の式によって、簡単に計算することができる。

又は、信頼度指標の第２の定義を用いて、以下の式によって計算することができる。

この信頼度指標を実際に用いて、特定の特徴のセグメンテーションにおいて、そのモデルが有する信頼度を評価するときには、解析中に、問題となる特徴の長さが考慮に入れられなければならない。信頼度指標の形式がいずれも特定のＥＣＧ特徴の長さとともに変化するので、これが必要とされる。

ＥＣＧ特徴の長さを信頼度指標解析に組み込む最も実効的な方法は、単に、多数の正常なＥＣＧ特徴（すなわち、異常と見なされないか、又は雑音によって汚染されていない特徴）の場合に、特徴長（ｘ軸上）に対する信頼度指標（ｙ軸上）のプロットを作成することである。図１０は、そのようなプロットを示しており、対象となる特定の特徴は、個々のＥＣＧ波形（すなわち、１つの心拍全体に対応するＥＣＧ信号の部分）であると見なされる。

ＥＣＧ特徴長に対する信頼度指標のプロットを与えるとき、標準的な密度モデリング技法を用いることによって、高い信頼度指標に関連付けられるこの２−Ｄ空間の領域を特定することができる。詳細には、正常な（すなわち高い信頼度の）ＥＣＧ特徴に対応する１組の２−Ｄデータ点は、ガウス分布又は確率的な線形回帰モデルのような、任意の適当な密度関数によってモデル化することができる。

所与の密度関数は自然に、この２−Ｄデータ空間内の高い信頼度の領域のための境界を定義する。ガウス分布の場合、この境界は、その平均から２つの標準偏差内に存在するデータ空間の領域であると見なすことができる。確率的な線形回帰モデルの場合、この境界は、適当な統計的有意水準にある平均回帰線のための下限信頼度よりも上に存在するデータ空間の領域であると見なすことができる（R. J. Freund及びW. J. Wilson著「Regression Analysis」（Academic Press, 1998）に記述される）。この下限（９９％の有意水準にある）は、図１０において破線によって示される。

セグメンテーションされた特徴のそれぞれの長さとともに、多数の新たなＥＣＧ特徴のセグメンテーションのための信頼度指標を与えられるとき、任意の信頼度指標が高い信頼度領域に対する境界の外側に該当するか否かを判定することは簡単なことである。その際、それらの信頼度指標がこの領域の外側に該当する特徴のためのセグメンテーションの結果は、最終的な解析から除去することができるか、又は別法では、これらの特徴は、訓練を受けたＥＣＧ解析者によってさらに手動で解析するためにハイライトすることができる。

図１１はＥＣＧ信号の一部を示し、それは、人為的な筋雑音によって汚染された２つのＥＣＧ波形を含む。またその図は、持続時間制約を組み入れた隠れマルコフモデルによる各波形のセグメンテーションも示す。その２つの波形は、下側のｘ軸上で「Ａ」及び「Ｂ」を付される。

残念なことに、人為的な筋雑音は、標準的なＥＣＧ信号の帯域幅に類似の帯域幅を有し、それゆえ、簡単なフィルタリング技法によって除去することはできない。その雑音は、ＥＣＧ信号の根底を成す波形特徴を隠し、波形の正確なセグメンテーションを非常に難しくする。このため、そのような信頼度のないセグメンテーションを検出できることが重要である。

図１２は、雑音のある２つのＥＣＧ波形の場合の、波形の長さに対する信頼度指標のプロットを示す。９９％下限信頼度（破線）も示される。雑音のある２つの波形のための信頼度スコアはいずれも、この線の下側にあり、それゆえ、これらの波形が、信号の自動解析中に自動的に検出されるようにする。

本発明の上記の実施形態は、ＥＣＧ信号解析との関連で記述されてきた。しかしながら、本発明はこの応用形態には限定されず、他の生物医学信号の解析のためにも用いることができる。適する生物医学信号の例には、脳波図（ＥＥＧ）信号（その場合には、解析は無呼吸(apneoa)を検出するために用いることができる）；筋肉の電気的活動を測定する筋電図（ＥＭＧ）；心臓電気生理現象の場合のような侵襲的に得られる心信号（その場合には、カテーテルが心臓に挿入され、各カテーテルは多数の電極を含む）；及びペースメーカの移植前に用いられる、人体に電極を挿入することによって心臓から導出される信号を含む心内(endocardial)電位図信号が含まれる。持続時間制約付きビタビアルゴリズム及び信頼度指標技法は、隠れマルコフモデルとともに用いることには限定されない。たとえば、持続時間制約付きビタビアルゴリズムは、任意の有限状態離散時間マルコフ過程とともに用いることができ、信頼度指標は確率的に動作する任意のセグメンテーションアルゴリズムとともに用いることができる。持続時間制約付きビタビアルゴリズムは、生物医学以外の信号とともに、たとえば誤り検出／訂正符号の分野において用いることもできる。

本発明の実施形態は、コンピュータシステム上で実行されるコンピュータプログラムによって実行される。そのコンピュータシステムには、任意のタイプのコンピュータシステムを用いることができるが、一般的には、任意の適当な言語において書かれるコンピュータプログラムを実行する従来のパーソナルコンピュータである。そのコンピュータプログラムはコンピュータ読取り可能媒体に格納することができ、そのコンピュータ読取り可能媒体は、任意のタイプであってよい。たとえば、コンピュータシステムのドライブに挿入することができ、磁気的、光学的又は光磁気的に情報を格納することができる、ディスク形媒体のような記録媒体；ハードドライブのような、コンピュータシステムの固定記録媒体；又は固体コンピュータメモリである。解析されるべき信号は、たとえばＥＣＧ装置から直にコンピュータシステムに入力することができるか、又はコンピュータシステムが、予め得られている信号の記憶装置から、信号を表す情報を読み出すことができる。

十分に標識が付けられたＥＣＧ波形を示す図である。 専門家であるＥＣＧ解析者によって求められるような波形境界、及び標準的な隠れマルコフモデル手法を用いて推測される波形境界とともにＥＣＧ波形を示す図である。 専門家であるＥＣＧ解析者によって求められるような波形境界、及び新規の持続時間制約手法を用いて推測される波形境界とともに、図２において用いられるのと同じＥＣＧ波形を示す図である。 専門家であるＥＣＧ解析者によって求められるような波形境界、及び標準的な隠れマルコフモデル手法を用いて推測される波形境界とともに、第２のＥＣＧ波形を示す図である。 専門家であるＥＣＧ解析者によって求められるような波形境界、及び新規の持続時間制約手法を用いて推測される波形境界とともに、図４において用いられるのと同じＥＣＧ波形を示す図である。 新規の持続時間制約付きビタビアルゴリズムの動作を詳述する流れ図である。 ＥＣＧ信号を、Ｐ波領域、基線領域、ＱＲＳ群領域、Ｔ波領域、任意選択のＵ波領域及び第２の基線領域にセグメンテーションするための標準的な隠れマルコフモデルのアーキテクチャを示し、Ｐ波状態のための遷移確率ａ及び観測モデル確率分布ｂも合わせて示される図である。 図７ａのモデルと同じようにＥＣＧ信号をセグメンテーションするが、任意選択のＵ波領域を用いることなくセグメンテーションするための標準的な隠れマルコフモデルのアーキテクチャを示す図である。 ＥＣＧ信号をＰＲ間隔領域、ＱＲＳ群領域、Ｔ波領域及び基線領域にセグメンテーションするための標準的な隠れマルコフモデルのアーキテクチャを示す図である。 ＥＣＧ信号をＰＲ間隔領域、ＱＴ間隔領域及び基線領域にセグメンテーションするための標準的な隠れマルコフモデルのアーキテクチャを示す図である。 ＥＣＧ信号をＱＴ間隔領域、及び信号の残りの部分に対応する第２の領域（Ｘ）にセグメンテーションするための標準的な隠れマルコフモデルのアーキテクチャを示す図である。 持続時間制約の隠れマルコフモデルへの組み込みを図式的に示す図である。 ＥＣＧセグメンテーションのための組み込み式の持続時間制約を用いる隠れマルコフモデルのアーキテクチャの一部を示す図である。 ９９％信頼度下限とともに、多数の正常なＥＣＧ波形のための、ＥＣＧ波形長に対する信頼度指標のプロットを示す図である。 雑音（詳細には、「筋肉のアーテファクト(muscle artefact)」）によって汚染された２つのＥＣＧ波形を含み、雑音のある２つの波形（Ａ及びＢを付される）のセグメンテーションも示される、ＥＣＧ信号の一部を示す図である。 図１１に示される、雑音のある２つのＥＣＧ波形のための信頼度指標を示す図である。

Claims (8)

1. 確率的なセグメンテーションアルゴリズムによってセグメンテーションされた信号を解析するための、コンピュータによって実施される方法であって、
   複数のセグメンテーションされた信号特徴毎に信頼度指標を計算することと、
   個々の信号特徴長に対して前記信頼度指標をプロットすることと、
   密度モデリング技法を適用することであって、それによって、高い信頼度特徴に関連付けられるデータ空間の適切な領域を求める、適用することと、
   特定の信号特徴のための前記信頼度指標がこの領域の外側に該当するか否かを判定することと
   を含む方法。
2. 前記信号特徴は、より長い信号の一部を含む、請求項１に記載の方法。
3. 前記信頼度指標は、
   前記モデルを与えられるときの前記信号特徴の対数尤度と、
   前記モデルを与えられるときの前記信号特徴と最も確からしい状態シーケンスとの結合対数尤度と
   から成るグループから選択される指標である、請求項１又は２に記載の方法。
4. 前記信号は、心電図信号と、脳波図信号と、筋電図信号と、心臓電気生理現象信号と、心内電位図信号とから成るグループから選択される生物医学信号である、請求項１〜３のいずれか一項に記載の方法。
5. 前記信号の解析の結果についての情報を表示することをさらに含む、請求項１〜４のいずれか一項に記載の方法。
6. 前記信号は心電図信号であり、前記信号のセグメンテーションからＱＴ間隔についての情報が得られる、請求項１〜５のいずれか一項に記載の方法。
7. コンピュータシステム上で実行されることができるコンピュータプログラムであって、そのように実行されるときに、前記コンピュータシステムに請求項１〜６のいずれか一項に記載の方法を実行させることができるコンピュータプログラム。
8. 請求項７に記載のコンピュータプログラムを格納するコンピュータ読取り可能媒体。

Images