+'"&%+'&  "&*-%$) #"(  !,2019年12月11日
横浜市立大学 生命医科学研究科
理化学研究所 MIH
池口 満徳

-• –  6:<>A6=?-#

6:<>A6=?*.

– 795?$%'1$• %/

 6:<>A6=?-@

– +-0&,



6:<>A6=?3(2-'B•  6:<>A6=?*(AI)-!
– LINC5?8A64;-– -

")-AI- =DFJM=GL1 9I<H?Eü MD-2#3212 fs O10-15s) 1ü 1>A@B1#0&&7.,7.#
0#86,400>A@B=86 ps+&-(/$Nü 1>A@B1DH--(7.#0#86,400,000>A@B=86 ns0571

1

ü *1-41C:;K1MD012)6$&&7(8%02#"/'!NHistory of Molecular Dynamics Simulations ofProteins
Ca2+-ATPase (more than 270,000 atoms including water)2005 Sugita, Ikeguchi, Kidera, Toyoshima, JACS, 127,6150

BPTI (3000 atoms)
Nature, 267, 585, 1977
in vacuum, 10 ps

Lysozyme(20,000 atoms)
1990s
in water,several ns

Flagship Supercomputer “K” スパコンの世界順位
13
57911131517192120196201811

2018620171120176

2016112016620151120156

20141120146

20131120136

20121120126

20111120116

Ranking of K-Computer  Top500

 K-computer
663,552 cores
1.327 PB memory
11.28 PFLOPSRpeak
10.51 PFLOPSRmax

Shutdown in Aug. 2019

Fugaku23:JQS=#9 5<
#;" 4'1" 796T
~  knz
 |s < n0 Y  s  w~

~ q>hjVdX}}

< 3Nn: S

3N &[n < }~n*aAEaA3NS-mop4Pn`'_*asw~7,!

A 3N &2"5(%r{7:D+?

wUaC>,

~ fgknMDd$
:D+?e2"5QM/MM Free Energy7.

/\-wUa^306]

IRxd?H?

KM ;1A
:D+?

O@uvTF7
2"5ts~1A
H?

Q9
 |s 

~ MDq3Wlijc)APTFBn0G
3N &[DB=LD#dJ$P |s]PbZx,8?
1A  w
,!
~yx,?

AL&*,URVW

AL&*-U R VW

AL&*.U R (VWMCI7:0/+ 8>P@ODN

#%$:)

#KF?HSE! BCGN:)

RIKEN Sugita Group: GENESIS

https://postk-web.r-ccs.riken.jp/perf.html

   

$%




* *!(+472*1+50/86)'&Folding, Ligand Binding,Conformational Changes

MD,+60

$%

5
2
11

MD1
MP-CAFEE

aggregation #(-36.8

100nsμsms"汎用スパコンでの動的機能解析
数十μs/day計算可能な専用計算機（ANTON2等）と違って、我々が使用可能なのは汎用スパコンである
一つの計算では、汎用計算機だと、現状からの飛躍的向上は、非常に困難である
しかし、汎用スパコンでは、ノード or CPUコアは沢山ある

レプリカを沢山用意して、統計的に計算を行うことで、実質的に、長時間ダイナミクス情報を得る

汎用スパコンでの動的機能解析
レプリカ系のMD計算
方法① 熱力学的方法
自由エネルギー分布を得る（取りやすい複数の状態と頻度比を得る）
e.g. レプリカ交換、アンブレラサンプリング等
方法② 速度論的方法
キネティックス（速度論）を得る（取りやすい複数の状態間遷移を得る）
e.g. マルコフ状態遷移モデル（MSM）
Weighted Ensemble

Hierarchy of Protein Dynamics
10-15
femto
bondvibration

10-12
pico

10-9
nano

methylrotation

10-6
micro

loop
motion

10-3
milli

100
seconddomainmotion103
hourslow dynamics

protein folding

side chain rotamers

106
day

drug binding
Conventional MD Simulations
Enhanced Sampling, Ensemble Simulations, etc.

MD;&;4F7<WX•)@CE5G%F– $(@=73E26/=7302Y• +@-G%F
– $(@=73*8:+8:0F2Y• G%F– RIUN<@A/=7302Y• – !

,GF"HOSJV/HULSPV/HUMTPVGF• RIUN!

QVKGF– RIUN@>!

•)

QVK

G%F

– =@D1?/$(A)

9F@2Y• .!<#9F)GF– NMRBSAXS>=@.!<G'MD•   – 

MD68:<>6;=.: RMSD, RMSF
Root Mean Square Deviation (RMSD)ü %1.*34ü .
ü ! .

42.-0ü 795=)+$2?
."!# / &',2
795="!# / ($RMSD % =∑( )( (%) − )-./,(21MD.0246.35&8RMSD
7RMSD from crystal structure)

RMSD (Å)

Root Mean Square Deviation (RMSD)

RMSD =

(2xref i )i xiNü )&"+,
ü MD&  ,*&%(ü /1-5!#*8
&RMSD' $*/1-5RMSD'
domain A Time (ns)domain B

RCA

Oroguchi, MI et al., Biophys. J., 96, 2808 (2009)

RCB

MDORTWZOUXG.$: RMSD, RMSF
Root Mean Square Deviation (RMSD)
ü '0;JGBLM-%ü 3G4
ü 65G

M,KGFI(ü PSNX1 AC:K^
G7685H#/*=?EKPSNX7685H@:RMSD % =∑( )( (%) − )-./,(21RMSF 4 =∑5 )( (%) − )67.,(81Root Mean Square Fluctuation (RMSF)
ü 0;JGBL[ü)Y")&\G4ü J:D:K+9:+Mü ]!02HJ><@:ü VZQ2HJ><=:Conformational Changes of β subunit
upon ATP bindingγγ

βEβE

βTPβTPProtein Fluctuation (RMSF)
Empty; Open)ATP bound; Closed)
ADP bound; Closed)

ü Fluctuation of C-terminal domain: βE>βTP>βDPü Despite fairly similar conformations of βTP and βDP,
fluctuations of C-terminal domain are different.

Y. Ito & M. Ikeguchi, J. Comp. Chem. 31, 2175 (2010)

βE

βTP

MD
•

#%& !"$ #– "'• –

$ # #'MD!#%"$
Correlation Matrix &ü ü ü ':, PCA
Cij =Δri ⋅ Δr j
Δri2 Δr 2j $

Cij = 1 Cij = 0 Cij = −1

Correlated Motion in ProteinCij = 1 Cij = 0 Cij = −1

Y. Ito & M. Ikeguchi, J. Comp. Chem. 31, 2175 (2010)

MD9;=?A9>@2#: &!, PCA
Correlation Matrix B&!C
ü $&2ü 1üCij =Δri ⋅ Δr j
Δri2 Δr 2j(-(57:<8@*6)5Principal Component Analysis (PCA)
B#Cü '$

ü $71#

Cij = 1 Cij = 0 Cij = −17ü 7 3,5+/.
$ 24)0%*6)5!7"∆"# ∆"$% = 1 ⋯ 3*

MD"#$&("%') PCA
 Principal Component Analysis (PCA)

PCA

ü PCA ! ü  !
 Computational Biochemistry and Biophysics, Becker et al. eds., 2001.

PCAF F1-ATPase β;@B?<=DE
ü Top mode 4$7*3 30 % :68
ü Top ten modes 4$7*3 70 % :68

ü 43AC>4 ATP3,4$# 2 )1 :.($7*/580$("79-235$7'/'8(&8%Y. Ito & M. Ikeguchi, Chem. Phys. Lett. 490, 80 (2010)+$0'!:

Principal Component Analysis: Mode 1-4üRelaxation from interactions withadjacent subunits

Mode 1: magenta
Mode 2: yellow
Mode 3: cyan

Mode 4: green
βE→βDP: red
Y. Ito & M. Ikeguchi, Chem. Phys. Lett. 490, 80 (2010)üIntrinsic flexibility is well correlatedwith structural transition

MD!!("./•  $&')(–  $# '#0• $ )(–  $#  (0• )(– ,*-+"$  %#0Interaction between ligand and protein

Interaction between CYP1A2 and EPA

Isobe, YI, MI, Arita, Sci. Rep. 8: 7951 (2018)

MD*

*$1'+@A• -/0%21– -,'# 0"& ,'#!"B• -21– -,'#()()!1"B• 21– ;4>8+-. ,'# !"B•  2 1– 39<5? 3>6<:? 3>7=:?2 1

MD*-/38:4;


<A=

i=0
λ=0




2#$%0&!
5769+

i=1
λ=0.01

<A=.<B=**+,

i=2
λ=0.02


<B=

i=N
λ=1

1/' +38:4;

2)3, − 2)38:4; ∆" = −%& ' ( ln exp −Free Energy Perturbation%& ')),62  Thermodynamic Integration=497
5 67 8"+*,Bennett Acceptance Ratio (BAR))(
ΔGbindΔGcomplex=ΔGligand=ΔGcomplex

ΔGligand

MDFEP, TI, BAR

ΔGbind = ΔGcomplex − ΔGligand

8 1 6(MD MD



816(MD

Temperature;;;'/7)<9<(73.24
-5,1&"':-90<!;;; ;;;'/7)<

High

Low





816( MD #$
+:16:*%Temperature

High

Low

Replica 3

Replica 2

Replica 2

Replica 2

Replica 3

Replica 0

Replica 1

Replica 0

Replica 3

Replica 0

Replica 1

Replica 1

…
…
…
…

.)-%MD(REMD)!"

Torsion angle [θ]MDTorsion angle []: Oseltamivir &*+.1&,0#REMD# .)-%MD !"   
$0'/(1


Time [ns]ΔSClausius [cal/mol/K]

Time [ns]18.4

MD: REMD: 18
17.6
17.2 0

Hikiri, MI, et al. J Chem Theory Comput. 2016 12: 5990.

50

100

150

Time [ns]20038:<>+Folding/Unfolding+=9<1

MD(REMD)
UnfoldFold

38:<>(10)

ΔSc

REMD
(255K-400K)
<>
16000

Unfold :Fold :Fold

Unfold

12000
8000
4000

Simulation time [ns]00 1 2 3 4 5 6 7#-+RMSD (Å)

Fold/Unfold

100
PFold(T) [%]*"$.#-+RMSD (Å)

7;4027<?6?5(T=305K)

80

Tm(REMD)

60
40Tm()200 250 280 310 340 370 400

Temperature [K]• REMD(38:<>*&!' +Folding/Unfolding+/• Folding/Unfolding+/ .%)* • REMD(+(Tm=305K),(315K)*!

 Folding  ΔSc(T)Unfold

Fold

2.0
ΔGFoldUnfold(T) [kcal/mol]

ΔSc-10Fold
ΔSUnfold:

ΔS

:

0.5
0.0-0.5

Temperature [K]0.0-20-30-40285

1.0-1.0 250 280 310 340 370 400

ΔSFoldUnfold(T) [cal/mol/K]

ΔSc(T) [cal/mol/K]

01.5295

305

315

Temperature [K]Hikiri, MI, et al. J Chem Theory Comput. 2016 12: 5990.

325

335

-5.0

-10.0
-15.0
-20.0
-25.0

-30.0 250 280 310 340 370 400

Temperature [K]MD..(6+/GH• 245)76– 20+'5&*$0+'%&I• !2#76– 20+'

,-!,-%6&I• 76– B9E>/23$0+'• %&I"76– 8?C:F$8E<C@F$8E=D@F76• B9E>AF;76– B9E>2

1AF;

(Generalized) Replica Exchange with Solute Tempering (REST / gREST)
REST: Liu et al. PNAS, 102, 13749 (2005).
gREST: M.Kamiya and Y.Sugita, JCP, 149, 072304 (2018).
highEnergy

450K
300K

Temperature(T)

600K

Conventional REST

replica3

replica3

replica1

replica1

replica2

replica2

replica2

replica1

replica3

replica2

replica1

replica3

replica1

replica2

replica2

replica3

replica3

replica1

MD(600K)

MD(450K)
MD(300K)
lowtime(t)Solute region

Solute region
ChargeLJ
BondAngle
Dihedral
CMAP

gREST

ChargeLJ
BondAngle
Dihedral
CMAP4;1
D

938. N-M

938. 3-1
938. 2-1
938. 1-1

ü 938.)+$%'8/:2
",( -ü 2&*ü 4;1- 

DN
 TM

DN-1 TM

DN-1
 TM

D3
 T1
D2
 T1
D1
 T1

D3
 T1
D1
 T1
D2
 T1

D1
 T1
D3
 T1
D2
 T1

-!(

938.*MD&→ #0569;07:MD22+<.3MN• !79;,=<– 74.*;)-&4.*')O• #7%=<– 74.*"/1#/1'<)O• =<– H?KD378&4.*')O• – $=<>EI@L&>KBIFL&>KCJFL=<• H?KD– H?KD7• !

GLA=<5GLA

=<– 47:(6&8!

0<7)O

Method: umbrella sampling simulations

Umbrella Sampling MDµs - mseq. MDeq. MDFree energy

Basin aroundopen structure

Basin aroundclosed structure
Conformations

Method: umbrella sampling simulations
0.5-13)
,5/25+((%MD#%!
$
4/2* MD"%&Replica Exchange UmbrellaSampling (REUS)Free energy)'

6Reaction Coordinate7
Conformations

Method: umbrella sampling simulations

/-2,01) $&(#% ))$!$*.1+3) 'Free energy)"'

4Reaction Coordinate5
ConformationsFree Energy Simulations
1. Initial Path
Nudged Elastic Band
openclose

2. Umbrella Sampling Simulations
For each intermediates along path, simulations
with restraint,

w j = K rmsd (ΔDrmsd − ΔDmin ) 2

ΔDrmsd = rmsd(X, X open ) − rmsd(X, X closed )was carried out.
Restraint potential is applied to both main chains
and side chains.
B. Roux JACS (2005) C.L. Brooks PNAS (2007)

3. WHAM
Weighted Histogram
Analysis Method was used
for removing restraints andcalculating free energy
profiles along structuraltransition between open
and closed conformations.

Free Energy Profiles for Open-Close Transition
blue: open
yellow: closed
green: ATPopen

close

opencloseü Meta-stable state is found:ATP-bound open conformation.
encounter complex

Ito, Oroguchi, Ikeguchi, JACS, 133, 3372 (2011)  transition: H-bond in P-loop

OpenH-bond partner of Asp256 is changedfrom Lys162 to Thr163 (P-loop).
Mutation of Asp256, Lys162, or Thr163 results inremaining open conformation even with ATP bound.

Closedマルコフ状態遷移モデルによるダイナミクス解析
A lot of short MDs

Definition of‘State’
Simulation time

Simulation time

…
…

…
…

…

…

…
…

…
…

Markov state modelTransition matrix
Population
Network

U#&*MDINORUIPS

All atom (MM)Slow motion
High precision

Coupling force Coarse-grained (CG)-Fast motion
Low precision

TMMBCG8*245"9GSMQSF
TRMQD:>@MM61;BAV %9)B3AJSLE';+ B,$;$);<.?KHCS3A; !

(MDClosed)7;
(@;$Open)7;
/9=?0

MM-CG
*MD

K. Moritsugu, T. Terada, A. Kidera, 2017?

)0%+,$/&"-0(*!-.#'MSES   
weighted ensemble  Moritsugu et al., JCP 2018

)0%+,$/&"-0(*!-.#'

MSES   
weighted ensemble  '&)8@4:;$",3?51<@7)!9/<=26MSES(+, !%weighted ensemble)*-#3-basin kinetic model:

3?5@6

open


WE2 0.011 ms-1

WE1 0.047 ms-1

PC 2

energy barrier

closed0.0062 ms-1
DG+- = -kBT ln (k+ / k- )0.30 ms-1

3?5@6

4.6 kJ/mol
closed/helicalPC 1

3?5@6(+,conformational switch

open/
helical

3.1 kJ/mol
open/
coiledk = 1.1 ms-1=>0@9@)!

.



MD;&;4F7<WX•)@CE5G%F– $(@=73E26/=7302Y• +@-G%F
– $(@=73*8:+8:0F2Y• G%F– RIUN<@A/=7302Y• – !

,GF"HOSJV/HULSPV/HUMTPVGF• RIUN!

QVKGF– RIUN@>!

•)

QVK

G%F

– =@D1?/$(A)

9F@2Y• .!<#9F)GF– NMRBSAXS>=@.!<G'89P#IT#',7
B'Q$&
人工知能（AI）技術の活用
ライフ系企業

IT企業
キックオフミーティングの様子
アカデミック

0@O8
理化学研究所
医薬プロセス最適化プラットフォーム推進グループ
医薬基盤・健康・栄養研究所

10%LO-.IO=
1. N
2.  N
3. 6OA4.  1BEJO1FM
5. D9/CN NADMET
6. 8GM3JO1F:INH0O5
7. ;)+K2-3NNK?6(-38. "N ! ND7(,IN'<*'O49.  >O3NNLP
10. AI 

池口満徳
WG4代表

本間光貴
WG5代表

LINC

LINCパンフレットより

ターゲット探索・検証

スクリーニング
ヒット to リード

疾病関連遺伝子

創薬標的タンパク質
X線結晶解析
NMR
PJ11

電子顕微鏡

化合物ライブラリ
PJ12

in silico

HTS

スクリーニング

疾病動物モデル
薬物動態
毒性（ADMET）
PJ18

ヒット化合物
PJ13
PJ14

タンパク質立体構造
WG4

前臨床試験

リード最適化

薬効
立体構造

候補化合物
化合物
最適化

PJ15

PJ16
PJ17

メドケム
WG5

立体構造決定法
X線結晶構造解析

NMR(核磁気共鳴)

電子顕微鏡

実験データ
人工知能
コンピュータ
A. I.

タンパク質立体構造

PJ11

AI化することで、熟練した技能が必要
これまで、熟練者のみが可能であった高度な解析が可能に

PJ12：AIによるドッキング計算高度化
これまでの問題点：精度の低さ（1~2%のヒット率）
黄：実験、緑：計算

正解ポーズ、活性化合物群

黄：実験、緑：計算

不正解ポーズ、非活性化合物群

学習
AIは、画像認識のように、正しい結合を間違った結合から見出すことができるか？
J. Chem. Inf. Model. 2017, 57, 942–957.

画像認識

結合計算

色（RGB）

原子種

２次元対象（猫等）

３次元

結合ポーズ

PJ13：分子動力学計算におけるAI活用
これまでの問題点
時間軸を含む大量のデータ（静止画像 → 動画）
熟練した技術者が長時間かけて解析

AIによって、重要な箇所を発見

膜タンパク質アクアポリンのMD
全原子数～16万数千万～数十億ステップの計算

動的なタンパク質の制御による創薬のブレークスルーを目指すPJ14：AIを用いた高精度分子力場
分子力場

分子軌道

MDでよく使われる分子運動◯

反応解析でよく使われる電子状態

分子運動

☓

電子状態を加味した
高精度力場を高速計算

☓

人工知能
A.I.
分子運動電子状態◯◯電子状態◯

電子状態を学習まとめ• 分子動力学シミュレーションの新展開‒ スパコン「京」から「富岳」へ‒ アンサンブル シミュレーションの有効性• 分子動力学シミュレーションの解析・応用‒ 分子動力学シミュレーションの活用例• 分子シミュレーションと人工知能(AI)の連携‒ LINCコンソーシアム：WG4‒ 分子子ミューレヨン分野でのAIの活用