第 4 章 病的状態が 2 足歩行モデルの motor module に与える影響 56
4.2 方法
4.2.1 骨格モデル
頭部,上肢,体幹 (HAT), 下腿,足部を有する2次元7リンクの骨格モデルを構築した (図4.1). 身体形状のそれぞれのセグメントは表 4.1に示した [48].
表4.1 骨格モデルパラメータ
HAT 大腿 下腿 足部
長さ (m) 0.800 0.4165 0.418 0.066
重さ (kg) 44.070 6.500 3.055 0.975
慣性モーメント (kgm2) 5.823 0.117 0.048 9.347×10−5
各関節は線形な粘性要素を持つピン構造の関節としてモデル化した. 股関節,膝関節,足 関節の粘性係数はそれぞれ1.09, 3.17, 0.943 Nms/radとした [49]. また,膝関節,足関節 は過屈曲,過伸展防止の制限をかけた. 膝関節,足関節のバネ・ダンパ係数は両関節ともに
2.0×103Nm/rad,5.0×102Nms/radとした. 地面に足が着いた際,踵またはつま先に床 反力が生じる. そのバネ・ダンパ係数は水平成分が5.0×103N/m,1.0×103Ns/mであ り, 垂直成分は2.5×104N/m,5.0×102Ns/mとした.
ニュートン・オイラー法によって運動方程式を導いた [44, 45]. 具体的な式を付録A.8に 示す. 一般的な式の形は以下のように記述される.
¨
x=P(x)F+Q(x,x,˙ Fg(x,x),˙ Tr(F’m)), (4.1) x は7リンクの位置と角度のベクトル(21×1)であり, P は(21×12)の行列, Fは並進 力を表すベクトル(21×1), Qは外的な力・トルクを表すベクトル(21×1),Fgは床反力 のベクトル(8×1)である. Trは関節トルクのベクトル(6×1)である. F’mは筋張力の ベクトル(18×1)であり,次の小節で詳しく説明をする. 関節トルクの式を付録A.9に示 す. 床反力Fgの要素は以下の式で表わされる.
Fgxrt =
{−kgx(xrt−xrt0)−bgxx˙rt {yrt <0}
0 それ以外,
Fgxrh =
{−kgx(xrh−xrh0)−bgxx˙rh {yrh <0}
0 それ以外,
Fgyrt =
{−kgy(yrt−0) +bgyf(−y˙rt) {yrt <0}
0 それ以外,
Fgyrh =
{−kgy(yrh−0) +bgyf(−y˙rh) {yrh <0}
0 それ以外,
Fgxlt =
{−kgx(xlt−xlt0)−bgxx˙lt {ylt <0}
0 それ以外,
Fgxlh =
{−kgx(xlh −xlh0)−bgxx˙lh {ylh <0}
0 それ以外,
Fgylt =
{−kgy(ylt−0) +bgyf(−y˙lt) {ylt <0}
0 それ以外,
Fgylh =
{−kgy(ylh −0) +bgyf(−y˙lh) {ylh <0}
0 それ以外,
(4.2) (xrt, yrt), (xrh, yrh), (xlt, ylt), (xlh, ylh) はそれぞれ右と左のつま先と踵の座標である. xrt0,xrh0, xlt0, xlh0はそれぞれ,一歩ごとに右と左のつま先または踵が地面に着いた瞬間
のxrt, xrh, xlt, xlh を保持する. Fgxrt, Fgxrh, Fgxlt, Fgxlh は右と左の水平方向の床反力, Fgyrt, Fgyrh, Fgylt, Fgylhは右と左の垂直方向の床反力をそれぞれ示す. kgx, bgxとkgy, bgyはそれぞれ水平方向と垂直方向のバネ・ダンパ係数である. 特に水平方向では足が地 面に着いてからの変位とその速度に応じて床反力が生成される.
関節の拘束力は以下の運動学的拘束を受ける.
C(x)¨x=D(x,x)˙ (4.3)
C は(12×21)の行列, Dは(12×1)のベクトルである. 具体的な式を付録A.10に示す. 式(4.1)を式(4.3)に代入すると,
F= [C(x)P(x)]−1[D(x,x)˙ −C(x)Q(x,x,˙ Fg(x,x),˙ Tr(F’m))] (4.4) が得られ, Fを取り除くために,式(4.4)を式(4.1)に代入すると,
¨
x=P(x)[C(x)P(x)]−1[D(x,x)˙
−C(x)Q(x,x,˙ Fg(x,x),˙ Tr(F’m))]
+Q(x,x,˙ Fg(x,x),˙ Tr(F’m)).
(4.5)
が得られる.この式に筋張力F’mを与えると, 2足歩行の筋骨格系システムの動きが得ら れる．各パラメータの値を付録A.11に示す.
Pattern Formation (PF) 筋骨格システム
フィードバック 反射ループ
HAT
大腿
下腿 足部
① ②
⑤ ⑥③ ④
⑦
⑧ ⑨
CPG
神経システム
CPGループ Rhythm Generator (RG)
Motor Modules Motor Modules
左 右
図4.1 提案した神経筋骨格モデルの概略図. 左側： 筋骨格モデル. 各筋は, 1,大臀 筋 (GM); 2,腸腰筋 (IL); 3,大腿二頭筋長頭 (BFL); 4,大腿直筋 (RF); 5,大腿二頭 筋短頭 (BFS); 6,中間広筋 (VA); 7, 腓腹筋 (GC); 8,ヒラメ筋 (SO); 9,前脛骨筋 (TA)を表す. 頭部,上肢,体幹をまとめてHATとした. 右側: 神経システム. Rhythm Generator (RG) ネットワークとPattern Formation (PF) ネットワークを有する階 層的CPGモデル. RGでは, 4つのニューロン内2つ (左側の破線円)が位相性,周期 性のある活動パターンを作り出す.その活動を基にして, PFでは左下肢筋を活動させ る5つのmotor moduleを生み出す (神経システムから筋骨格システムへの矢印). ま た, RGでの他の2つのニューロン (右側の破線円)活動をによって, PFでは右下肢筋 を活動させる5つのmotor moduleを生み出す(矢印は記載していない).
4.2.2 筋モデル
左右の下肢に各 9 筋, 計 18 筋をモデル化した (図4.1A). モデル化した筋は大臀筋 (GM), 腸腰筋 (IL),大腿二頭筋長頭 (BFL),大腿直筋 (RF), 大腿二頭筋短頭 (BFS), 中 間広筋 (VA),腓腹筋 (GC),ヒラメ筋 (SO),前脛骨筋 (TA) である.筋はα運動ニューロ ンの活動に対応し,筋の張力–長さ,筋の張力–速度に依存した収縮を起こす.式(2.5)の筋 モデルを用いた [46].
Fm′ = ¯FmCE·k(ξm)·h(ηm)·αm+FmPD+FmPE,
k(ξm) = 0.32 + 0.71 exp[−1.112(ξm−1)] sin[3.722(ξm−0.656)], h(ηm) = 1 + tanh(3.0ηm),
FmPD =cPDm L˙m,
FmPE =kPEm {exp[15(Lm−L¯m)]−1},
(4.6) Fm′ は左右の筋張力であり, mは左右の筋を示す (m= 1, ...,18). ¯FmCEは収縮要素 (CE) による最大筋張力, k(ξm)は筋の張力-長さ関係, h(ηm)は筋の張力-速度関係, αmはα 運 動ニューロンに対応する刺激信号 (0≤αm ≤1), FmPD とFmPE はそれぞれ粘性と弾性力, ξm と ηm はそれぞれ筋の最適長L¯m によって正規化された筋の長さと収縮速度である. ξm=Lm/L¯m, ηm= ˙Lm/L¯˙m であり, Lm, ˙Lm, L¯˙mはそれぞれ筋の長さ, 収縮速度,最大 の筋収縮速度である. cPDm は粘性係数, kmPE は弾性係数である. 各パラメータの値を付録 A.1に示す.
4.2.3 神経モデル
Matsuokaニューロンモデル(式(2.6))を使用し,上位であるRhythm Generator (RG) ネットワークと下位であるPattern Formation (PF)ネットワークで構成される14ニュー ロンの階層型CPGモデル [25]を構築した(図4.2). 具体的には以下で表される．
τiu˙i =−ui+
∑n i′=1
wi,iCPG′ yi′−βvi+uθi+ Feedi, τi′v˙i =−vi+yi,
yi = max(0, ui).
(4.7)
ui はi番目ニューロンの内部状態 (1≤i≤14), viはi番目ニューロン自身への抑制であ り, τi, τi′ は時定数, β はvi の係数, wCPGi,i′ はi′ 番目ニューロンからi番目ニューロンへ の結合重みを示す. uθi は外部からの入力, Feediは筋骨格系からのフィードバック情報で ある. RGネットワークは互いに抑制し合う4つのMatsuokaニューロンで構成され, リ ズムを生み出す. 一方, PFネットワークは互いにかつ自身を抑制する5つのMatsuoka ニューロンを左右の下肢に配置し,それらはそれぞれ左右5つのmotor moduleに対応す る. 5 つのmotor module は外部入力として, RGのリズミックな活動に反応し, 異なっ
たフェーズで一つ一つ順番に活動する (図4.2B). 各パラメータを付録A.12に示す. PF ニューロンは運動指令をα運動ニューロンへ送り,筋を活動させる. α運動ニューロンは 姿勢反射,交叉性伸展反射などの様々な反射からもフィードバック情報を受け取る. α 運 動ニューロンの出力αm (1≤m≤18)と反射の出力Reflexmは以下の式で与えられる.
αm =wconditionm 1
1 + exp(−4(∑5
n=1wm,nα PFn+ Reflexm)), (4.8) Reflexm =
{∑
n′(cm,n′θn′ +c′m,n′θ˙n′) + POSm Fg>0,
0 otherwise, (4.9)
wconditionm , wm,nα , cm,n′, c′m,n′ は重み係数であり, PFn はPFニューロンの出力である. θn′ は各関節角度(n′ ∈ {hip, knee, ankle} ), POSmは姿勢制御, Fgは床反力を示す. 各 パラメータの値を付録A.13に示す.
これまでにも, 2足歩行リズムを生み出す階層型 CPGモデルは提案されているが, あ らかじめ決められた位相のみのパターンを生み出すものや各ニューロンの結合が少なく, 生み出されるパターンがある程度決まっているものであった [48, 49, 107]. 本モデルでは, RGは一定のリズムを保ちながら, PFニューロン同士の結合が複数あるため, 生み出され るパターンは多彩である. そのため,特異的な歩行パターン形成が可能となる.
(B) RG
PF
4 3 1
2
5 6
7
8 9 10
1
2
3 4
1 2 3 4
Time [sec]
1 2 3 4
Time [sec]
RG
PF
Motor modules
(A)
右脚 左脚
図4.2 構築した階層型CPGモデル. (A) Rhythm Generator (RG)ネットワークと Pattern Formation (PF) ネットワークの概略. 黒い円は抑制性、矢印は興奮性の結合 を示す。破線は RGからPFへの出力を示す. (B) RGとPF間での階層性. RGの 2つのニューロンはPFにおいて異なる位相を持った5つの活動パターンを生み出す. RGからの出力が特定の値 (破線円)に達すると, PFのニューロンは活動する. 1つの PFニューロンは1つのmotor moduleを表す. 片側下肢の神経制御のみ記載.
4.2.4 パラメータ探索
構築したモデルは未知なる 56個のパラメータ(wm,nα , cm,n′, c′m,n′, POSm, Feedi) を 有しており,安定した歩行させるために調整が必要であった. そこで, GA [108]を用いて これらのパラメータを調整した. 評価関数は以下の式を使用した.
wdD+wsS →max (4.10)
wd, wsは重み係数,Dは本モデルが転倒するまでに歩いた距離を表す. Sは歩行中の歩数 である. 初期に 50の個体をランダムに作成した. 0.7 の確率で1点交差, 0.1の確率で突 然変異が生じるように設定した. 構築したモデルは安定した5秒間の2足歩行を達成する のを確認した. 本来,ヒトは試行錯誤的に歩行を獲得するとされ,それをGAで模擬した.
4.2.5 健常と病的歩行モデル
健常歩行は GAを用いて, 56 個のパラメータを調整しシミュレートした. その際, 式 (4.8) において, wconditionm = 1.0 と設定した. 健常歩行パターン獲得後に, 脳梗塞のよ うな病的状態でもシミュレートした. 特に, m = 1,· · · ,6に対して, wconditionm = 0.8, m = 7,8,9に対して, wconditionm = 0.6と設定し,脳梗塞患者で観察されるようなMEPs
の弱化 [106]をシミュレートした. その後,歩行モデルは即座に転倒してしまったので, 2
つの仮想的なリハビリテーション過程を検討した. 1 つ目の過程は患者が反射機構のみ による歩行の再獲得を想定し, 28個の反射パラメータ(cm,n′, c′m,n′)をGAを用い,再調 整した. この過程を反射代償モデルとした. もうひとつの過程は患者が制御器自体の順 応により歩行を再獲得することを想定し, 48個のCPG へのフィードバックパラメータ (Feedi)をGAを用い,再調整した. この過程をCPG 代償モデルとした. 両過程ともに, 歩行は可能となったが, motor moduleはそれぞれ異なった変化を示した.
4.2.6 シミュレーションの実装
GAは並列計算ライブラリである Message Passing Interface (MPI) を用いて実装し た. 全てのプログラムはC言語を用い,微分方程式は4次のRunge-Kutta法で求解した. 時間刻み幅は0.1ミリ秒とした. これは先行研究 [44–46, 48] を参考にして,精度よく求解 できる値に設定した. 加えて,構築したモデルが直立を保てるようにHATを筋張力と比 例微分制御法,両方で制御をした. 比例微分制御法はHATの関節角度と角速度を用いた. 関節角度と角速度係数はそれぞれ, 1.0×103, 1.0×102 とした. これは上肢などを使用し, 支持物に掴まり歩行していることを意味している. 本研究はこのような一定の条件で,歩 行中の下肢筋活動パターンであるmotor moduleの変化を検討した.