야마다 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_{max}$는 근육의 등척성 최대 수축력으로, 실험 데이터에 기초하여 100N으로 되어 있다.

근육의 활동도 다이나믹스는, 등척성 근육 수축에 의한 근력을 조사한 것에 의해 모델화되어 있다. 등척성 근육 수축에 의한 근력은, 신경근 접합부에 의한 $Ca^{2+}$농도와 가교형성에 의존한다. 이 모델에 대해 활성도 다이너믹스 $a$는, 칼슘 방출과 재흡수를 반영한 Low-pass 필터에 의해 칼슘 다이나믹스 $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 )}$으로 치환하는 것으로 행하였다. 단 $N$은 근방추의 활동도의 변화율을 결정하는 파라미터이다. 최대 변위 $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}$으로 했다.

3.3.5 태아의 나이에 기초한 태아 물리 모델

작성한 태아 모델에 기초하여, 원하는 태아의 나이를 시뮬레이션 하기 위한 방법에 대해서 설명한다. 태아의 나이의 변화에 기초하여 변화시키는 파라미터는 체중과 근력 두 가지이다. 구체적으로는 이하의 두 개의 수순으로 행한다. 첫 번째로, Archie 등에 의한 각 태아의 나이에 대한 태아의 체중의 분포의 평균치 [14]를 사용하여 체중을 결정하고, 태아 모델을 스케일링한다. 두 번째로, Stickland에 의한 사람 태아의 봉공근 (역주: 넙다리빗근, Sartorius)의 단면적과 체중의 관계성에 대해서 데이터 [468]을 기초로 각 근육의 근단면적을 구해, 근단면적과 근력의 관계성 [213]으로부터 대상이 되는 태아의 나이에 대한 각 근육의 최대 수의 수축력을 결정한다. 이 모델화에서는 각 신체 부위의 비율은 일정하다고 한 것에 유의해 둔다. 작성한 태아 모델의 예시를 Fig. 3.15에 보인다.

3.3.6 태아 모델을 위한 동역학 계산과 충돌 계산

태아 모델의 물리 시뮬레이션을 ODE 위에서 안정적으로 행하기 위해 위에서의 문제점과 강구한 대책에 대해서 설명한다. 생겨난 문제는 크기가 작은 것, 폴리곤 모델인 것 두 가지에 기인한 것이다.

첫 번째로 크기가 작기 때문에 생겨난 문제와 대책에 대해서 설명한다. 이번의 모델화하고 권동을 확인한 태아 모델은, 가장 작을 때는 12주일 때 두둔 길이(역주: CRL, 머리에서 엉덩이까지의 길이)는 76.1mm이다. 이 12주 태아 모델의 두부는 질량 $5.1\times 10^{-3}\mathrm {kg}$, 두부의 관성 텐서가 관절 축에 대해서 약 $2\times 10^{-7}\mathrm {kgm^{2}}$이다. 이 정도의 order가 되면, ODE의 내부에서 행해지는 가동역 한계를 갖는 관절이나 충돌 계산의 불안정성이 현저하게 된다. 이것은 시뮬레이션의 시간 스텝을 작게 하는 것이 하나의 해결책이지만, 그것만으로는 충분하지 않다. 행한 대책은 두 개가 있어, 하나는 관절에의 점성 저항의 도입, 한 개는 단위계의 변경이다. 관절에의 점성 저항은 근사 후퇴 오일러법에 의해 모델화했다. 단위계의 변경은 길이와 시간의 단위계를 변경하는 것에 의해 행했다.

두 번째로, 태아 모델이 삼각 메시에 의해 표현된 폴리곤 모델이기 때문에, 충돌 계산이 적절하게 수행되지 않는 문제점이 있었다. ODE에서는 충돌 계산은 Geometry의 클래스 전 조합에 대해 개별으로 정의되어 있지만, 폴리곤 끼리의 충돌 계산은 매우 불안정하다고 말하고 있어, 실제 태아 시뮬레이션에서도 충돌 계산이 적절하게 행하여지지 않고 통과해버렸다. 한편, 폴리곤과 '캡슐 또는 구' 와의 충돌 계산에 대해서는 안정적인 경향이 있어, 태아 시뮬레이션의 범위에서는 통과하는 문제와 발생하지 않았다. 여기서 충돌 계산에 쓰인 기하학 형태에 대해서는, 사지를 캡슐 및 구로 근사했다 (Fig. 3.16). 얼굴에 대해서는, 눈과 귀, 입 등의 부분에 대해서 Re-topology를 행하는 것으로 한층 더 안정화와 고속화를 보였다.

3.4 태아의 척수 및 피질 모델

3.4.1 개요

척수 및 피질 모델에 대해서 설명한다. 척수 모델은, 근육의 고유 감각 피드백을 입력으로 하고, 근육에의 운동 지령을 출력으로 한다. 이 모델은 발달 초기의 자발 운동의 모델을 상정한 것이고, 이 가정의 타당성에 대해서는 제 4장의 실험에서 검증을 행한다. 피질 모델은 감각 피드백을 입력으로 하고, 가소성에 기초한 신경계 및 인지발달을 시뮬레이션하기 위해 사용한다. 구체적으로는, 제 5장에서 피질의 일차 영역인 체성 감각 영역, 제 6장에서는 뇌 전 영역을 시뮬레이션 할 때, 척수 모델에 피질 모델을 추가하는 형태로 시뮬레이션이 행해진다. 이 두 개의 모델에 대해서는, 기본적으로는 같은 뉴런 모델 및 가소성을 따르고, 체성 감각 영역 모델의 확장이 뇌 전 영역의 시뮬레이션이다고 말할 수 있다. 거기서 첫 번째로, 피질 모델의 기본 요소인 뉴런 모델과 가소성의 규칙에 대해서 설명을 행한다. 뉴런 간의 시냅스 결합에 대해서는 동일 영역 내를 포함한 근방에서는 랜덤을 기본으로서 생선된다. 대뇌 피질 전 영역의 시뮬레이션에 대한 영역 간의 결합에 대해서는, 신생아의 뇌 이미지 해석에 의해 도출을 행하여, 이것에 기초하여 결정했다. 거기서 두 번째로는, 이 영역 간의 결합 관계 도출에 대해서, 그 방법과 결과에 대해서 설명을 행한다. 또한 피질 모델에 대해서는, 대뇌 피질 출력과 운동의 관계성은 본 연구의 대상 외이므로, 피질으로부터 척수 및 근육에의 출력은 모델화되어 있지 않다.

3.4.2 척수신경회로 모델

척수 신경 회로 모델은, 태생기의 전반의 매우 빠른 초기 발달의 시기로부터 시작하는 자발 운동에 대해서, 자발 운동의 기전 및 신체성과의 관계성에 대해서 찾아보는 것을 가능하게 하는 모델일 필요가 있다. 그를 위해, 명시적인 운동 표현이나 근육 간 협조 관계를 상정하지 않고, 생물학적 지식에 기초하여 최소한의 가정만으로부터 구성되는 모델이 원해진다. 하지만, 자발 운동과의 관계성이 많은 문헌에서 지적되어 있는 신경 진동자에 대해서는, 단일의 근육에의 출력계로서 모델에 포함한 경우는, 본 연구의 목적으로부터 타당한 가정이라고 생각된다. 상기의 이유로부터, 척수 신경 회로 모델에 대해서는, Kuniyoshi and Sangawa에 의한 모델 [263]을 사용했다. 이 모델은, 근육 장력, 근육 길이의 감각 정보를 입력으로 하여, 근육의 활동도를 결정하는 운동 지력을 출력으로 한다. 각 근육은 신경진동자를 포함한 단위 유닛에 의해 독립으로 제어되고, 동일한 근육 레벨 이외의 반사 회로는 들어있지 않고, 신경계 레벨에서 근육 간의 결합은 일절 갖지 않는다. 즉, 명시적인 운동 표현이나 근육 간의 관계성을 결정하는 요소는 일절 들어있지 않다. 이 모델에 의한 연구에서는, 운동 창발의 예비 실험만이 보고되어 있어, 평가는 정성적인 것으로 생물의 자발 운동과의 관계성에 대해서는 취급되어 있지 않다. 여기서 본 연구에서는 이 모델을 사용하여, 창발한 운동에 대해서 정량화 및 사람의 자발 운동 특징과의 비교를 행하는 것으로부터 시작하는 것으로 했다.

모델의 단위 유닛은, 근육, $\alpha$, $\gamma$운동 뉴런, 척수 내 개재성 감각 뉴런 S0, 그리고 신경진동자 모델으로부터 구성된다 (Fig. 3.17). 이 모델 간의 결합에 대해서는 생물학적 지식을 고려한 지연과 이득이 설정되어 있다. 이하, 각 모델에 대해서 설명을 행한다.

척수 내 뉴런의 다이나믹스는, 고양이나 사람의 $\alpha$운동 뉴련에 대한 실험적인 데이터 [26, 84, 305]에 기초하여, Band-pass 필터로서 모델화 된 He에 의한 운동 뉴런 모델 [181]이 사용되고 있다.

• ${\frac {m(s)}{i(s)}}=k_{\text{m}}{\frac {(1+s/33+(s/33)^{2})}{1+2(s/58)+(s/58)^{2}}}\qquad (3.29)$

여기서, $i$는 입력신호 $[\mathrm {nA} ]$, $m$은 출력 신호 $[\mathrm {p/s} ]$이다. 입출력 신호 둘 다 0에서 1까지의 범위이고, 범위에서 나간 경우에 대해서는, 그 범위에 들어가도록 역치처리를 하고 있다. 또한, $k_{\text{m}}$상수 파라미터이며, 입력 전류 $1\mathrm {nA}$에서의 발화율을 나타낸다. 실험 데이터 [26]으로부터 $k_{\text{m}}=1.5\mathrm {p/s/nA}$으로 설정되어 있다. 척수 개재 뉴런의 다이나믹스에 대해서는 연구가 적기 때문에, 운동 뉴런의 다이나믹스에 의해 대체된 경우가 많다 [181, 85]. 본 연구에서도, 척수 개재 뉴런에 대해서는 상기의 운동 뉴런 다이나믹스를 사용한 모델을 채용하고 있다.

신경 진동자 모델에 대해서 설명한다. 본 연구에서는, 생리학적 상세한 모델이 아닌 비선형 진동자에 의한 모델을 사용했다. 이것은, 신경 세포 모델에서부터 구성된 신경 회로망과 리듬적인 활동 패턴의 관계성을 조사하는 것을 주로 하고 있는 것, 거기에, 전신의 규모에서 시뮬레이션을 행하는 것은 계산 비용이 높은 것이 있다 [159, 63, 420, 508]. 한편, 비선형 진동자는 신경 세포 레벨보다 추상도가 높은 단순한 모델이기는 하지만, 리듬적인 활동 패턴과 운동의 관계성을 조사하는 이론 연구의 많은 경우에 실제로 사용되고 있어, 계산 비용도 낮다 [172, 95, 148]. 게다가, 이러한 진동적인 활동 패턴은 발달 초기부터 척수에 대해서 이미 관찰되어지는 것이 닭이나 쥐의 배를 조사한 연구에 의해 보고되어 있다 [176, 399, 536, 13, 177]. 본 연구에서는, 어떻게 진동적인 활동이 신경 회로망으로부터 창발하는가에 대해서는 다루지 않고, 그 활동을 베이스로서 환경과 상호작용하는 것에 의해 발생하는 감각/운동정보구조나 그 발달적 변화에 대해서 다루는 것을 목적으로 하고 있기 때문에, 신경 진동자 모델으로 충분하다고 판단했다.

신경 진동자 모델으로서는, Asai 등에 의한 모델 [18, 1, 17]에 기초한 Kuniyoshi and Sangawa [263]에 의한 것을 사용했다. 먼저, Asai 등에 의한 신경 진동자 모델을 설명한다. 이 모델은, 두 개의 계가 상호로 제어하여 맞추는 것에 의해 상호적 활동 패턴을 생성하는 하프 센터 형태 [60, 224, 311]의 신경 진동자 모델이다. 또한, 모델의 권동은 임상적으로 관찰되어 있는 활동 패턴의 경향을 재현 가능하다고 되어 있다 [18, 17]. 각각의 계는, Fitzhugh-Nagumo (Bonhoeffer-van der Pol: BVP) 방정식에 의해 모델화 되어 있다 [132, 340]. Fitzhugh-Nagumo 방정식은, 상세하게 생물 물리학적 타당성이 높은 Hodgkin-Huxley 방정식 [194]와 비교하면 단순화 되어 있지만, 신경 세포의 흥분성 및 진동막에서 보여지는 전형적 성질을 표현 가능하다고 되어 있다 [355, 535]. 두 개의 Fitzhugh-Nagumo 방정식으로부터 구성된 신경 진동자 모델은 다음의 식으로부터 성립한다 [18].

• ${\frac {\mathrm {d} x_{1}}{\mathrm {d} t}}=c{\Bigl (}x_{1}-{\frac {1}{3}}x_{1}^{3}-y_{1}+z_{1}{\Bigr )}+\delta (x_{2}-x_{1})\qquad (3.30)$
• ${\frac {\mathrm {d} y_{1}}{\mathrm {d} t}}={\frac {1}{c}}(x_{1}-by_{1}+a)+\epsilon x_{2}\qquad (3.31)$
• ${\frac {\mathrm {d} x_{2}}{\mathrm {d} t}}=c{\Bigl (}x_{2}-{\frac {1}{3}}x_{2}^{3}-y_{2}+z_{2}{\Bigr )}+\delta (x_{1}-x_{2})\qquad (3.32)$
• ${\frac {\mathrm {d} y_{2}}{\mathrm {d} t}}={\frac {1}{c}}(x_{2}-by_{2}+a)+\epsilon x_{1}\qquad (3.33)$

여기서, $x$,$y$는 각각 막 전위, 불응성을 나타낸다. $z$는 외부 전류 자극에 해당하고, $z_{1}$,$z_{2}$의 관계성을 변화시키는 것으로 출력의 성질을 변화시키는 것이 가능하다. 구체적으로는, 동일의 경우 교대성의 주기적인 출력 $x_{1}$,$x_{2}$가 관찰되고, 다른 경우에는 그것이 복잡한 권동을 나타내는 것이 보고되어 있다 [17]. 또한, 파라미터는 $a=0.7$, $b=0.675$, $c=1.75$, $\delta =0.013$, $\epsilon =0.022$의 값으로 설정되어 있다.

외부 입력에 의해 주기적인 출력과 복잡한 출력의 전환이 가능한 점에 착목하여, Kuniyoshi and Sangawa는 초기의 자발 운동의 모델으로서, 신경 진동자 모델을 다음으로 정의했다.

• $\tau {\frac {\mathrm {d} x}{\mathrm {d} t}}=c{\Bigl (}x-{\frac {1}{3}}x^{3}-y+z{\Bigr )}+\delta (s_{S0}-x)\qquad (3.34)$
• $\tau {\frac {\mathrm {d} y}{\mathrm {d} t}}={\frac {1}{c}}(x-by+a)+\epsilon s_{S0}\qquad (3.35)$

Asai 등에 의한 모델과의 차이는 한 쪽의 Fitzhugh-Nagumo 방정식만을 사용하고, 그 한 쪽으로부터 입력을 입력 피드백 $s_{S0}$으로 한 것이다. Asai 등의 모델을 포함한 하프 센터 형태 모델은, 굴근(역주: 굽어지는 근육)과 신근(역주: 펴지는 근육), 또는 좌우의 다리의 교대성의 진동적 패턴을 대상으로 하고 있다. 그를 위해, 로코모션 중의 근육 활동을 두 개의 그룹으로 크게 나눠 설명하려는 문제점이 있다 [158, 311]. 게다가, 전신 운동을 대상으로 하면 문제가 더욱 심각해져서, 어떻게 근육 군을 분리할까도, 그러한 근육 간 또는 좌우의 협조 관계를 앞서 가정하는 것이 초기 발달을 대상으로 한 모델으로서는 타당성이 있는가도 문제가 되어 온다. Kuniyoshi and Sangawa는 이 문제는, 한 쪽의 반대가 되는 Fitzhugh-Nagumo 방정식을 신경계가 아닌 신체/환경 측에서 구하는 것으로 해결을 시험해 보았다. 신경 진동자 모델에 있는 단일의 Fitzhugh-Nagumo 방정식은, 그 출력에 의해 신체를 움직이고, 환경에서 움직여서, 그 결과로서 감각 피드백을 얻는다. 감각 피드백이 랜덤이 아닌, 신체/환경에 의해 구조화 되어, 그 관계성이 일정하지 않다는 것을 보이면, 그것은 Asai 등에 의한 상호로 일정한 관계성으로 결합된 두 개의 계와 동질한 권동을 하고 있다, 라는 가정에 기초하고 있다. 이 Kuniyoshi and Sangawa에 의한 신경 진동자 모델은, 다양한 운동의 창발을 확인했다, 라는 정성적인 보고에 멈춰 있기 때문에, 본 연구에서는 이 모델을 사용하여, 정량적인 평가와 상기의 가정이 성립하고 있는가의 이론적 검증을 행하는 것으로부터 시작할 필요가 있다.