야마다 2015 박사 논문

3.3 동역학 시뮬레이션을 위한 태아 물리 모델과 태내 환경 모델

3.3.1 개요

여기서는, 작성한 태아의 근골격 신체 모델을 강체의 동역학 시뮬레이션에서 돌리기 위한, 태아 모델의 작성에 대해 설명한다. 먼저, 신체 모델을 21개의 신체 부위로 분할하고, 골격에 기초한 관절을 설정하고, 강체의 링크계로서 모델화를 행하였다. 여기에 수반하여, 단일 신체 부위 모델 내에 붙어 있는 근육에 대해서는, 강체 시뮬레이션이므로 제거했다 (역주: 관절과 연관되지 않은 근육은 제거했다). 다음으로 신체를 구동하고, 환경과의 상호작용의 결과를 피드백으로서 얻기 위한 운동기와 감각기의 각각에 대해서 모델화를 행한다. 구체적으로는, 운동기로서의 근육을, 감각기로서는 근육의 고유 감각인 근방추, 골지근 기관 및 촉각의 감각수용기에 대해서 모델화를 행했으므로, 각각에 대해서 설명을 행한다. 또한, 자궁막과 양수로부터 구성되는 태내 환경 모델에 대해서 설명한다. 작성한 32주에 해당하는 태아 모델을 원하는 태아의 나이의 모델로 변환하는 방법에 대해서도 설명한다. 마지막으로 ODE를 사용하여 태아 모델으로 강체 시뮬레이션을 행하기 위해 행한 개선점에 대해서 설명한다. 작성한 태아 모델에 의한 운동학 시뮬레이션의 외관을 Fig. 3.10에, 시뮬레이션의 스냅샷을 Fig. 3.11에 나타낸다.

3.3.2 태아의 강체 링크 모델과 근육 배치

태아의 근골격 신체 모델에 기초한 태아의 신체 모델을 21개의 신체부위와, 그것을 연결하는 20개의 관절, 36자유도의 강체 링크 모델으로서 모델화했다 (Table 3.4 and Fig. 3.12). 신체의 분리는 골격에 기초하여 행하고, 피부를 거기에 응하여 분리하고, 관절은 골격에 대해서 설정된 것을 사용했다. 각 강체의 질량이나 관성 모멘트에 대해서는, 밀도가 일정하다 가정하고 계산하는 것으로 구했다. 관절 가동 영역에 대해서는, 태아에 대해서 조사한 것은 찾을 수 없었다. 따라서, 사지 또는 어캐에 대해서는 와타나베 등이 보고한 신생아의 것 [548]을, 몸통과 머리에 대해서는 어른의 것 [546]으로 대용했다. 이것들의 관절 가동영을 참고하여, 이 범위를 넘지 않도록 태아 모델의 관절 가동각 한계를 설정했다 (부록 C).

근육에 대해서 설명한다. 각 신체 부위 모델에 대해서는 단일 또는 복수의 뼈가 할당되어 있어, 이것에 기초하여 근육은 고정되어 있는 뼈가 소속하는 신체 부위에 대해 시뮬레이션 속에서는 고정된다. 신체 부위의 분리에 수반하여 같은 신체 부위 내에 닫혀 있는 근육에 대해서는 제외했다. 그 결과로서 최대 494개의 근육이 강체 시뮬레이션 위에서 구동 가능하게 되었다 (Fig. 3.13). 이 근육 중에서는, 복수의 선분으로 모델화된 근육의 일부가 시뮬레이션 가능하게 된 것도 있다는 것에 주의한다. 일부를 포함한 근육을 정리하여 지운 경우 또는, 관절 가동각 한계까지 움직였을 때에 길이가 변화하지 않는 근육을 포함한 경우, 이 근육 개수로부터 제한 것이 된다.

3.3.3 근육 및 근육의 고유감각 모델

근육 모델

본 연구에서 사용하는 근육 모델은, Hill의 특성 방정식을 사용하여 모델화[193]된 He에 의한 근육 모델[181]을 베이스로 한다. 근육 모델은, 0부터 1의 연속 값인 무차원수의 운동 지령을 입력으로서 힘을 발생하는 것이 가능하다. 이 모델은, 자연 길이 0.036m의 근육을 베이스로 한 모델이기 때문에, 자연장이 다른 근육에 대해서도 사용하는 것이 가능하도록 확장해서 사용했다. 전신에 배치되어 있는 각 근육은, 자연 길아와 최대 수의 수축력을 파라미터로서 동일의 다이나믹스를 따른다.

먼저, He에 의한 근육 모델의 설명을 행한다. 각 근육이 발생하는 힘은 능동 수축력 $F_{act}$와 수동 수축력 $F_{p}$의 힘의 합에 의해 부여된다. 근육의 능동 수축력 $F_{act}$는 근육의 활성도 $a$, 길이 $l$, 속도 $v$에 의존한다. 이 모델은 이하의 식으로서 모델화된다.

• $F_{act}=F_{max}\;a\;f_{l}(l)\;f_{v}(v)\qquad (3.2)$

다음으로 활성도 $a$, 수의 수축력의 근육 길이 의존 성분 $f_{l}(l)$ 및 속도도 의존성분 $f_{v}(v)$의 다이나믹스에 대해 각각 설명한다. $F_{act}$는 근육의 등척성 최대 수축력으로, 실험 데이터에 기초하여 100N으로 되어 있다.

근육의 활동도 다이나믹스는, 등척성 근육 수축에 의한 근력을 조사한 것에 의해 모델화되어 있다. 등척성 근육 수축에 의한 근력은, 신경근 접합부에 의한 $Ca^{2+}$농도와 가교형성에 의존한다. 이 모델에 대해 활성도 다이너믹스 $a$는, 칼슘 방출과 재흡수를 반영한 로우패스 필터에 의해 칼슘 다이나믹스 $q$ (Eq. (3.3))과 토로포닌에 의한 칼슘 결합에 의한 가교 사이클을 표현한 비선형 필터 (Eq. (3.5))에 의해, 하기와 같이 모델화 되어 있다 [463, 47].

• ${\frac {\mathrm {d} q}{\mathrm {d} t}}=18.8496(m-q)\qquad (3.3)$
• $u={\frac {q^{2}}{k^{2}+q^{2}}}\qquad (3.4)$
• ${\frac {\mathrm {d} a}{\mathrm {d} t}}=b(u-a)\qquad (3.5)$
• $b=22(1-0.51a)^{2}\qquad (3.6)$

여기서, $m$은 운동 뉴런으로부터의 신호로, $m$, $q$ 둘 다 단위는 [p/s]이다. 등척성 근육 수축과 운동 뉴런의 출력과의 사이에 시그모이드의 관계성이 있는 것이 모델 연구에 의해 보여지고 있다 [185, 135]. 이 지식에 기초하여 파라미터 $k$는, 운동 뉴런의 발화율이 30[p/s]에서 안정 상태의 등척성 근육 수축의 반절의 값을 출력하는 시그모이드가 되도록 $k=30[p/s]$로 정해져 있다. 상태의존 파라미터인 $b$는, 근육의 히스테리시스와 큰 근력으로부터의 부드러운 이완을 반영한 파라미터이다. 근육 활동도 $a$는 0에서 1의 범위를 갖는 무차원 수이다.

다음으로, 근육 길이 의존 성분과 속도도 의존성분의 다이나믹스에 대해 설명한다. 이 모델은, Zajac에 의한 근육의 일반화 모델 [537]을 He 등이 확장한 것이다 [182, 183]. 근육 길이 의존성분 $f_{l}(l)$은 Otten에 의해, 최적 근육 길이 $L_{opt}$에 의해 근육 길이 l의 상대 길이에 의존하는 형태로 하기에 의해 모델화 되어 있다 [363].

• $f_{l}(l)=\mathrm {exp} {\biggl \{}-{\biggl (}{\frac {|(l/L_{opt})^{2.8286}-1|}{0.6182}}{\biggr )}^{2.3680}{\biggr \}}\qquad (3.7)$

$L_{opt}$는 Spector 등의 실험 데이터 [454]와 일치하도록, 0.05m의 값으로 설정되어 있다. 속도도 의존성분 $f_{v}(v)$는, 0부터 1.8의 값을 갖고, 다이나믹스는 수축시와 신장시에 각각 다른 쌍곡선 함수에 의해 이하와 같이 모델화 되어 있다.

• $f_{v}(v)={\begin{cases}(1+7.3125v)/(1+4.0625v),&v\geq 0{\text{: Stretching}}\\(1+v)/(1-2.25v)&v\leq 0{\text{: Contracting}}\end{cases}}\qquad (3.8)$

여기서, $v$는 최대 근육 수축 속도 $V_{max}$에 의해 정규화된 속도이고, 여기서의 $V_{max}=0.25m/s$로 되어 있다.

다음으로 근육의 수동 수축력 $F_{p}$의 다이나믹스에 대해 설명한다. 수동 힘의 다이나믹스는, 근육 길이와 그 속도에 의존하는 형태로 모델화 되어있다. 근육 길이에 의존하는 부분은, 자연 길이보다 길게 되면 지수함수로 수동 수축력은 증가하고, 어떤 길이를 넘으면 선형으로 모델화 되어있다. 근육의 속도에 의존하는 부분에 대해서는, 확률적 섭동을 사용해서 근육의 수동 수축력에 대한 점탄성의 관계성에 대해 서술한 Kirsch 들에 의한 연구 [248]에 기초해 모델화 되어 있고, 이하를 따른다.

• $F_{p}=f_{p}(l)+B_{p}(l)v\qquad (3.9)$
• $f_{p}(l)={\begin{cases}0&(l\leq 0.036)\\0.1663(\mathrm {exp} (208.2(l-0.036))-1)&(0.036<l<0.057)\\13.0076+2743(l-0.057)&(l\geq 0.057)\end{cases}}\qquad (3.10)$
• $B_{p}(l)={\begin{cases}0.34624\ \mathrm {exp} (208.2(l-0.036))&(l<0.057)\\27.4&(l\geq 0.057)\end{cases}}\qquad (3.11)$

마지막으로, 자연 길이 및 최대 수의 수축력이 서로 다른 근육에 대해서도 사용하는 경우에 대해 설명한다. 자연 길이가 서로 다른 점에 대해서는, 대상이 되는 근육의 자연 길이 및 근육 길이를 각각 $l'$,$l_{0}'$이라 하면, $l\leftarrow 0.036l'/l_{0}'$으로 치환하는 것으로 구했다. 등척성 최대 근육 수축력 $F_{max}$에 대해서는, 대상의 그것을 $F_{max}'$으로서 수동 힘과 능동 힘의 합력을 구한 뒤에 $F_{max}'/F_{max}$을 곱하는 것으로 구했다.

근육의 고유 감각기 모델

근육의 감각기로서는, 근육 길이와 근육 장력을 반영한 감각 피드백을 반환하는 근방추와 골지건 기관의 모델이 조합되어 들어있다. 이 두개에 대해서도, 정합성의 관점으로부터, 근육 모델과 같은 방식으로 He에 의한 모델을 채용했다 [181]. 이하, 각각의 모델에 대해서 설명한다.

골지건 기관의 모델은, Houk and Simon에 의한 근력에 대한 응답을 모델화한 모델 [202]에 대해, Prochazka and Gorassini가 고양이의 로코모션을 예시로서 통상 운동 상황하에서의 응답을 설명 가능한 레벨으로 확충한 모델 [389]에 기초한다. 골지건 기관의 출력 $r_{Ib}$는, 전달함수의 형태로 Band-pass 필터로서 이하의 식으로 구해진다.

• $r_{Ib}(s)=K_{\text{GTO}}{\frac {(1+s/0.15)(1+s/1.5)(1+s/16)}{(1+s/0.2)(1+s/2)(1+s/37)}}\qquad (3.12)$

이 모델에서는 골지건 기관이 $\alpha$운동 뉴런과의 결합은 직접 결합하는 것을 가정하여 모델화되어 있고, 전달 지연과 이득 $K_{\text{GTO}}$는 상수로서 다루어지고 있다 [181]. 여기서는, 근육 장력이 $F_{max}$가 될 때에, 출력이 $1{\text{nA}}$가 되도록 이득 $K_{\text{GTO}}$는 $0.01{\text{nA/N}}$으로 설정되어 있다.

다음으로, 근방추의 모델에 대해 설명한다 [178, 181]. 먼저 $x$의 근육 길이, $y$를 방추 속 근육 길이, $z$를 근방추의 센서 영역 길이라고 하면, 해부학적 제약에 의해 이하의 관계성이 성립한다.

• $x=y+z\qquad (3.13)$
• ${\dot {z}}={\dot {x}}-{\dot {y}}\qquad (3.14)$

근방추 출력 $r_{Ia}\ \mathrm {p/s}$는 이 $z$를 사용하여, 이하의 전달함수로 표현되는 다이나믹스를 통하는 것으로 얻을 수 있다.

${\frac {r_{Ia}(s)}{l_{sp}(s)}}=K_{\text{MS}}{\frac {(1+s/7.23)(1+s/74.07)}{(1+s/12.46)(1+s/123.28)(1+s/250)}}\qquad (3.15)$

여기서, $K_{\text{MS}}=3.2\times 10.5^{5}\mathrm {p/s/m}$의 상수로서 취급되어 있다. 근 길이 $x$는 시뮬레이션에서 직접 얻는 것이 가능하므로, Eq. (3.14)에 의해 ${\dot {y}}$가 얻어지면 수치적분에 의해 근방추의 센서 영역 길이 $z$가 얻어져서, Eq. (3.15)에 의해 근방추 출력을 얻는 것이 가능하다. 아래에서, ${\dot {y}}$를 얻는 방법을 설명한다.

방추 속 근육 길이 $y$의 미분 방정식은 다음과 같이 모델화 되어 있다.

• ${\dot {y}}=f_{\text{vi}}^{-1}{\Bigl (}{\frac {{\frac {K}{E}}\ z}{f(l)}}{\Bigr )}\qquad (3.16)$
• $f(l)={\begin{cases}0.0009197\mathrm {exp} (400l),&(l<0.0025)\\l&(l\geq 0.0025)\end{cases}}\qquad (3.17)$
• $l=x-z-c\qquad (3.18)$

여기서, $c$는 근방추의 느슨함이 시작되는 길이이며, $c=0.035\mathrm {m}$의 상수로서 취급되고 있다. 또한, $K$는 센서 영역의 강성(剛性), $E$는 수축 영역의 강성이다. 센서 영역의 방추 속 근육에 대한 강성 비율 $K/E$는, 근방추에의 입력 신호 $i\ \mathrm {p/s}$에 의존하고, 이하의 식에 의해 부여된다.

• ${\frac {K}{E_{0}(i)}}={\frac {99}{1+0.08i}}\qquad (3.19)$

거기에, 근육 길이의 속도가 양수일 때의 힘과 속도의 관계 $f_{\text{vi}}({\dot {y}})$는 다음의 형식으로 일반화 가능하다.

• $f_{\text{vi}}({\dot {y}})={\begin{cases}(1+C_{1}{\dot {y}})/(1+C_{2}{\dot {y}})&({\dot {y}}>0)\\(1+C_{3}{\dot {y}})/(1-C_{4}{\dot {y}})&({\dot {y}}<0)\end{cases}}\qquad (3.20)$

여기서, $C_{1}=85.5\mathrm {s/m}$, $C_{2}=22.7\mathrm {s/m}$, $C_{3}=1\mathrm {s/m}$, $C_{4}=62.8\mathrm {s/m}$의 상수 값으로서 모델화 되어 있다. $f_{\text{vi}}({\dot {y}})$의 역함수를 구하는 것은 용이하므로, 이상의 식에 의해 ${\dot {y}}$는 계산 가능하다.

상기의 근방추 모델은, 근육 모델과 같은 방식으로 자연 길이 $0.036\mathrm {m}$의 특정의 근육에 대한 모델이기 때문에, 전신의 근육에 적용하기 위해서 확장할 필요가 있다. 각 근육에 대해 자연 길이 뿐만이 아니라 근육 길이의 변화하는 범위가 각각 다르기 따문에, 이 두 개의 값을 기초로 변환을 행하였다. 구체적으로는, 근육 길이 및 자연 길이를 $l'$, $l_{0}'$, 근육의 자연 길이로부터의 최대 변위가 $D$일때, $l\leftarrow 0.036{\bigl (}{\frac {N-1}{D}}(l'-l_{0}')+1{\bigr )}$으로 치환하는 것으로 행하였다. 단 $D$에 대해서는, 각 관절에 대해 관절이 움직일 수 있는 한계(관절 가동각)의 상한 하한을 랜덤으로 정했을 때에 취하는 자세를 충분히 많이 뽑아, 그 때의 최대 변위를 기초로 결정했다.

3.3.4 태내 환경과 촉각 모델

태내 환경은 자궁막과 양수로부터 모델화되어 있다. 태아 모델은, 자궁막에 의한 중심에의 구심력과, 양수에 의한 유체저항 및 부력을 받는다. 자궁막에 의한 구심력은 각 강체에 할당된 매 접촉점마다 정의된다. 접촉점 $i$의 법선 벡터를 ${\boldsymbol {n}}_{i}$, 표면적을 $A_{i}$, 위치 벡터를 ${\boldsymbol {p}}_{i}$라고 하면, 자궁막에 의한 구심력 ${\boldsymbol {F}}_{u,i}$는 이하의 식에 의해 모델화된다.

• ${\boldsymbol {F}}_{u,i}=-f_{u,i}A_{i}{\boldsymbol {n}}_{i}\qquad (3.21)$
• $f_{u,i}={\begin{cases}(K_{u}r_{i}+C_{u}{\dot {r}}_{i})+{\frac {({\boldsymbol {n}}_{i}\cdot {\boldsymbol {q}}_{i})_{+}}{\lVert {\boldsymbol {q}}_{i}\rVert }}&(r_{i}>0)\\0&(r_{i}\leq 0)\end{cases}}\qquad (3.22)$
• $r_{i}=\lVert {\boldsymbol {q}}_{i}\rVert -R_{u}\qquad (3.23)$
• ${\boldsymbol {q}}_{i}={\boldsymbol {p}}_{i}-{\boldsymbol {p}}_{u}\qquad (3.24)$

여기서, $f_{u,i}$는 접촉점 $i$에서의 압력을 나타내고, ${\boldsymbol {p}}_{u}$는 자궁막의 중심 위치 벡터, $R_{u}$는 자궁막 저항이 발생하는 반경을 나타낸다. $K_{u}$, $C_{u}$는 각각 탄성과 점성의 상수 파라미터이다. $(x)_{+}$는, $x$가 양수인 경우에 $x$를 반환하고, 음수인 경우 0을 반환한는 사상이다. 본 실험에서는, $R_{u}=0.05\ \mathrm {m}$, $K_{u}=1.0\ \mathrm {Pa/m}$, $C_{u}=0.01\ \mathrm {Pa\cdot s/m}$으로 했다. 또한, 강체에 속한 접촉점에 분배된 표면적 $A_{i}$의 합계는, 그 강체의 표면적에 일치한다.

유체 저항에 대해서도 같은 방식으로 매 접촉점에 계산이 행해진다. 각 접촉점의 유체저항 $F_{a,i}$는 이하의 식으로 정의된다.

• ${\boldsymbol {F}}_{a,i}=-f_{a,i}A_{i}{\boldsymbol {n}}_{i}\qquad (3.25)$
• $f_{a,i}={\frac {1}{2}}\rho C_{D}({\boldsymbol {\nu }}\cdot {\boldsymbol {n}}_{i})^{2}\qquad (3.26)$

여기서, $f_{u,i}$는 접촉점 $i$에서의 유체 저항에 의한 압력을 나타내고, ${\boldsymbol {\nu }}$는 소속된 강체의 속도, $\rho$는 양수의 밀도, $C_{D}$는 항력 계수를 나타낸다. 본 연구에서는, $\rho =1010\mathrm {kg/m^{3}}$[204], $C_{D}=1.0\ \mathrm {Pa\cdot m/kg\cdot s^{2}}$으로 설정했다. 부력은 각 강체의 체적을 밀도 일정으로서 계산되어 있다.

태아 시뮬레이션에서는, 촉각 자극의 감각 피드백을 모델화 하기 위해서, 전신에 800개의 촉각감각기 모델을 배치했다 (Table 3.5 and Fig. 3.14). 이 분포는 Weinstein에 의한, 사람의 각 신체 부위의 이점변별역 (역주: 두 점 변별 경계) [518]에 기초하여 결정되어 있다. 각 신체 부위에 대해서는, 그 분배 수를 랜덤하게 배치하는 것으로 촉각의 배치를 행했다. 촉각기의 모델으로서는, 사람의 촉각의 기계 수용체의 하나이며, 주로 물체의 접촉 검지에 기여하고 있다고 하는 메르켈 세포를 모델화했다. 본 연구에서 취급하는 것은 강체 시뮬레이션이므로, 피부의 진동이나 비틀림을 상정하고 있지 않기 때문에, 남은 세 개의 기계수용기의 모델화는 이번에는 생략했다. Freeman 등이 보인 사람의 감각 수용기가 고유의 주파수 응답 특성을 갖고 있다고 하는 지견 [133]에 기초하여, Shirado 등은 메르켈 세포를 컷오프 주파수 50Hz의 Low-pass 필터로서 모델화 하고 있다 [442]. 본 연구에서도, 이 모델을 메르켈 세포의 모델으로서 채용했다.

감각 세포 모델의 입력 전압은, 태내 환경에서는 자궁막에 의한 압력, 유체 저항, 강체끼리의 접촉에 의한 압력 세 가지가 입력된다. 자궁막에 의한 압력, 유체 저항에 대해서는, 전신에 배치된 촉각점마다 계산된 $f_{u}$, $f_{a}$가 사용된다. 접촉에 의한 압력에 대해서는, 물리 시뮬레이션에 의해 계산된 강체끼리의 접촉력을 접촉점마다 분배하는 형태로 계산한다. 접촉점 $i$가 속한 강체 $j$가, 강체 $k$와 충돌한 것에 의해 강체 $j$가 받은 힘을 ${\boldsymbol {F}}_{p,jk}$로 하고, 접촉점 $i$와 강체 $k$와의 거리를 $l_{i,k}$라고 하면, 접촉에 의한 압력 $f_{p,i}$는, 이하의 식에 의해 계산된다.

• $f_{p,i}={\frac {1}{A_{i}}}{\frac {d_{i,k}}{\textstyle \sum _{n}\displaystyle d_{n,k}}}\lVert {\boldsymbol {F}}_{p,jk}\rVert \qquad (3.27)$

• $d_{i,k}=A_{i}(L_{c}-l_{i,k})_{+}\ (-{\boldsymbol {F}}_{p,jk}\cdot {\boldsymbol {n}}_{i})_{+}\qquad (3.28)$

여기서, $L_{c}$는 접촉력을 분배하는 최소의 거리를 나타내는 상수 파라미터이다. 본 연구에서는, $L_{c}=5\mathrm {mm}$으로 했다.