計測自動制御学会東北支部 40 周年記念学術講演会
2004 年 12 月 22 日，23 日 日本大学工学部

講演番号：A07

補助人工心臓における差圧・流量推定○花岡 哲文，田中 明 † ，吉澤 誠 ∗ ，ポール オレガリオ，小川 大祐，阿部 健一，
白石 泰之 ‡ ，山家 智之 ‡ ，仁田 新一 ‡東北大学大学院工学研究科, † 福島大学共生システム理工学類, ∗ 東北大学情報シナジーセンター,‡
東北大学加齢医学研究所

Estimation of Pressure Difference and Flow Rate for a Ventricular Assist Device*Tetsubumi Hanaoka, Akira Tanaka† , Makoto Yoshizawa∗ , Paul Olegario, Daisuke Ogawa,
Kenichi Abe, Yasuyuki Shiraishi‡ , Tomoyuki Yambe‡ , Shinichi Nitta‡
Graduate School of Engineering, Tohoku University,†College of Symbiotic Systems Science, Fukushima University,∗
Information Synagy Center, Tohoku University,‡Aging and Cancer, Tohoku University
Abstract– The authors have previously proposed a new method of estimating pressure difference and flow rate
without direct measurement. This method uses an ARX model from supplied current and rotational speed, whichrepresents the dynamic state where rotational speed changes with time. This study indicated that it is possible toreduce the number of input parameters without reducing estimation accuracy.
Key Words: ventricular assist device, pressure head and flow estimation, ARX model1

はじめに補助人工心臓 (ventricular assist device; VAD) は，急性・慢性心不全に陥った心臓に対して，血液循環のためのポンプ機能を補う治療用装置である．症例の増加，埋め込み期間の長期化，機械的信頼性の高さなどが考慮され，最近，VAD に対して定常流型ポンプを採用する例が増えている．米国では，MicroMed Technology 社の DeBakey VAD が FDA の承認を受け [1]，日本でも，TERUMO 社の DuraHeart が臨床試験に供されるようになった [2]．
人工心臓の制御や監視を行う際，血圧や血流量の計測が必要となる．しかし，血圧センサや血流計の使用は，• 血液と接触する部分での血栓形成• ドリフト誤差による精度低下• 装着の侵襲性を行っている．一方で，著者らの時系列法 [7, 8] は，7
入力の ARX モデルを用い，動特性を考慮した推定を行うことを目的としており，血液の粘性変化の補償が可能である．
本研究では，時系列法による推定法，Funakubo [3] らによる定常特性から構築されたモデルによる推定法との比較を行うとともに，時系列法の前提となる静特性モデルの導出法を改めることで，推定モデルに必要な入力パラメータ数を減らすことを検討した．

2 方法
2.1

従来法

静特性試験の結果から，ポンプの差圧または流量 yは，b j ( j = 1 ∼ 6) を定数として以下のように表すことができる [7]．

y = b1 N 2 I + b2 NI + b3 I + b4 N 2 + b5 N + b6• 電磁ノイズ発生 (電磁流量計)• 再キャリブレーションが困難(1)

この式に基づき，動特性を表現するために，(1) 式を拡張した ARX モデル (2) を考える．のような問題があるため，それらの直接的使用はなるべく避ける必要がある．
一方で，定常流型人工心臓の血液ポンプを駆動するモータの回転数や供給電流は，差圧 (揚程) やポンプ拍出量と高い相関をもつ．このことを利用し，供給電流，回転数からポンプの差圧または流量を推定する手法が報告されている [3–6]．これらの手法は，回転数の時間変化を考慮しない定常状態において差圧・流量の推定69y(k) = −L
Xi=1

ai y(k − i) +Mj
7 X
X
j=1 i=0bi j u j (k − i) + w(k) (2)

ここで，k は t = k∆t (k = 1, 2, 3, . . . ) を満たす離散時間(∆t[s] はサンプリング周期を表す)，w(k) は残差，y(k)は離散時間 k における流量 Q または差圧 P，M j ，L は入力，出力項の次数である．u j (k) は ARX モデルの入力であり，以下に示されるような 7 種の入力である．

u1 (k) = N 2 (k) · I(k)また，差圧 P についても同様の操作を行い，以下の推定式を得た．

(3)

0(4)

u3 (k) = I(k)(5)

2u4 (k) = N (k)(6)

ここで，a ，b ，α0 ，β0 は定数である．

u5 (k) = N(k)

(7)

2.3

u6 (k) = 1(8)u7 (k) = K(k)(9)

2.1 節の考え方を用い，(15) 式を拡張し，(17) 式で表される ARX モデルを考える．

P̂ =0

ここで，u7 (k) = K(k) は，(1) 式には含まれないが，係数 ai ，bi j を同定するときと推定システム運用時での環境変化を補償するために追加されたものである [8]．ここでは，(10) 式のように，過去 n 個の電流 I と回転数N の比として計算している．ここで，使用するデータ数 n は 200 としている．

00

N(k) + N(k − 1) + · · · + N(k − n + 1)=
I(k) + I(k − 1) + · · · + I(k − n + 1)

Q(k) = −ŷ(k) = −
2.2

i=1

ai ŷ(k − i) +j=1 i=0

bi j u j (k − i)LX
i=1

ai Q(k − i) +(10)Q̂(k) = −(12)α，β は回転数 N の関数として，以下のように表した．
(13)(14)β(N) = β0 N

ここで，a，b，α0 ，β0 は実験的に求められた定数である．
(13)，(14) 式は α(N) と β(N) を N の指数関数で表しているのに対し，(1) 式では Q = γ(N)I + δ(N) とおいて，γ(N) と δ(N) を N の 2 次関数で近似している点が異なる．(13)，(14) 式を (12) 式に代入することにより，推定流量 Q̂ を求める式が (15) 式のように得られる．

Q̂ =β01
IN −a − N b−aα0α0

bi j u j (k − i) + w(k)(17)

u1 (k) = N −a (k) · I(k)(18)

u2 (k) = N(19)b−a(k)(20)

オフラインで未定係数 ai ，bi j を算出し，以下の式で推定値 Q̂ を求める．

(11) 式では，推定に必要な入力数が 7 個必要である．
入力パラメータ数の少ないモデルを得るために，Funakubo [3] の手法を参照し，静特性モデルの構築を行った．流量 Q と電流 I の関係は，傾きを α，切片を β として，以下のように近似することができる．

bj=1 i=0

ここで，(17) 式で新たに追加された入力 u3 (k) = K(k)は，(10) 式と同一である．

(11)α(N) = α0 N aMj3 X
X

u3 (k) = K(k)

入力数を減らした静特性モデルの構築

I = α(N)Q + β(N)0u j (k) は ARX モデルの入力で，以下に示されるような3 種の入力とする．

推定システムの運用の際には，y(k)，u j (k) の時系列データをあらかじめ測定し，(2) 式の ARX モデルの未定係数 ai ，bi j をオフラインで一括型最小 2 乗法を用いて同定しておく．係数同定後，オンラインで計測される測定値 u j (k)，K(k) を使用し，(11) 式のように y(k) の推定値 ŷ(k) を算出する．
Mj
7 X
X(16)

本研究での提案法

u3 (k) = K(k)L
X000
β1−a− 00 N b −a0 IN
α0α0

u2 (k) = N(k) · I(k)

L
X
i=1

ai Q̂(k − i) +Mj
3 X
X
j=1 i=0bi j u j (k − i)

(21)

また，前節と同様，差圧推定式 (16) についても ARXモデルを適用し，推定式を得た．

P̂(k) = −Mj
L
3 X
X
X00
ai Q̂(k − i) +bi j u j (k − i)i=1(22)j=1 i=0

3 実験
実験は Fig. 1 のような模擬循環系を用いた．定常流
人工心臓ポンプには，Terumo 製 CAPIOX 遠心ポンプを用いた．専用ドライバを使用し，モータ回転数の制御を PC を介して行い，モータ回転数 N ，電流 I ，流量Q，差圧 P (Pout − Pin ) をサンプリング周期 50ms にて計測した．
静特性を計測する際には，Fig. 1 において，ポンプ流入部にあるクランプにより前負荷を変化させた．
また，動特性を検証するデータとして，回転数 N を方形波状に変化させたものを用いた．ARX モデルの未定係数を同定するための同定用データと，推定用データの 2 種類を用意した．
得られたデータを用いて，オフラインで推定された各推定法の比較を行った．推定精度の評価の指標とし

(15)

70

て，真値と推定値間の相関係数 R ならびに推定誤差 を (24) 式のように算出した．

PnD

PnD

0.12

¯k=1 {y(k) − ȳ}{ŷ(k) − ŷ}2
k=1 {y(k) − ȳ}PnD(23)¯2k=1 {ŷ(k) − ŷ}

vtn
DX
1
={y(k) − ŷ(k)}2nD k=1,

R= q0.16



0.08


0.04(24)0500

ここで，y は流量 Q，または差圧 P を表すこととする．
nD は使用データのサンプル数，ȳ，ならびに ŷ¯ は，それぞれ y，ŷ の平均値を示す．

1000
1500
Rotational Speed [rpm]2000

Fig. 3: α, β の回転数特性

Venous6

Pump:B

R=0.97^
Q[L/min]

54
3

PoutPin(N)(N)2
Pump

Flowmater

1

Arterial tank

2345
Q [L/min]6Fig. 1: 模擬循環系の構成
Fig. 4: (15) 式による静特性の流量推定結果

4

結果・考察

4.1

静特性モデルによる推定0Fig. 2 に実験で得られた capiox ポンプの流量–電流静特性を示す．
(3) 式における α，β と回転数との関係を Fig. 3 に示す．(15) 式における a，b，α0 ，β0 は，a = 1.0，b = 0.87，
α0 = 9.0 × 10−6 [A·min/L]，β0 = 3.0 × 10−4 [A] となった．
(15) 式を用いた静特性の推定結果を Fig. 4 に示す．実測値，推定値間の相関係数 R は R = 0.97 となった．

P [mmHg]0.35

R=0.99200150100^

0.3Current [A]

0300
250

0.25

5000.2

-50

0.15

0.05
012345Flow [L/min]

6050

100 150 200 250 300 350
P [mmHg]

800rpm1000rpm
1200rpm1400rpm
1600rpm
1800rpm0.100

また，差圧においては，(16) 式における a ，b ，α0 ，0000
β0 はそれぞれ，a = 1.0，b = 2.7，α0 = −2.0 ×
010−6 [A/mmHg]，β0 = 3.0 × 10−9 [A] と求められた．(16)式を用いた静特性の推定結果を Fig. 5 に示す．実測値，推定値間の相関係数 R は R = 0.99 となった．以上から，(15)，(16) 式による推定モデルは静特性をよく記述できているといえる．

Fig. 5: (15) 式による静特性の差圧推定結果

74.2

動特性モデルによる推定

4.2.1

静特性モデルとの比較

回転数が時間的に変化する動特性のもとで，推定精度の比較を行うために，静特性モデル (15) 式，(16) 式

Fig. 2: 電流–流量静特性71と，時系列法による動特性モデル (21) 式，(22) 式を用いて推定を行った．

数 R をそれぞれ Table 1 に示す．
本提案法では入力パラメータ数を減らしたことによる大幅な推定精度の低下は認められず，ほぼ同程度とみなすことができる．誤差は流量において 2.1%，差圧において 2.6%であり，実用に十分な推定精度であると考えられる．

まず，同定用データを用い，(17) 式における未知の係数を求めた．ただし，入出力次数 M j ，L はすべて 3としている．次に，推定用データを用い，(15) 式，(21)式で流量を推定した．

(15) 式での推定結果を Fig. 6 に，(21) 式による推定結果を Fig. 7 に示す．推定誤差 r，相関係数  は，Fig. 6で r = 0.56， = 2.3L/min，Fig. 7 で R = 0.95， =0.09L/min となった．Fig. 6 において，過渡状態の推定精度が低いが，Fig. 7 においては，定常状態だけでなく，過渡状態においても高い精度を保つ推定が行えることがわかる．

Table 1: 従来の提案法との流量・差圧推定精度の比較
入力数 7 の (11) 式による推定 入力数 3 の (21)(22) 式による推定流量
 0.08L/min (1.9% )(平均 4.2L/min) R
0.97
差圧 2.2mmHg(2.5% )(平均 88mmHg) R
0.98

Q[L/min]

また，差圧についても同様の手順で推定を行ったところ，推定誤差 r，相関係数  は，(16) 式で r = 0.66，
 = 23mmHg，(22) 式で R = 0.95， = 2.3L/min となり，流量の場合と同様の結果が得られた．

 

84001020304050607080

Time[s]Fig. 6: 静的モデル (15) 式による動特性の流量推定結果

Q[L/min]

本研究では，定常流型人工心臓において，回転数が動的に変化する状況での差圧・流量推定に関し，時系列法が有効であることを示した．また，ARX モデルのベースとなる静特性モデルの指数項部分を改良することによって，従来法とほぼ同等の推定精度を保ちながら，モデルの入力パラメータ数を減らすことが可能となった．このことにより，測定機器の異常や外乱が入力に加わった場合でも，推定値に大きな誤差が生じにくくなることも期待できる．
今後は，実用化に向け，動物実験での検討，他の遠心ポンプを用いた場合での本推定手法の汎用性や，運用時に未知の係数を同定する手法を検討していく必要がある．
[1] http://www.fda.gov/bbs/topics/ANSWERS/2004/ANS01280.html[2] 日経産業新聞，2004 年 1 月 21 日号[3] A.Funakubo, S.Ahmed, I.Sakuma, Y.Furuki, Flow rate andpressure head estimation in a centrifugal blood pump, ArtifOrgans, 26, pp. 985–990, 2002

405 おわりに

参考文献

 

80.09L/min (2.1% )
0.95
2.3mmHg(2.6% )
0.9701020304050607080[4] 築谷朋典, 赤松映明, 西村和修, 磁気浮上式遠心ポンプの流量測定法, 人工臓器, 26, pp. 98–102, 1997

Time[s][5] T.Kitamura, Y.Matsushima, T.Tokuyama, S.Kono,
K.Nishimura, M.Komeda, M.Yanai, T.Kijima, C.Nojiri,
Phisical Model-Based Indirect Measurements of BloodPressure and Flow Using a Centrifugal Pump, Artif Organs,24, pp. 589–593, 2000

Fig. 7: 入力数を減らした ARX モデル (21) 式による動特性の流量推定結果

4.2.2[6] P J Ayre, N H Lovell, J C Woodard, Non-invasive flow estimation in an implantable rotary blood pump: a study considering non-pulsatile and pulsatile flows, Phisiol. Meas,
24, pp. 179–189, 2003

従来法との比較

ARX モデルの入力パラメータ数を減らしたことによる推定精度への影響をみるために，(21) および (22) 式を用いた本提案手法と，(11) 式を用いたこれまでの提案手法について，推定精度の比較を行った．

[7] M.Yoshizawa, T.Sato, A.Tanaka, K.Abe, H.Takeda,T.Yambe, S.Nitta, Y.Nose, Sensorless estimation of pressure head and flow of a continuous flow artifical heart basedon input power and rotational speed, ASAIO Journal, 48,
pp. 443–448, 2002

まず，共通の同定用データを用い，双方の ARX モデルの未知の係数を決定した．次に，共通の推定用データを用いて，(21) 式 (または (22) 式) および (11) 式で推定し，精度を比較した．ただし，ARX モデルの入出力次数 M j ，L は全て 3 に統一した．

[8] A.Tanaka, M.Yoshizawa, Y.Aizawa, P.Olegario, K.Abe,
T.Yambe, S.Nitta, Estimation of Pressure Head and Flow
Rate in a Continuous-Flow Artificial Heart –In VIVO
Evaluation–, sice02, 2002

このときの流量ならびに差圧の推定誤差  ，相関係72