予測医療に向けた階層統合シミュレーション高木周

予測医療に向けた階層統合シミュレーション
東京大学
機械工学専攻
&
大学院工学系研究科
バイオエンジニアリング専攻
高木
周
2
階層統合生体力学シミュレーションとその医療応用
シミュレーションによる病態予測と治療支援
方針：・ミクロ（分子生物学，細胞動態）とマクロ（病態予測，治療法
検討）を結びつける．
・マクロ(病態）からミクロ（分子）へトップダウンのアプローチ．
（病態をシミュレートしながら，ミクロ因子の影響を調べる）
階層統合シミュレーションによる新しい予測医療の構築
分子
細胞
血小板
血小板
vWF
vWF
組織
血小板
GP1ｂ
vWF
コラーゲンの露出した血管壁
臓器
心筋梗塞
For the Design and Development
of High Intensity Focused Ultrasound
Therapy (HIFU) Device
High Intensity Focused Ultrasound therapy
•Prostate cancer
•Breast cancer
http://www.prostatecancercentre.co.uk/treatments/hifu.html
HIFU therapy has been developed for the treatment of deeply-placed cancer.
Brain cancer
Liver cancer
http://www.imasonic.com/
Displacement of focal point due to the
reflection and refraction of ultrasound at
the interfaces of bones
Control of the focal point
by an array transducer
Recent Development of HIFU Device
Products
(1) InSightec(GE) : Exablate 4000
Brain Treatment, Phase Control by Array Transducer
(1024 elements),
MR-Guided, (MR-ARFI under development)
(1) InSightec: Exablate 4000
(2) Phillips: Treatment of uterus myoma,
MR-Guided
(3) Siemens: collaborating with
Chongqing(重慶) HIFU, MR-Guided
(2)Phillips;Temperature
Monitoring by MRI
(4) Supersonic Imaging: collaborating
with Prof. Fink (Inserm), MR-ARFI
(4) 256 Channels
--------------------------------------------------------------(5) University of Michigan, Profs, C. Cain, B. Fowlkes, Z. Xu
Histotripsy ( Cavitation induced Treatment)
(5) Plaque: Utilization of Cavitation, to Control the Focused Area
blown away
HIFU simulator
Comparison of the focal point
Without phase delay
Pmax=1.7MPa
0
2
Pressure [MPa]
With phase delay
Pmax=2.2MPa
0
2
Pressure [MPa]
Experimental Apparatus
Stage
Stage controller
40 mm
Trigger
56ch PZT
Transducer

z
y
x
o
LabVIEW™
Trigger
61ch Amplifier
Needle hydrophone
Water
Oscilloscope
Acrylic plate (3 mm)
c=2660m/s
Ultrasound through acrylic plate (20deg)
Acrylic plate is arranged with the angle of 20 degrees.
6
0.3
5
0.25
4
w/o phase shift
y [mm]
3
0.2
2
1
0.15
0
0.1
-1
-2
-3
-15
0.05
-12.5
-10
-7.5
-5
x [mm]
-2.5
0
1
6
0.3
5
0.25
4
w/ phase shift
y [mm]
3
0.2
2
1
0.15
0
0.1
-1
-2
-3
-15
0.05
-12.5
-10
-7.5
-5
x [mm]
-2.5
0
1
HIFU Ｓｉｍｕｌａｔｉｏｎ for Focus Control for Breast Cancer using MRI data
MRI image
Numerical model
Fat
Water
•
•
•
56ch transducer
Focal distance 100 mm
Frequency 2 MH
•
•
•
•
3D orthogonal mesh
Size 120x110x110 mm
1200x1100x1100 grid
8x4x4 sub domains
Parenchymal
Speed of
[m/s]
Density
[kg/m3]
Fat
1465[1]
985[1]
lacteal gland
1547[2]
1032[2]
connective
tissue
1615[1]
1090[1]
Sound
[1] ICU REPORT 61, [2] T. D. Mast,, 2000
w/o control
w/ control
Parkinson’s Disease Simulation
Close collaboration among teams of
K. Doya(OIST),
Y. Nakamura(Univ. Tokyo),
T. Nomura(Osaka Univ.),
and S. Takagi (Univ. Tokyo)
12
研究実施体制
(筋骨格‐脳神経系階層統合シミュレーション)
パーキンソン病の階層統合シミュレーションによる病態理解と
治療応答予測および運動と脳の関係のモデル化
大脳基底核-皮質
-小脳神経回路
大規模シミュレーション
（OIST－銅谷）
脊髄反射の
神経回路モデル
（東大-中村，高木）
筋繊維レベルからの
筋収縮の有限要素解析
（東大―高木）
サルの実験データ
に基づくパーキン
病脳モデル
(ミクロモデル）
（OIST―銅谷，
生理研－橘）
患者データに基
づくパーキン病モ
デル
（マクロモデル）
（阪大―野村
刀根山病院ー
佐古田）
全身の分布質量筋モデルによる姿勢制御予測（東大―中村）
13
ミクロスケール実験データ：
パーキンソン病モデルサルの脳シグナル (Tachibana et aｌ. (2011)）
マクロスケール実験データ：
パーキンソン病における姿勢反射障害発生機序の仮説モデル (野村（阪大））
‐ PLoS One (2009)
‐ Hum Mov Sci (2008)
‐ J Theor Biol (2012)
‐ Math Biosci (2013)
従来仮説
連続的制御モデル
前後(mm)
前後(mm)
左右(mm)
左右(mm)
間欠制御の
崩壊？
PSD (rad2/Hz)
（シミュレーション）
適切なタイミングで
制御OFF
柔かな運動
UPDRS=53(OFF)
PSD (rad2/Hz)
新しい仮説
間欠制御モデル
パーキンソン病患者
若年健常者
Freq (Hz)
Freq (Hz)
（シミュレーション）
常に制御ON
高ゲインで安定化
硬い運動
今までの実施状況1 (脳神経系シミュレーション) (2011-2014年度)
Diesmann (ユーリッヒ研究所）, 五十嵐（OIST銅谷チーム）ら
2013年7月，ISLiMで開発された脳神経系シミュレータ NESTを
用いて，「京」上で， 10兆個の結合を持つ大脳皮
質局所神経回路のシミュレーションに成功
－世界最大の脳神経シミュレーション－
・ 17億3,000万個の神経細胞
・ 10兆4,000億個のシナプス結合
マーモセットなどの小型霊長類の全脳規模に匹敵
ヒトの脳機能解明に向けた第一歩
プレスリリース http://www.riken.jp/pr/topics/2013/20130802_2/
5
16
16
今までの実施状況2 (2011-2014年度)
筋骨格-神経系階層統合シミュレーション（中村，高木，銅谷）
脊髄神経モデルと筋骨格ワイヤモデルの接続
(2) K-Body+HI-Muscle
STIMULATION
(1) PyNEST + MUSIC + K-Body
全身骨格モデルと
筋有限要素モデルの接続
外力
(刺激)
伸張反射
シミュレーション
実験データを用いた筋・神経の
モデルパラメータの同定
Excitatory Synapses
Inhibitory Synapses
下肢の筋活動度推定計算
二頭筋SN発火→IN発火→三頭筋の筋活動抑制
ニューロン番号
(サイズ順)
IN
三頭筋SN
二頭筋SN
三頭筋MN
二頭筋MN
運動(MN)・感覚(SN)・介在ニューロン(IN)のスパイク列
(サイズの小さいものから順に発火する)
Numerical Method
Suitable for Medical Image Data
(Full Eulerian FSI simulations)
Gaehtgens et al.(1980) Blood Cells., 6, 799.
Background
Fluid-Structure coupling analysis of living body
Diagnostic image (CT, MRI)
Voxel data (Volume Fraction of Constituents)
representing multi-component geometry
without Mesh Generation
(Finite-Difference or Finite-Element) Simulation
on
Eulerian frame
http://www.rehab.research.va.gov/
jour/02/39/3/ledouxf06.jpg
Full Eulerian approach for Fluid-Structure Interactions
vs.
Lagrangian
Eulerian
How is the two-phase distinguished ?
s
1
Solid
Solid
0.5
Fluid
Fluid
0
by the boundary of mesh
by solid volume fraction
How is the solid deformation described?
by the displacement of
material points themselves
by left Cauchy-Green
deformation tensor
・ Sugiyama, Ii et al. (2011) J. Comput . Phys., 230, 596.
・ Ii, Sugiyama et al. (2011) Int. J. Numer. Meth. Fluids, 65, 150.
Basic equations for conservation and constitutive law
  v  0,  m   t v  ( v   ) v    p    σ m ,
Cauchy 's stress tensor in mixture form
σ m  (1   s )σ f   s σ s ,
solid volume fraction
solid (visco-hyperelasto):
fluid (Newtonian):
σ f  2  f D ,
Mooney-Rivlin model:
1
D  v  vT  ,
2
strain rate
σ s  2 c1B   2 c2 (tr(Β )Β  B  B )  2  s D .
left Cauchy-Green deformation
 xi
B  F  F , Fij 
,
X j
T
・ Sugiyama, Ii et al. (2011) J. Comput . Phys., 230, 596.
Basic equations for kinematics of solid motion
solid volume fraction
 t s  ( v   ) s  0,
left Cauchy-Green deformation tensor
 t B  (v   ) B  L  B  B  L ,
T
where
L  v .
T
・ Sugiyama, Ii et al. (2011) J. Comput . Phys., 230, 596.
Validation I
Comparison with available numerical data
Gao & Hu (2009)
J. Comput. Phys. 228, 2132.
particle-particle interactions
in a shear flow
full Lagrangian
Particle-particle interactions in a shear flow
● material points (tracers)
(roll-over and bounce-back modes)
Movie
s=f=1, Lx=8, Ly=4, =20, C1=0, C0=80
contour: vorticity (spin)
Vtop=1
-3
Fluid
Solid
Solid
Vbottom=1
1
Comparison with Gao & Hu（JCP2009)'s results

Gao & Hu (full Lagrangian)
 Present (full Eulerian)
NxxNy=1024x512
t=0.2
t=0.8
t=4.0
t=5.6
t=30.0
t=33.2
Extension to 10-particle system
● material points (tracers)
Movie
s=f=1, Lx=8, Ly=4, =20, C1=0, C0=80
contour: vorticity (spin)
Vtop=1
-3
Vbottom=1
1
Feature: full Eulerian approach
・makes it easily possible to perform FSI simulations
with complex geometry
Initial voxel data
FSI simulation
without mesh generation/reconstruction
・ Sugiyama, Ii et al. (2010) Comput . Mech., 46, 147.
・ Nagano, Sugiyama et al. (2010) J. Fluid Sci. Tech., 5, 475.
・ Takagi, Sugiyama et al. (2012) J. Appl. Mech., 79, 010911.
The Super-Fast FSI Solver: ZZ-EFSI
ZZ-EFSI achieved the actual
speed of 4.5 PETA FLOPS!!!
Software availablle at
http://www.islim.org/islim-dl_e.html
Development of Multiscale
Thrombosis Simulator
Multiscale modeling of initial stage of thromobosis
Blood flow
(continuum mechanics)
Protein(wWF)-Protein
(GP1b) binding
(stochastic process)
•Substance diffusion
•Metabolic reaction and activation
•Morphology change
Molecular interaction
(molecular dynamics)
Blood flow simulation
• [0,44]×[0,22]×[0,22] [m], D=20 [m],
• 30 RBCs, 10 platelets (Ht ≈ 20 [%])
• <u> ≈ 2 [mm/s]
320x160x160 mesh, MPC-RICC at RIKEN (256 cores)
・ Ii, Sugiyama et al. (2012) J. Biomech. Sci. Eng.,7, 72.
Stochastic Monte Carlo Simulation of
Cell adhesion on vessel wall
R+L
kf
kr
C
R： Receptor
L： Ligand
C： R-L complex
kf： Binding React. Rate Cnst.
kr： Dissosc. React. Rate Cnst.
Receptor－Ligand Bind. Force
f b  Function of  l  l0 
Jianrong Li, et. al(2009)
Platelet
Vessel Wall
l
l0
Ligand-Receptor Binding Model
Stochastic model with energetic elasticity
Eyring (1935) J. Chem. Phys., 3, 107.
Bell (1978) Science, 200, 618.
Dembo (1988) Proc. R. Soc. Lond. B, 234, 55.
Hammer & Apte (1992) Biophys. J., 63, 35.
P f  1  exp(  k f  t )  R f  formation
Pr  1  exp(  k r  t )  Rr  breakage
( R f , Rr  [0,1] : Random numbers )
Luo et al. (2007) Blood, 109, 603.
forward reaction rate
Formation

( l  l0 ) 2
k f ( l )  k f 0 exp   ts

2 k bT

reverse reaction rate
f   p ( l  l0 )



Model Parameter
 p  10  4 [N/m],  ts  0.9 p [N/m]
Breakage

( l  l0 ) 2
k r ( l )  k r 0 ex p  ( p   ts )

2 k bT




l0  60 [nm], k r 0  3 [s  1 ]
Fox et al. (1988) J. Biol. Chem., 263, 4882.
Arya et al. (2005) Biophys. J., 88, 4391.
Kim et al. (2010) Nature, 466, 992.
Motion of an adhering platelet
dimensionless time t* =  t
 = 1000 (1/s)
 = 3000 (1/s)
 = 6000 (1/s)
Simulated Example : Wall-Adhering Platelet moving in Shear Flow
Experiment
Reininger et al., blood, 107
(2006) 3537
Simulation
(= 3000 s-1)
シミュレーションモデル構築のための実験
(後藤信哉(東海大）より提供）
Ht value
© Goto
(Tokai Univ.)
Number of adhering platelets
シミュレーション用データ取得のための
フローチャンバーを用いた実験
血小板粘着に対する赤血球の影響に関する実験データ
（後藤（東海大）ら 2012)
Experimental data by Goto & Tamura.
Adhesion of Platelets in a Shear Flow
Ht = 5.5 %
 = 800 (s-1) t  200 (ms)
Ht = 15.3 %
Ht = 21.9 %
Numerical Simulations for the Experiment by Goto et al.
Ht = 0%
Ht = 20%
Domain Size: 144mx72mx72m
# of grid points: 480x240x24
本年度：実験系と同じ流路サイズのシミュレーション
Ht = 20%
Comp. extent:
400mx100mx100m
Num. grid points: 2,048x512x512
2,048 nodes (16,384 cores) on the K computer
injured
wall
platelet adhesion
41
Effect of arteriosclerosis size ( Hematocrit Value: 20% )
Effect of arteriosclerosis size ( Hematocrit Value: 0% )
抗血小板薬のモデリング
循環器系
・抗血小板薬の薬効再現のためのモデリングとシミュレーション
45
実測データを用いた循環器系統合シミュレーション
開発済み
開発中
局所３次元計算
術後予側シミュレーション
Left
Left
Left
Left
Blood pressure(mmHg)
140
external illiac artery (Preoperation)
posterior tibial artery (Preoperation)
external illiac artery (Postoperation)
posterior tibial artery (Postoperation)
0D-1D-3D連成モデル
120
神経系モデルの
導入
100
80
60
4
5
Time(s)
Femoral artery (preoperation)
Femoral artery (Postoperation)
Graft
25
Blood flow(ml/s)
20
15
10
5
0
4
-5
5
Time(s)
患者データの活用
およびデータ同化
0D-1D全身血管シミュレーション
血圧計開発
マルチスケール血栓症シミュレーション
課題３．予測医療に向けた階層統合シミュレーション
まとめ
現状
・医療応用に向けてソフトウェアの準備ができ，本格計算が
始まったところ
・真の意味での医療応用はこれからが勝負．
医学・医療さらには看護に対して何を創出できるか．
今後の予定，
・超音波診断・治療器設計に向けた計算
・血栓シミュレータによる血小板粘着過程の再現と
抗血小板薬の薬効評価（P2Y12阻害薬他）
・脳神経系と筋骨格系の統合シミュレーション（パーキン
ソン病の再現）
シミュレーションが支援する新しい医学・医療の創出ヘ…..
Research Achievement (2015)
Won a Distinguished Simulation Award!!!
UT-Heart （Hisada et al.)
謝辞
Principal Investigators
後藤信哉（東海大・医），姫野龍太郎（理研），横田秀夫（理研），
松本洋一郎（東大・工），久田俊明（東大・新領域），
大島まり（東大・情報学環），劉浩（千葉大・工），
天野晃（京大・情），和田成生（阪大・基礎工），
小野謙二（東大・生研），松澤照男（北陸先端大・情報セ），
山口隆美（東北大・工），野村泰伸（阪大・基礎工），
岡澤重信（広島大・工），中村仁彦（東大・情報理工）
銅谷賢治（沖縄科技大）, Huaxiong Huang (York Univ.)，
長谷部光泉（東海大・医）
Main contributors for software development):
杉山和靖（東大），伊井仁志（東大），沖田浩平（日大），塩崎聖治（理研），野
田茂穂（理研），石川顕一（東大），梁夫友（理研），山村直人（理研），片岡博
之（理研），石峯康浩(理研），K.Lee(理研），X. Gong（上海交通大）, 清水
和弥（東大），島本憲夫（東大），世良俊博（阪大），古田琢哉(理研）

Download Report