JP2016067490A - 末梢血管抵抗推定方法 - Google Patents

Abstract

【課題】推定精度を向上することができる末梢血管抵抗推定方法を提供する。【解決手段】電気回路モデルを用いて末梢血管抵抗を推定する末梢血管抵抗推定方法において、電気回路モデルの式の減衰定数α、角周波数ω、振幅係数Ｂ、時間係数ｔｂの値を変化させて前記式を用いた計算により得た橈骨動脈圧ｖｐの値と、測定により得た橈骨動脈圧ｖｐの値と、を比較し、前記計算により得た橈骨動脈圧の値ｖｐと前記測定により得た橈骨動脈圧の値ｖｐとの平均二乗誤差を最小化する減衰定数α、角周波数ω、振幅係数Ｂ、時間係数ｔｂの最適解を決定するにあたって、（ｉ）前記α、ω、Ｂ、ｔｂのそれぞれについて複数回の探索を行う、（ｉｉ）各回の探索では複数のサンプリング値を用いる、（ｉｉｉ）次回の探索では前回の探索よりも探索範囲を狭める、の処理を行う。【選択図】図６

Description

本発明は、生体の末梢血管抵抗を非侵襲にて推定する方法に関する。

末梢血管抵抗（ＴＰＲ）は動脈硬化及び高血圧の基本的な指標である。血液が細動脈を通り抜けることに抵抗するＴＰＲは、血管の硬さ及び直径に関連がある。血流に対するＴＰＲは、主として細動脈にある。その血管の直径がわずかに変化しても抵抗が著しく変わり、循環系全抵抗に大きく響く。

従来、このようなＴＰＲを非侵襲にて推定する方法及び装置が、例えば非特許文献１〜３に記載されている。これらの文献では、中枢側の血管抵抗、末梢側の血管抵抗、血液による慣性、及び、血管のコンプライアンス（粘弾性）の４要素を電気回路モデルに模してＴＰＲの推定値を算出する方法が記載されている。さらに、上述の４要素に加えて、血管のコンプライアンスを中枢側と末梢側に分けることで要素数を１つ増やし、５要素を電気回路モデルに模してＴＰＲの推定値を算出する方法が記載されている。

石山仁，天野和彦，笠原宏，上馬場和夫，許鳳浩: 橈骨動脈圧波形から構築した五要素電気回路モデルによる循環動態の非侵襲的評価法, 電子情報通信学会技術研究報告, Vol.97, No.215, MBE97-68, (1997-07). 石山仁，天野和彦，笠原宏，上馬場和夫: 五要素電気回路モデルによる大動脈圧波形の非侵襲的推定法, 電子情報通信学会技術研究報告, MBE98-55, (1998-07). 石山仁，笠原宏，児玉和夫，許鳳浩: 電気回路動脈系モデルを用いた橈骨動脈波形による循環動態パラメータの非観血的測定法, 電学論C, 113巻4号, 1993.

ところで、非特許文献１〜３に記載されているような、血管の要素を電気回路モデルに当てはめることでＴＰＲを推定する方法は、非侵襲にてＴＰＲを測定できるので、患者への負担も小さく、ＴＰＲを容易に得ることができるので大変便利である。

しかしながら、従来のこの種のＴＰＲ推定方法は、演算過程において各パラメータの最適解を求めることができなくなる場合が生じる欠点があった。また、推定精度の点で未だ不十分であった。

本発明は、以上の点を考慮してなされたものであり、血管の要素を電気回路モデルに当てはめることで末梢血管抵抗を推定する方法において、演算アルゴリズムを改良することで、より的確な演算を行って、推定精度を向上することができる末梢血管抵抗推定方法を提供する。

本発明の末梢血管抵抗推定方法の一つの態様は、
大動脈圧と、橈骨動脈圧と、動脈系中枢部での血管抵抗と、動脈系末梢部での血管抵抗と、血管のコンプライアンスと、血液の慣性と、を電気回路モデルに当てはめることで末梢血管抵抗を推定する末梢血管抵抗推定方法であって、
時間をｔ、大動脈圧の立ち上がりからその圧力が最低血圧値になるまでの時間をｔ_Ｐ１、橈骨動脈の圧力すなわち橈骨動脈波をｖ_ｐ、最低血圧値をＥ_ｍｉｎ、橈骨動脈の圧力における振幅係数をＢ、橈骨動脈の圧力における時間係数をｔ_ｂ、橈骨動脈の圧力における振動成分の第１の振幅係数をＤ_１、指数関数すなわちｅｘｐをｅ、減衰定数をα、角周波数をω、第１の位相をθ_１、１拍の時間をｔ_ｐ、橈骨動脈の圧力における振動成分の第２の振幅係数をＤ_２、第２の位相をθ_２としたときに、これらの血管の要素を下記の数式４及び数式５のように電気回路モデルによって表し、
数式４及び数式５における減衰定数α、角周波数ω、振幅係数Ｂ、時間係数ｔ_ｂの値を変化させて数式４及び数式５を用いた計算により得た橈骨動脈圧ｖ_ｐの値と、測定により得た橈骨動脈圧ｖ_ｐの値と、を比較し、前記計算により得た橈骨動脈圧の値ｖ_ｐと前記測定により得た橈骨動脈圧の値ｖ_ｐとの平均二乗誤差を最小化する減衰定数α、角周波数ω、振幅係数Ｂ、時間係数ｔ_ｂの最適解を決定するにあたって、次の（ｉ）、（ｉｉ）、（ｉｉｉ）の処理を含む。
（ｉ）前記α、ω、Ｂ、ｔ_ｂのそれぞれについて複数回の探索を行う。
（ｉｉ）各回の探索では、複数のサンプリング値を用いる。
（ｉｉｉ）次回の探索では前回の探索よりも探索範囲を狭める。

本発明によれば、より的確な演算を行うことができるので、末梢血管抵抗の推定精度を向上させることができるようになる。

実施の形態の末梢血管抵抗推定方法を実行するための測定システムの概略構成を示す図 デジタルデータに変換した血圧波形モデルを示す図 動脈系の四要素集中定数モデルに等価の、電気回路モデルを示す図 図４（ａ）は大動脈起始部の圧力波形を示す図であり、図４（ｂ）はそれを近似した三角波を示す図 若年健常者の平均波形を示す図であり、図５（ａ）は１０拍の脈波からの平均動脈圧波形を示す図、図５（ｂ）は計算波形と平均した実測波形とを示す図 減衰定数αに関する探索例を示す図 ４種のパラメータ（α，ω，Ｂ，ｔ_ｂ）の解の全てが各々の中域に含まれる場合の探索処理を示す図 ４種のパラメータ（α，ω，Ｂ，ｔ_ｂ）の解の１つ以上が各々の中域に含まれない場合の探索処理を示す図 減衰定数αに関する探索終了例を示す図 橈骨動脈圧波形の加算平均実測波形と、近似式から計算した波形との比較を示す図 数式３１の電気回路モデルを示す図 五要素集中定数動脈系電気回路モデルを示す図 正弦半波で近似した左心室圧波形を示す図 電気回路モデルにより推定算出されたＴＰＲと、観血的測定によるＴＰＲとの相関性を示す図 数式３６により推定算出されたＴＰＲと、観血的測定によるＴＰＲとの相関性を示す図

以下、本発明の実施の形態について、図面を参照して詳細に説明する。

＜１＞測定システム
まず、本実施の形態の末梢血管抵抗推定方法を実行するための測定システムについて説明する。図１は、その測定システムの概略構成を示す。また、図１の測定システムを用いて、本実施の形態の末梢血管抵抗推定方法の推定結果の検証も行った。

心電計は、ベッド上であお向けの状態の患者の両手首と両足首との４箇所に貼り付けた表面電極を介して心電図波形を得る。本実施の形態の場合、心電図波形を約１０［ｍｉｎ］間測定する。また同時に、右橈骨動脈脈波に非侵襲なトノメトリー法血圧計のセンサを設置し、橈骨動脈圧波形を非観血的連続的に採取する。

データレコーダを用いて得られたアナログ信号は、オシロスコープにより波形確認を行った。そして、その信号を１２［ｂｉｔ］Ａ／Ｄコンバータを用いてデジタルデータに変換する。この血圧波形モデルを、図２に示す。波形１周期において、血圧値が最も急激に変化するＤＢＰとＳＢＰとの時間間隔は約０．１［ｓ］である。この間隔内において５０点サンプリングできればほぼ元の波形に再現できると考えたため、また心電図波形においてＲ波の時間位置を精度良く決めるため、サンプリング周期はΔｔ＝０．００２［ｓ］とした。また、血圧（電圧）振幅値は１２［ｂｉｔ］で量子化した。そのデジタルデータを計算機に取り込んで、計算機によって演算及び解析を実行した。

＜２＞電気回路モデル解析による循環動態の推定
循環器系の状態を診断する場合、最も一般的に測定されるのが血圧や心拍数である。さらに詳しい診断には、血管の粘性抵抗やコンプライアンス(粘弾性)を測定することが必要である。血管の粘性抵抗やコンプライアンスを表す時、動脈系の振る舞いを記述する代表的なモデルの1つである四要素集中定数モデルがよく用いられる。

このモデルは、循環器系を電気的な等価回路とみなすものであり、動脈系中枢部における血液による慣性、中枢部における血液粘性による血管抵抗(粘性抵抗)、中枢部における血管のコンプライアンス(粘弾性)、及び末梢部における血管抵抗(粘性抵抗)の４つのパラメータから構成されている。

しかしながら、これらのパラメータを算出するには、大動脈起始部と切痕部の圧力波形及び血流量を測定する必要がある。その測定法には、動脈にカテーテルを挿入し直接測定する方法や、超音波などで間接的に測定する方法があるが、前者の方法は侵襲的な大がかりな方法となり、後者では圧力波形までは求められない。

近年、非侵襲的に橈骨動脈の脈波を検出する装置や、血圧と１回拍出量を収縮期面積法により同時に測定できる装置（脈波・コロトコフ音記録計）が開発された。これらの装置は、小型で測定法が非侵襲的である点に特徴がある。現在、手軽で非侵襲的にしかも安価な装置により、循環動態パラメータの評価ができるシステムを開発することを目的として上記の装置を使用し、橈骨動脈の脈波形と１回拍出量を計測するだけで、四要素集中定数モデルのパラメータ値を循環動態のパラメータとして近似的に計算する方法が見い出されている（非特許文献１〜３参照）。

＜２．１＞四要素電気回路モデル
＜２．１．１＞動脈系モデルと橈骨動脈波の近似式
血管中の血流Ｆ、血圧差Ｐ及び流れに対する抵抗Ｒ_ｂ、これら３つの量の関係は、次式に示すように、ちょうど電気回路における電流Ｉ、電圧（電位差）Ｅと電気抵抗Ｒ_ｅとの間におけるオームの法則と類同である。

１つの器官について言えば、その動脈と静脈における血圧の差をその器官内の血管の全抵抗で除すと、その器官の血流量が定まる。

動脈系の代表的なモデルの１つである四要素集中定数モデルは、動脈系を図３のような電気的な等価回路とみなすことにより仮定されるモデルであり、表１のように対応づけられる。なお、静電容量Ｃに対応づけたコンプライアンスは血管の軟度を表す量であり、粘弾性のことである。

ここで、大動脈起始部の圧力波形は一般に図４（ａ）のような波形であるから、圧力波形を図４（ｂ）の三角波で近似することにする。同図において近似波形の振幅と時間をＥ_ｏ、Ｅ_ｍ、ｔ_Ｐ、ｔ_Ｐ１とすると、任意の時間ｔにおける大動脈圧ｅは次式で表される。ただし、Ｅ_ｏは最低血圧、（Ｅ_ｏ＋Ｅ_ｍ）は最高血圧であり、ｔ_Ｐは１拍の時間、ｔ_Ｐ１は大動脈圧の立ち上がりからその圧力が最低血圧値になるまでの時間である。

数式２と数式３の大動脈圧ｅを、図３の等価回路に入力して、橈骨動脈の圧力（橈骨動脈波）ｖ_Ｐを求めると、次式となる。

数式４と数式５が、動脈系を四要素集中定数モデルで表した時の橈骨動脈波の近似式である。ここで、次式

とおけば、数式４と数式５の各定数α、ω、Ｂ、Ｄ_１、…などは次式のようになる。

ここで、ｖ_０１とｉ_０１はｔ＝０におけるｖ_Ｐとｉ_Ｃの初期値であり、ｖ_０２とｉ_０２はｔ＝ｔ_Ｐ１における初期値である。

＜２．１．２＞ 循環動態パラメータ
数式７の角周波数ωより、中枢部における血管抵抗Ｒ_Ｃは、次式のようになる。

Ｒ_Ｃが実数でかつ正となる条件は、次式の通りである。

ここで、一般にＲ_Ｐのオーダは１０^３［ｄｙｎ・ｓ／ｃｍ^５］程度、Ｃは１０^−４［ｃｍ^５／ｄｙｎ］程度であり、また、ωは脈波に重畳している振動成分の角周波数であるから１０［ｒａｄ／ｓ］以上であるとみてよい。よって、数式１９の上限はほぼ１／（ω^２Ｃ）とみなせることから、Ｌを簡略化のため近似的に、次式

とおくと、Ｒ_Ｃは次式となる。

また、数式２０と数式２１の関係により数式７の減衰定数αは、次式の通りとなる。

数式２０〜数式２２の関係を用いて、αとω及び４つのパラメータのいずれか１つ、例えば、血液による慣性Ｌで残りのパラメータを表すと、次式のように表すことができる。

上式より、モデルのパラメータはαとω及びＬが定まれば求めることができる。αとωは橈骨動脈波の実測波形から、また、Ｌは図３の等価回路から計算した１回拍出量（数式２６）を、測定した１回拍出量と一致するような値に設定することにより算出することができる。

図３の等価回路より、１回拍出量ＳＶは、次式となる。

そして、この式の単位を、ＭＫＳＡ単位系からＣＧＳ単位系に変換する。

＜２．１．３＞ ｉ_Ｃ初期値の算出法
循環動態パラメータの算出にあたって、まずｉ_Ｃの初期値を定める。健常者の実測より、手首の橈骨動脈は半径ｒ＝１．３［ｍｍ］、最高速血流量（収縮期）は ｖ_０１＝０．９０４［ｍ／ｓ］、最低速血流量（拡張期）は ｖ_０２＝０．２５９［ｍ／ｓ］であった。血流量は、ｆ_０ｊ＝πｒ^２ｖ_０ｊ （ｊ＝１，２）である。図３の回路において、末梢側のＣとＲ_Ｐの電圧が等しいことから、次式が成立する。

またｉは、次式のように表される。

数式２７及び数式２８よりｉ_Ｃは、次式により表される。

これらとｉ_０ｊ／Ｃにおいて、［ｄｙｎ／ｓ／ｃｍ^２］単位を［ｍｍＨｇ／ｓ］に変換する係数ａ＝１３３２［ｄｙｎ／ｃｍ^２／ｍｍＨｇ］を用いると、ｉ_０ｊ／Ｃは、次式により表される。

ここで、ｉ_０１／Ｃ及びｉ_０２／Ｃの算出段階においてはＣを決定できないため、非特許文献３などでも説明されている標準的な算出値Ｃを用いた。

＜２．１．４＞ 平均波形における特徴点の検出法
記録した橈骨動脈圧波形を１拍ごとに重ね合わせた。従来法では、１［ｍｉｎ］間における１拍当たりの平均波形を求めていた。しかし、被検者によっては心拍数が異なり、心拍数が多い場合は平均波形が滑らかになる可能性がある。また、標準的な１［ｍｉｎ］間の心拍数である７０拍の平均波形では、血圧変動が激しい場合は、平均波形があまりにも滑らかになり、各パラメータの算出が困難なケースが生じた。近年の手法では、１０拍の平均波形が利用されている点、また、この平均波形は元の波形と大きな相違がない点を考慮して、本実施の形態では、１０拍の加算平均化波形を用いることにした。

図５は若年健常者の平均波形を示すものであり、図５（ａ）は１０拍の脈波からの平均動脈圧波形を示し、図５（ｂ）は計算波形と平均した実測波形とを示す。若年健常者の平均波形（図５（ａ））において、最大値（第１ポイント）の時間ｔ_１と血圧値ｙ_１、凹部（第２ポイント）の時間ｔ_２と血圧値ｙ_２、及び２番目のピーク（第３ポイント）の時間ｔ_３と血圧値ｙ_３を検出した。また、１拍の時間ｔ_Ｐ、最低血圧値Ｅ_ｍｉｎ（数式４と数式５の第１項目）を読み取った。抽出した値、及び１回拍出量ＳＶの算出例を表２に示す。

１回拍出量ＳＶは、一般的に循環動態には影響が小さく、また簡易な測定法の開発という観点から、被検者に依存しない標準的な固定値とした。なお、凹部と２番目のピークがはっきりしないなだらかな脈波の場合には、第２と第３ポイントの時間をｔ_２＝２ｔ_１、ｔ_３＝３ｔ_１としてその点の血圧値を読み取った。

＜２．１．５＞ α，ω，Ｂ，ｔ_ｂの算出アルゴリズム
読み取ったデータを基にして、数式４のα，ω，Ｂ及びｔ_ｂを決定する。α，ω，Ｂ及びｔ_ｂを適当に変化させ、ｔ＝ｔ_１，ｔ_２，ｔ_３における橈骨動脈波の圧力ｖ_Ｐ１，ｖ_Ｐ２，ｖ_Ｐ３を計算し、（ｙ_１−ｖ_Ｐ１），（ｙ_２−ｖ_Ｐ２），（ｙ_３−ｖ_Ｐ３）が最小となるようなα，ω，Ｂ及びｔ_ｂを算出する。従来法では、第１ポイント〜第３ポイントの３点だけで実測の平均波形とモデルによる計算波形とを比較し、その平均二乗誤差を最小化していたが、次の２つの問題点があった。

１． ３点だけの比較では、平均波形と計算波形との一致が困難であり、大きな誤差が生じてしまう。
２． 高齢高血圧患者の場合では反射波が存在しないケースが多く、第２及び第３ポイントを的確に判定することができない。

そこで、本実施の形態では、計算量を考慮して Δｔ＝０．０１［ｓ］ごとに平均波形と計算波形を比較し、その平均二乗誤差を最小化する手法を提案する。

また、α，ω，Ｂ及びｔ_ｂの従来の算出法では、４種のパラメータにおいて各々の指定範囲を設定し、その領域内における全ての場合について計算し、最適解を求めていた。この手法では、計算量が膨大であり、かつ設定した指定範囲に依存するため、循環器系疾患や刺激後の動脈圧波形の場合は、最適解を得ることができるとは限らない。

そこで、本実施の形態では、効率が良く、かつ安静状態の健常者以外の場合でも利用できる可変領域を用いるアルゴリズムを提案する。

本実施の形態による、α，ω，Ｂ，ｔ_ｂの算出アルゴリズムは以下の通りである。
１） まず、α，ω，Ｂ，ｔ_ｂの初期設定探索領域、最終探索サンプリング間隔Δα，Δω，ΔＢ，Δｔ_ｂ、及び領域縮小率を表３のように定めた。表３に示す各パラメータの初期条件は、被検者データの分析、生理学的に適切な範囲内、算出誤差の低減、計算時間量の抑制などを考慮して、経験則に基づいて決定した。なお、初期設定探索領域は８回縮小した時に、探索サンプリング間隔が最終値に到達するように設定した。図６は、αに関する探索例を示す図である。

２） １ブロックの探索において、１種のパラメータの比較可能ポイント数は７点とした。

３） ７^４通りの各々の組（α，ω，Ｂ，ｔ_ｂ）について、計算波形と平均波形との平均二乗誤差を求め、平均二乗誤差が最小となる解を算出する。

４） ４種のパラメータにおいて１ブロックの探索が完了した時、全てのパラメータにおいて各々の設定探索領域を左域・中域・右域と３等分に分割する。

５） その時に算出した解の全てが各々の中域に含まれる場合は、最適解の局所領域に順調に探索していると判断し、領域縮小率に従って設定探索領域を縮小し、次のブロックの探索をする。また、その時に算出した解の１つ以上が各々の中域に含まれない場合は、最適解の局所領域に到達していないと判断し、設定探索領域は算出した解を中心とする領域にそれぞれ移動させ、次の探索をする。図７は４種のパラメータ（α，ω，Ｂ，ｔ_ｂ）の解の全てが各々の中域に含まれる場合の探索処理を示す図であり、図８は４種のパラメータ（α，ω，Ｂ，ｔ_ｂ）の解の１つ以上が各々の中域に含まれない場合の探索処理を示す図である。

６） 前回の探索よりも次回の探索のサンプリング間隔を小さくするものとし、次回の探索のサンプリング間隔が既定の最終値に到達し、かつ次回の探索で求めた解の組が前回の探索で求めた解の組と一致した場合に、その解を最適解として、探索を終了する。図９は、αに関する探索終了例を示す図である。

このようなアルゴリズムを利用して得られた各パラメータの算出例を下記に示す。
α＝４．７１［１／ｓ］，ω＝１８．３２８［ｒａｄ／ｓ］，Ｂ＝３４．８７［ｍｍＨｇ］，ｔ_ｂ＝０．９７４［ｓ］

ｔ_Ｐ１，Ｅ_ｍ，Ｅ_ｏは数式８、数式９、数式２３〜数式２５から求めることができる。その算出例は、次の通りである。
ｔ_Ｐ１＝０．９４８［ｓ］，Ｅ_ｍ＝３６．１７［ｍｍＨｇ］，Ｅ_ｏ＝５３．５３［ｍｍＨｇ］

高齢高血圧患者は若年健常者に比べて、最低血圧Ｄとその直後の最高血圧ＳとのＤ−Ｓ時間間隔ｔ_１が長いこと、そして、血圧単位のパラメータＢ及びＥ_ｍが大きい傾向が示された。ｔ_１が長くなる原因としては、血管弾性力（瞬発力）の低下が考えられる。またＢ及びＥ_ｍが大きい原因としては、抵抗Ｒ_Ｐが増大しているためである。

＜２．１．６＞ 循環動態パラメータの算出例
次に、数式２６から計算する１回拍出量が、固定値と一致するような血液による慣性Ｌの値を算出した後、残りのパラメータ値を数式２２により求める。そのパラメータの算出例は、以下の通りである。

Ｃ＝２．４０７×１０^−４［ｃｍ^５／ｄｙｎ］，Ｌ＝１２．３６８［ｄｙｎ・ｓ^２／ｃｍ^５］，Ｒ_Ｃ＝５８．２９［ｄｙｎ・ｓ／ｃｍ^５］，Ｒ_Ｐ＝８８１．５０［ｄｙｎ・ｓ／ｃｍ^５］

また、直流的な総末梢血管抵抗（中枢部と末梢部の血管抵抗の和）ＴＰＲは、
ＴＰＲ ＝ Ｒ_Ｃ＋Ｒ_Ｐ ＝ ９３９．７９［ｄｙｎ・ｓ／ｃｍ^５］
である。確認のため、算出したパラメータで数式１９を計算すると、１２．１６７≦Ｌ≦１２．３６８ となり、数式２０の近似は妥当であるといえる。

図５（ｂ）に、算出したパラメータを用いて計算した橈骨動脈波形（実線）と１０拍分を平均した実測波形（破線）との比較を示す。この例では、計算した血圧波形は、実測波形とほぼ類似しているといえる。

３種類の特徴点だけを比較する従来法では、計算波形と実測波形との平均二乗誤差は、１０．４５９［ｍｍＨｇ^２］であった。これに対し、本実施の形態によるΔｔ＝０．０１［ｓ］ごとに計算波形と実測波形を比較し、その平均二乗誤差を最小化する方法による平均二乗誤差の算出結果は、２．６６４［ｍｍＨｇ^２］であり、大幅に誤差を低減させることができた。

＜２．２＞ 五要素電気回路モデル
＜２．２．１＞ 電気回路モデルの改善
前節では、動脈系の振る舞いを表す代表的なモデルの１つである動脈中枢側の血管抵抗、末梢側の血管抵抗、血液による慣性、血管のコンプライアンス（粘弾性）の４つの要素で構成された四要素集中定数モデルについて述べた。四要素集中定数モデルを電気回路に置き換え、入力電圧を三角波で近似した大動脈圧波形に、ＣＲ並列回路部の電圧を橈骨動脈圧波形に、平均電流を１［ｓ］間当たりの心拍出量に対応させて４つのパラメータを非侵襲的に算出する方法を示した。

ところで、血管のコンプライアンスは、中枢側に向かうほど血管壁の組成がより弾性的になって大きくなるため、中枢側のコンプライアンスは、動脈系全体のコンプライアンスの中で大きな割合を占めている。四要素集中定数モデルでは、血管抵抗については動脈系の中枢側と末梢側に分離して示しているため両者を評価できるが、コンプライアンスについては血管抵抗のように分離されていないため、中枢側と末梢側に分離して評価することができない。

本実施の形態では、まず、コンプライアンスを中枢側と末梢側に分離して配置した動脈系電気回路モデルの回路構成を考える手がかりとして、橈骨動脈圧波形の近似式を実測波形より導出した。その近似式をもとに回路構成を考察した結果、中枢側コンプライアンスに相当する静電容量を電気回路モデルの入力端子に並列に、末梢側コンプライアンスに相当する静電容量を電気回路モデルの出力端子に並列に配置した五要素集中定数動脈系電気回路モデルを構築した。

＜２．２．２＞ 橈骨動脈圧波形の近似式とその電気回路モデル
図１０は、橈骨動脈圧波形の加算平均実測波形と、近似式から計算した波形との比較を示す図である。図１０の点線は、実測した１０拍の加算平均橈骨動脈圧波形である。同図より、橈骨動脈圧波形の前半では減衰を伴う振動成分が含まれているが、後半では振動成分がほとんど認められない単調で緩慢な減衰を示していることが分かる。そこで、橈骨動脈圧波形の近似式を次式のように表すことにする。

波形の後半において、振動成分がほとんど認められない単調な減衰を示すためには、数式３１の減衰定数α_５，βと角周波数ω_５は、それぞれα_５＜＜β，α_５＜＜（ω_５／２π）でなければならない。この条件下で数式３１の各定数を適当に決定して計算した波形が図１０の実線であり、実測波形とほぼ一致した波形が得られたため、数式３１は橈骨動脈圧波形の近似式として妥当であると考えられる。この時の各定数は次の通りである。
α_５＝０．５４８［ｓ^−１］，β＝５．１５９［ｓ^−１］，ω_５＝１７．２７［ｒａｄ／ｓ］，Ｋ_１＝−４０．７４［ｍｍＨｇ］，Ｋ_２＝８７．０６［ｍｍＨｇ］

数式３１に基づいて電気回路モデルを考えると、まず、数式３１の第１項は、図１１（ａ）の抵抗Ｒ_１（＝Ｒ_Ｃ）、インダクタンスＬ、静電容量Ｃ_２の直列回路で表すことができる。第２項は、図１１（ｂ）の静電容量Ｃ_１と抵抗Ｒ_２（＝Ｒ_Ｐ）の放電回路で表すことができる。また、数式３１が繰り返し持続するためには電源が必要である。動脈系で電源に相当するものは心臓であるから、本電気回路モデルでは左心室の圧力を電源とした。

電源を考慮して、図１１（ａ）、図１１（ｂ）を重ね合わせて構成した回路が図１２の五要素集中定数動脈系電気回路モデルである。Ｃ_２Ｒ_２並列回路部の端子電圧を橈骨動脈圧波形に対応づけると、Ｃ_２とＲ_２は、末梢側におけるコンプライアンスと血管抵抗に相当し、電源側に並列に配置したＣ_１と電源と末梢側の間に直列に配置したＲ_１は、中枢側におけるコンプライアンスと血管抵抗に相当する。血管抵抗は末梢側にいくほど血管径が細くなるため、中枢側よりも末梢側の方が大きい。逆に、血管のコンプライアンスは中枢側にいくほど血管壁がより弾性的になるため、末梢側よりも中枢側の方が大きい。故に図１２において、Ｒ_１＜＜Ｒ_２，Ｃ_１＞＞Ｃ_２の関係が成り立つ。

図１０の点線の前半部分は収縮期に当たり、波形の変化が大きい。この区間では図１２のダイオードＤは導通状態であり、Ｃ_１は電源電圧で充電される。また、末梢側に流れる電流は変化が大きいから、主に電源→Ｒ_１→Ｌ→Ｃ_２の経路で流れると考えられる。図１０の点線の後半部分は拡張期に当たり、波形の変化は緩慢である。この区間ではダイオードＤは阻止状態であり、Ｃ_１は収縮期に充電した電荷を末梢側に放電するが変化が緩慢であるから、Ｌの効果は無視することができ、さらにＲ_１＜＜Ｒ_２，Ｃ_１＞＞Ｃ_２であるから、ほぼＣ_１Ｒ_２の時定数で放電すると考えられる。

以上のことから、図１２は数式３１の電気回路モデルとして妥当であると考えられる。ここで、図１２の電気的要素と動脈系の要素との対応づけをまとめると表４のようになる。

左心室圧波形を容易に、しかも安価な装置で非侵襲的に測定することは困難である。左心室圧波形は一般に図４（ａ）に示すような形状であるから、左心室圧波形を図１３に示す正弦半波で近似することにした。ここで、図１３は、正弦半波で近似した左心室圧波形を示すものであり、Ｅ_ｍは左心室の最高血圧、ｔ_ｓは左心室に圧力が印加されている間の時間（以下、加圧時間と呼ぶ）、ｔ_１は収縮期開始時間、ｔ_２は収縮期終了時間、ｔ_ｐは１拍の時間である。図１２より収縮期（ｔ_１≦ｔ≦ｔ_２）では、（左心室圧）≧（大動脈圧）であるからダイオードＤは導通状態であり、次式が成立する。

ここで、ω_ｓ＝π／ｔ_ｓである。拡張期（ｔ_２≦ｔ≦ｔ_ｐ＋ｔ_１）では、（左心室圧）≦（大動脈圧）であるからダイオードＤは阻止状態であり、次式が成立する。

数式３２、数式３３を解くと、収縮期と拡張期における橈骨動脈圧ｖ_２は、次式となる。

ここで、数式３４、数式３５の角周波数ω_１，ω_２、減衰定数α，γ、係数Ｂ_０，Ｂ_１，Ｂ_２，Ｂ_ｍ、位相角θ_ｓ，θ_１，θ_２はＲ_１，Ｒ_２，Ｌ，Ｃ_１，Ｃ_２，Ｅ_ｍで決定される定数であり、ｔ’＝ｔ−ｔ_１，ｔ”＝ｔ−ｔ_２である。

＜２．２．３＞ 循環動態パラメータの算出例
本実施の形態による、電気回路モデルにおけるパラメータを取得するための算出法を以下に示す。

α_５，β，ω_５，Ｋ_１，Ｋ_２の初期設定探索領域、最終探索サンプリング間隔Δα_５，Δβ，Δω_５，ΔＫ_１，ΔＫ_２及び領域縮小率を表５に示す。これらの初期条件と数式３１を用いることにより、計算波形を求めた。なお、その算出にあたっては前章で述べたα，ω，Ｂ，ｔ_ｂの算出アルゴリズムと同様な手法を採用した。

本実施の形態の電気回路モデルで決定しなければならないパラメータは、Ｒ_１，Ｒ_２，Ｌ，Ｃ_１，Ｃ_２，Ｅ_ｍ，ｔ_ｓの７つである。これを一度に決定するのは困難であるから、Ｒ_１，Ｒ_２，Ｌ，Ｃ_２，Ｅ_ｍについては前節で提案した算出方法で求めた。

Ｃ_１とｔ_ｓの算出にあたっては、従来、まず、ｔ_ｓを（ＱＴ＋０．１）［ｓ］〜（ＱＴ＋０．２）［ｓ］の範囲で変化させて、各ｔ_ｓに対して計算した橈骨動脈圧の最低血圧値が測定値と一致するようなＣ_１を求めていた。ここで、ＱＴは心電図におけるＱＴ間隔（電気的収縮時間）である。次に、各ｔ_ｓとＣ_１の組の中で、橈骨動脈圧の計算波形と加算平均波形の平均誤差が最小となるｔ_ｓとＣ_１を本実施の形態の電気回路モデルのパラメータとして採用していた。しかしこの手法では、ｔ_ｓとＣ_１との理論的関係が曖昧であり、解が通常の範囲を越える場合があった。

そこで、本実施の形態では、数式３３から一脈拍内におけるＣ_１の経時的変動を調べた。その結果、波形の前半ではＣ_１は激しく振動し、波形の後半ではＣ_１はほぼ中間に収束する傾向があることを示した。よって、その収束値をＣ_１の最適解とすることにした。またｔ_ｓは求めないこととした。

以下に算出した循環動態パラメータの例を示す。
Ｃ_１＝１．６９４×１０^−３［ｃｍ^５／ｄｙｎ］，Ｃ_２＝２．４０７×１０^−４［ｃｍ^５／ｄｙｎ］，Ｌ＝１２．３６８［ｄｙｎ・ｓ^２／ｃｍ^５］，Ｒ_１＝５８．２８９［ｄｙｎ・ｓ／ｃｍ^５］，Ｒ_２＝８８１．４９８［ｄｙｎ・ｓ／ｃｍ^５］

図１０は、橈骨動脈圧波形の加算平均実測波形（点線）と、算出したパラメータを用いて計算した波形（実線）とを比較したものである。両波形は、凹凸部を除いてほぼ類似している。Ｃ_１の算出手法を改良した五要素モデルの平均誤差は、１．０４１［ｍｍＨｇ］であった。この数値は、前節で改良した四要素モデルの平均誤差（１．６１０［ｍｍＨｇ］）及び、従来法の五要素モデルの平均誤差（２．４９０［ｍｍＨｇ］）に比べて、大幅に改善されている。

このようにして、血管のコンプライアンスを中枢側と末梢側に分離して配置した動脈系電気回路モデルの回路構成を橈骨動脈圧波形の近似式に基づいて考察した結果、中枢側コンプライアンスに相当する静電容量を電気回路モデルの入力端子に並列に、末梢側コンプライアンスに相当する静電容量を電気回路モデルの出力端子に並列に配置した五要素集中定数動脈系電気回路モデルを構築することができた。

＜２．３＞ 安静状態検査における循環動態評価
＜２．３．１＞ 解析結果
まず、若年者と高齢者の各グループにおいて、前処理パラメータの平均値を算出した。その結果を表６に示す。

これらの表における比較から、高齢者では以下のような特徴が考えられる。
１．数式４及び数式３１におけるα，βが低いことから、ＳＢＰ直後における減衰が緩やかである。
２．α_５が高いことから、波形の後半における減衰が急である。
３．ω，ω_５が低いことから、圧の振動性が小さく、弾力性に乏しい。
４．Ｂ，Ｅ_ｍ，Ｅ_０が高いことから、血圧が高い。
５．数式３１におけるＫ_１が低くＫ_２が高いことから、圧の振動性の影響がわずかに作用するだけであり弾力性に乏しい。

次に、定常状態において、若年者と高齢者からそれぞれの循環動態パラメータの平均値を算出した。その結果を表７に示す。

高齢高血圧患者は若年健常者に比べて、コンプライアンスＣ_１，Ｃ_２が小さい、インダクタンスＬが大きい、中枢部血管抵抗Ｒ_１（Ｒ_Ｃ）・末梢部血管抵抗Ｒ_２（Ｒ_Ｐ）・総末梢血管抵抗ＴＰＲが大きい傾向にあることがいえる。加齢に伴い血圧が上昇して動脈が硬くなると、血管弾性力Ｃが低減し、血液による慣性Ｌ及び血管抵抗Ｒが増大すると考えられる。従って、この電気回路動脈系モデルから得られる循環動態パラメータは、動脈硬化度の有効な定量的評価法につながる可能性があるといえる。また血圧波形測定は、非観血的かつ約１［ｍｉｎ］間だけの短時間測定であるため、患者に負担を強いることがない。ただ、問題点としては、１回拍出量を測定しないため、１回拍出量の測定値ＳＶを一定値にしたことである。

＜２．３．２＞ 観血的ＴＰＲ測定検査における循環動態評価
本実施の形態による電気回路モデルにより推定算出されたＴＰＲの正確性を検証するため、観血的測定によるＴＰＲとの相関性を調べた。その結果を図１４に示す。モデル推定ＴＰＲは、観血的測定ＴＰＲより小さく推定する傾向があり、相関度が低いといえる（ｒ＝０．５７， ｐ＜０．０１）。モデル推定ＴＰＲの推定誤差は、次のような原因が考えられる。

１．血圧波形の測定時には同時に観血的カテーテル検査を施行していることから、波形に雑音が重畳されることが多い。
２．モデルに入力する血圧波形は１０拍分の平均であることから、１対１の対応ではない。
３．このモデルは健常者を対象としたモデルであることから、血圧波形の形状に異常がある患者の場合には、循環動態に関わる理論式に適合するとは限らない。
４．図１４より、このモデルにＴＰＲの上限（飽和状態）があることから、ＴＰＲが極めて高い患者の場合でも、モデルＴＰＲが２０００［ｄｙｎ・ｓ／ｃｍ^５］を超えることは困難である。

このような観点から、モデルＴＰＲの絶対値は的確であるとはいえない。

モデルＴＰＲを観血的ＴＰＲに近づけるためには、モデルに入力する血圧波形においてＴＰＲが高いと判断される異常な局所領域またはＲ_１，Ｒ_２を検知し、その情報を何らかの基準から数値化してモデル理論式に付加する必要がある。そこで、本実施の形態では、既存の各種データを分析し、血圧波形の異常がＲ_１（動脈系中枢部での血管抵抗），Ｒ_２（動脈系末梢部での血管抵抗）に反映することを考慮し、モデルＴＰＲ（総末梢血管抵抗推定値）の推定式を、次式とすることを提案する。

ここで、数式３６の係数ｗ_１，ｗ_２は表８に従う。

電気回路モデルの数式３６により推定算出されたＴＰＲの正確性を検証するため、観血的測定によるＴＰＲとの相関性を調べた。その結果を図１５に示す。モデル推定ＴＰＲは、観血的測定ＴＰＲと良好な直線的正相関を示す（ｒ＝０．７７， ｐ＜０．００１）。

＜３＞． 結論
上述の実施の形態では、心血行動態の各指標を分析する手法を提案し、循環動態指標を算出する電気回路モデル解析について述べた。

橈骨動脈の脈波形を計測するだけで、四要素集中定数モデルのパラメータ値を循環動態の指標として近似的に算出する方法について調べた。次に、橈骨動脈圧波形の近似式から、血管のコンプライアンスを動脈系の中枢側と末梢側に分離して配置した五要素集中定数動脈系電気回路モデルを、先に提案された方式と組み合わせて構築した。

このモデルの妥当性を検証するため、応用例として血管拡張剤投与検査等を行った。投与前後における橈骨動脈圧波形の測定及び記録は、非侵襲的な装置を用いた。次に、入力電圧を正弦半波で近似した左心室圧波形に、出力電圧を橈骨動脈圧波形に、平均電流を１［ｓ］間当たりの心拍出量に対応させ、ＴＰＲ誤差及び計算量の削減に重点を置いて、モデルに対するパラメータを算出した。その結果、橈骨動脈の平均血圧と脈圧及び１回拍出量の経時的変化をモデルにおける５種の循環動態パラメータの変化として、動脈系を中枢側と末梢側に分離して、定性的かつ定量的な評価に利用できることが確認された。従って本方式は、血管抵抗、コンプライアンス及び血液の慣性を明確に評価でき、特に血管抵抗を中枢部と末梢部に分離して評価できるという利点を最大限に発揮した。

従来の装置においてもＴＰＲや動脈硬化度等を算出することができるが、本法は非侵襲的であり、苦痛が少なく容易に繰り返し測定可能な橈骨動脈圧波形だけで算出できるため、肉体的、精神的にも負担をかけない状態で、長時間にわたる循環動態の評価を可能にするものである。

この簡易な方法は、
（１）被検者に負担を与えない非侵襲的方法であること、
（２）一心拍ごとの値が得られる連続的方法であること、
等の特徴を有し、今後、血行動態力学に関する臨床研究の場で広範に応用可能であるといえる。

上述の実施の形態は、本発明を実施するにあたっての具体化の一例を示したものに過ぎず、これらによって本発明の技術的範囲が限定的に解釈されてはならないものである。すなわち、本発明はその要旨、またはその主要な特徴から逸脱することなく、様々な形で実施することができる。

本発明は、生体の末梢血管抵抗を非侵襲にて推定する場合に広く適用し得る。

Claims (5)

1. 大動脈圧と、橈骨動脈圧と、動脈系中枢部での血管抵抗と、動脈系末梢部での血管抵抗と、血管のコンプライアンスと、血液の慣性と、を電気回路モデルに当てはめることで末梢血管抵抗を推定する末梢血管抵抗推定方法であって、
   時間をｔ、大動脈圧の立ち上がりからその圧力が最低血圧値になるまでの時間をｔ_Ｐ１、橈骨動脈の圧力すなわち橈骨動脈波をｖ_ｐ、最低血圧値をＥ_ｍｉｎ、橈骨動脈の圧力における振幅係数をＢ、橈骨動脈の圧力における時間係数をｔ_ｂ、橈骨動脈の圧力における振動成分の第１の振幅係数をＤ_１、指数関数すなわちｅｘｐをｅ、減衰定数をα、角周波数をω、第１の位相をθ_１、１拍の時間をｔ_ｐ、橈骨動脈の圧力における振動成分の第２の振幅係数をＤ_２、第２の位相をθ_２としたときに、これらの血管の要素を次の数式１のように電気回路モデルによって表し、
   

   

   前記数式１における減衰定数α、角周波数ω、振幅係数Ｂ、時間係数ｔ_ｂの値を変化させて前記数式１を用いた計算により得た橈骨動脈圧ｖ_ｐの値と、測定により得た橈骨動脈圧ｖ_ｐの値と、を比較し、前記計算により得た橈骨動脈圧の値ｖ_ｐと前記測定により得た橈骨動脈圧の値ｖ_ｐとの平均二乗誤差を最小化する減衰定数α、角周波数ω、振幅係数Ｂ、時間係数ｔ_ｂの最適解を決定するにあたって、次の（ｉ）、（ｉｉ）、（ｉｉｉ）の処理を含む、
   （ｉ）前記α、ω、Ｂ、ｔ_ｂのそれぞれについて複数回の探索を行う、
   （ｉｉ）各回の探索では、複数のサンプリング値を用いる、
   （ｉｉｉ）次回の探索では前回の探索よりも探索範囲を狭める、
   方法。
2. さらに、次の（ｉｖ）、（ｖ）の処理を含む、
   （ｉｖ）前記複数のサンプリング値を、値が小さいものから大きいものへと少なくとも左域、中域、右域の３つの領域に領域分けし、
   （ｖ）ある回の探索において、前記α、ω、Ｂ、ｔ_ｂの全ての解が中域に含まれる場合には、次回の探索をその中域を中心に狭めた探索範囲で行う一方、ある回の探索において、前記α、ω、Ｂ、ｔ_ｂの解の１つ以上が中域に含まれない場合には、その解を中心とする領域に探索範囲を移動させて再度解の探索を行う、
   請求項１に記載の方法。
3. さらに、次の（ｖｉ）の処理を含む、
   （ｖｉ）前回の探索よりも次回の探索のサンプリング間隔を小さくするものとし、前記次回の探索のサンプリング間隔が既定の最終値に到達し、かつ前記次回の探索で求めた解の組が前記前回の探索で求めた解の組と一致した場合に、その解を最適解として、探索を終了する、
   請求項２に記載の方法。
4. 前記最適解を用いて推定した動脈系中枢部での血管抵抗をＲ_１、前記最適解を用いて推定した動脈系末梢部での血管抵抗をＲ_２、総末梢血管抵抗推定値をＴＰＲとしたとき、前記総末梢血管抵抗推定値ＴＰＲを、係数ｗ_１，ｗ_２を用いて、次の数式２により求める、
   

   

   請求項１から請求項３のいずれか一項に記載の方法。
5. 前記係数ｗ_１，ｗ_２を、前記（Ｒ_１＋Ｒ_２）の値の範囲に応じて変える、
   請求項４に記載の方法。
   

Images