3D-RISMを中心とした生体機能解析：
理論と京への実装そして応用〜京スパコンで可能になること〜
吉田紀生
九州大学大学院理学研究院化学部門
平田文男（分子研・立命館大），丸山豊（大阪大），
Saree Phongphanphanee(分子研)，Daniel J. Sindhikara(立命館大),清田泰臣（北里大）
SS研HPCフォーラム2012

次世代ナノ統合シミュレーション
ソフトウェアの研究開発
・目的
開発中の次世代スーパーコンピュータを分子・物質の計算科学のシミュレーションで最大限に活用できるよう、計算科学理論・方法論の確立と、それらに基づくソフトウェアの開発を行う。電子・原子・分子レベルからの精緻な大規模計算に基づいて、ナノスケールの物質で発現する特有の現象・特性を解明し、予測することを可能とする。具体的に以下の３つのグランドチャレンジ課題を設置し、それぞれの分野の解明にブレイクスルーをもたらす計算科学の飛躍知の発見・発明に繋げる。研究開発拠点に産学連携体制を構築し、研究成果の産業界への展開を図る。
（ナノ統合拠点ホームページより）

2

中核アプリケーション
・6つのアプリケーション
・RISM/3D-RISM （代表：分子科学研究所 平田文男*）
液体の統計力学理論に基づき蛋白質などナノ分子の水和構造や水和自由エネルギーを計算するプログラム。 分布関数（相関関数）に関する積分方程式を解くために、フーリエ変換を多用する。 そのため、高速フーリエ変換（FFT）の速度が重要である。
（ナノ統合拠点ホームページより）*現・立命館大学

3

実在系の分子理論

理想系（単一の分子・孤立系）
量子化学，Schrödinger方程式

実在系(アボガドロ定数個の分子）10個，20個の積み上げではなく，パラダイムシフトが必要
統計力学，液体論

液体の統計力学理論：3D-RISM理論

液体の構造を記述する
統計力学理論

液体の積分方程式理論：
3D-RISM理論

5液体の構造液体の構造

気体の構造
特定の構造を持たない“確率論的構造 を持つが、一義的には決まらない

ある粒子からの平均的な存在確率で記述する

固体の構造
特定の構造を持ち、格子定数で決まる分布関数

液体の構造(分布関数)を求める
・実験的にはX線や中性子で求められる・分子シミュレーション
サンプリングにより求める(MD: ニュートン方程式，MC:乱数)

液体の積分方程式理論“解析的”に“直接”得られる

液体の積分方程式理論
MOZ

RISM

3D-RISM

モデル化なし分子としてそのまま扱う
６次元

相互作用点モデル
１次元

溶媒のみモデル化
溶質はモデル化なし
３次元

全自由度を扱える
角度に依存した分布など

１次元（サイト間距離）の分布関数

３次元の空間分布関数

厳密な方程式で理論的精度が高い
計算コストが高く，大きい分子は難しい

簡便で高速
適用範囲が広い
近似が大きく得られる情報は限定的

生体高分子の溶媒和を計算可能

3D-RISM理論の魅力
連続体モデル分子シミュレーションなど，他にも液体を扱う理論はあるが･･･溶媒を 分子 として扱う
水分子の酸素原子の分布確率

水素結合などの溶媒の分子としての特徴 を記述できる。
水分子の水素原子の分布確率
連続体モデルでは溶媒の分子的特徴は記述できない。

熱力学極限を満たす得られた分布は無限個の溶媒分子，無限時間の統計平均により得られるものに対応。
計算の初期条件に依存しない

分子シミュレーションでは溶媒の個数，サンプリングは有限。

熱力学量を(簡単に)計算できる

生体分子の分子認識への展開

分子認識〜最も重要な生体内素過程•酵素反応の初期過程•チャネルタンパクの選択性•抗原･抗体反応など系全体におよぶエントロピー変化を含む自由エネルギー変化を駆動力とする

統計的 , 確率論的 過程

液体の理論で分子認識？？
分子認識

高い親和性
すなわち高い分布確率

溶媒･リガンドの分布

従来の分子認識解析の難点(1)
内部に閉じ込められた分子• 分子シミュレーションでは分布は初期条件に依存

分子シミュレーション– あとからは，入れない– 事前に内部の分子分布の予備情報が必要(水やイオンの有無，個数）• 3D-RISM理論は初期条件に依存しない

3D-RISM– 完全な統計平均量としての分布が得られる– タンパク質の構造以外の予備情報は不要
16従来の分子認識解析の難点(2)“レアイベント”の解析
分子シミュレーション• 分子シミュレーションの時間スケールでは十分な統計量が得られない– 良く起こること は容易だが，起きにくいこと の解析は困難– 反応座標を設定しサンプルを限定
J. S. Hub, B. L. de Groot, Proc. Natl. Acad. Sci. U.S.A. 2008, 105, 1198.

– 人工的に圧力や電場を加えるなどして，加速して解析する
M. Ø . Jensen, E. Tajkhorshid, K. Schulten, Biophys. J. 2003, 85, 2884.

• 3D-RISMでは，完全な統計平均量としての分布が得られる– 外場や反応座標の設定は不要
3D-RISM17例：リゾチームの選択的イオン認識
N. Yoshida, et al., J. Am. Chem. Soc. Communication, 2006, 128, 12042
X-Ray

3D-RISM

リゾチームによるCa2+認識部位を再現
Kuroki, R.; Yutani, K.
J. Biol. Chem. 1998, 273, 34310野生型，ミュータントに対して，Na+, Ca2+, Li+の親和性を求め，実験を再現した
分子シミュレーションでは，結合部位の予備情報無しでは不可能

タンパク質の構造さえあれば，イオンの選択性，認識部位を非経験的・分子論的に決定できることを示した

3D-RISM理論に基づく
大きなリガンド分子の分子認識解析手法の開発と応用
清田泰臣，吉田紀生, 平田文男, JCTC, 2011

清田泰臣 助教（北里大学）

RISM理論による分子認識解析• RISM理論⇒溶媒分子の分布関数（＝確率密度）が得られる
溶媒としてリガンド分子を含めれば・・・
RISM理論によりリガンド分子の分布が求まる＝分子認識問題について応用できる• 従来法の問題– リガンドを溶媒として扱う– 溶媒（リガンドを含む）のバルク構造に関する計算が難しい
20

•従来

Solute-Solvent

LigandProtein
SoluteSolvent•本研究
Solute-Solute

N solventuuuuh(1, 2)  c
(1, 2)   c  cuvc (1, 3)h vuc (3, 2)d(3)cLigandProtein
SoluteSolute (ρ=0)
21•Procedure to obtain uu 3D-DFs.
•Calculate correlation function of solvent(water, ion, ...)
by vv 1D-RISM.
vv v v'v"hab(r)    acv  ccdvv  (r)   db
   v    acv  ccdvv ' (rcd )  hdb(rdb )
c,d

•

v'c,d

Calculate correlation funciton of protein by uv 3D-RISM.

v'vv'v
ha (r)     cb (r  )  ba( r   r )  v 'h ba( r   r ) dr v'•bCalculate correlation function of ligand by uv 1D-RISM.
uuhuvb (r)     cuvc (r)   cbv    v     cuvc (r)  hcbv' v (r)

 ,c

•

 ,cv'

Calculate correlation function between protein-ligand
by uu 3D-RISM.
uuh (r)    c (r  ) ( r   r )dr      ca (r  )v ' h va' u( r   r )dr v'a22分布関数からの認識部位決定
Information of ligand:Distribution function

Protein

Where?

How?

Location

Orientation

Two different types of score function

23

ApproachFor determining the location and the orientationst screening
search1BoxDistributionfunction

Orientation
search2nd screening

Location

Orientation

24

System

Phospholipase A2 - Aspirin

Phospholipase A2 and Aspirin form the complex for inducing
Anti-inflammatory effects through the specific binding site. [1]X-ray structure (PDB: 1OXR)

Protein and ions: Amber99
Water: TIP3P
ligand: Antechamber (GAFF)(a) 2-acetoxybenzoic acid(aspirin)(b)4-acetoxybenzoic acid(p-aspirin)

aspirin[1] R. K. Singh, A.S. Ethayathulla, T. Jabeen, S. Sharma, P. Kaur, and T. P. Singh., JDT. 13, 2005, 133.

25

• 3D‑DFs of Aspirin
Red: COOH sites
Blue: OCOCH3 sites
Yellow: Aromatic sites

aspirin

Asp49

Phe5

(gγ >2)26

• 3D‑DFs of p‑Aspirin
Red: COOH sites
Blue: OCOCH3 sites
Yellow: Aromatic sites

p-aspirin(gγ >2.5)

27

Box SearchAffinity based on the probability of ligand:1 NN1g (r, x)dr for Vbox  V protein (x)  Vligand  V  V
Vboxprotein ( x)f loc (x)    box
 0for Vbox  V protein (x)  Vligand
Vboxx・x

Vprotein28• Result of Box Search

aspirinCenter of molecule

His48
Asp49
Phe5floc(x) >1.45

29

• Result of Box Searchp-aspirin

Center of molecule

His48
Asp49
Phe5floc(x) >1.4530

• Result of Box Search
aspirin

Red:Overlap
p-aspirin

Blue:(Volume:138.1 Å3)(Volume:174.3 Å3)

His48
Asp49
Phe5floc(x) >1.45

31

Orientation SearchAffinity based on overlap:
Nˆ (Ω))
f ori (x, Ω)   g (x  l   RDistributions+Center of molecule
Euler angle

Protein32

• Result of Orientation Search
order

αβγ

overlap(a)11033050

2.33E+05(b)

240

260

40

2.02E+05(c)

360

22050

1.83E+05(a)(b)(c)

aspirin

Red: Predicted
Blue: X-ray

His48
Phe5

Asp49

33

• Result of Orientation Search
order

αβγ

overlap(a)11033050

2.33E+05(b)24026040

2.02E+05(c)

360

22050

1.83E+05(a)

(b)

(c)

aspirin

Red: Predicted
Blue: X-ray

His48
Phe5

Asp49

34

• Result of Orientation Searchorderαβγ

overlap(a)

1250

3030

9.26E+07(b)260300150

7.72E+07(c)

3300

250

110

5.90E+07(a)

(b)

(c)

p-aspirin

Red: Predicted
Blue: X-ray

His48
Asp49
Phe535

• Binding modeaspirin

p-aspirin

His48

His48
Asp49

Phe5NNOO
Ca2+OH

Phe5N

Asp49

NO
O
Ca2+OOOCH3

Tyr64

Tyr64

36

Summary•これまで困難だったサイズのリガンドの分子認識解析が可能になった•Phospholipase A2によるアスピリンの認識へ応用•課題•リガンドの構造揺らぎ•タンパク質の構造揺らぎ

37

３D-RISMを用いた溶媒配置予測アルゴリズムの開発
Daniel J. Sindhikara, Norio Yoshida, Fumio Hirata, J. Comput. Chem., 2012

Daniel J. Sindhikara
Jan 9th, 2012

http://goo.gl/uLohs

Placement algorithm

http://goo.gl/uLohs

39

Examples
HIV-1 Protease

F-ATP Synthase

http://goo.gl/uLohs

40

KNI-275 bound form of HIV-1 Protease• HIV-1 protease– Key step in replication• KNI-275 is a TS-analog with pM
activity• Combined x-ray, neutron
crystallography reveal six key watermolecules (PDB ID: 2ZYE )Adachi, M.; Ohhara, T.; Kurihara, K.; Tamada, T.; Honjo, E.; Okazaki, N.; Arai,S.; Shoyama,Kimura, K.; Matsumura, H.; Sugiyama, S.; Adachi, H.; Takano,K.; Mori, Y.; Hidaka, K.; Kimura, Hayashi, Y.; Kiso, Y.; Kuroki, R. P Natl Acad
Sci USA 2009, 106(12), 4641-4646.
http://goo.gl/uLohs

41

NKI275近くに見つかった水分子

http://goo.gl/uLohs

42

What didn’t crystallography find?

89

110

http://goo.gl/uLohs

130

43

Rotor Ring of F-Type Na+-ATPase• Part of F-ATPase machinery(animation)• Rotor driven by electrochemicalgradient• Crystallography reveals sodium and(after refinement) coordinatingwater atomic structure• 11-fold symmetry

Meier, T.; Polzer, P.; Diederichs, K.; Welte, W.; Dimroth, P. Science 2005, 308(5722), 659-662.
43.
Meier, T.; Krah, A.; Bond, P. J.; Pogoryelov, D.; Diederichs, K.; Faraldo-Gomez, J. D. J Mol Biol 2009, 391(2),498-507. (PDB:2wgm)
http://goo.gl/uLohs

44http://goo.gl/uLohs

45http://goo.gl/uLohs

46

Summary•

Discrete distributions can be coerced out of 3D-RISM distributions forboth water and ions•This technique matches quantitatively (about 0.5 Å error) with experiment
while providing more information•

Ion concentration required manual adjustment•

Can be used to add perspective to 3D-RISM, get insight into specific water
and/or ion interactions, provide an initial snapshot for MD simulations,give insights into drug design•

Info and tutorial on website, work done using free software (Ambertoolsand scipy) http://goo.gl/uLohs

http://goo.gl/uLohs

47

3D-RISMの京への実装

丸山豊博士（大阪大）

筑波大学 佐藤グループ
多田野寛人助教，高橋大介准教授，佐藤三久教授
48

3D-RISMのアルゴリズム
溶質—溶媒相互作用
ポテンシャルClosure方程式

3D-FFT

3D-FFT

3D-RISM方程式

3D-DF
各種物理量
49

（従来型）3D-FFTの並列化の問題点• 三次元空間の1方向（例えばｚ軸）を分散して並列化• 最大で1軸分の離散点数までの並列化しか出来ない– 3D-RISMでは次世代スパコンで20483のグリッド点を想定– 従来型では最大で2048並列しかできない• 軸交換時に全対全通信を必要とする– 通信律速50

Volumetric並列三次元ＦＦＴ• 近年、筑波大学高橋大介准教授らによって開発された並列化手法• 三次元空間の2軸を分散して並列化– 20483なら20482並列まで可能– 軸交換時の全対全通信の通信量を削減可能• 3D-RISMに実装することで高効率並列計算が可能になる

51

従来の3D-FFTの並列化• Z軸方向へ一次元ブロック分割

zyxx軸方向にFFTy軸方向にFFTz軸方向にFFT

全対全通信による軸交換

全対全通信による軸交換
52

Volumetric並列3D-FFT• y,z軸方向へ二次元ブロック分割

zyx
x軸方向にFFTxy面での軸交換

y軸方向にFFTz軸方向にFFTyz面での軸交換

xz面での軸交換
53

通信の半減←通信←通信

実空間とフーリエ空間で分割のレイアウトが変わる事にプログラム側で対応する事によりFFT部分の通信量を半減させた。

54

ナノ統合中核アプリとしての目標•

1万ノード程度でスケール結果– 10243グリッド
2048ノードまでテスト済み（右図）– 20483グリッド
8192ノードでスケール4096ノード比で88%の効率高速化率（1024コア基準）

3D-RISMの京への実装
181614

3D-RISM全体

1210

FFTのみ

8642005000

10000

15000

20000

core

55

京で可能になること・やりたいこと• ドラッグスクリーニング• 各種連成計算– 構造サンプリング– 電子状態計算

56

ドラッグスクリーニング
標的たんぱく質とリガンドの認識を予測・解析• どこに？ 認識部位探索• なにが？ リガンドスクリーニング
3D-RISMのアドバンテージ• たんぱく質を溶質として，リガンドを溶媒（あるいは第2の溶質）として計算することで上記2つを同時に解析可能
京での展開• 1つのたんぱく質（構造）・１つのリガンドに対して1回の3D-RISM収束計算が必要• 3D-RISMを高速に解くことで多くの構造・リガンドへの解析が可能
課題• たんぱく質の構造揺らぎ• リガンドの構造揺らぎ

57

連成計算
各種連成計算でできること• たんぱく質構造サンプリング– ３D-RISMとMD, MCとの連成計算•化学反応解析

3D-RISM– 3D-RISMとQM/MMの連成計算

3D-RISMのアドバンテージ• 溶媒の温度依存性，密度依存性などが簡単に扱える• 共溶媒効果や希薄電解質を扱える

QM/MM

京での展開
• セルロース分解酵素反応解析• リグニン構造の温度依存性，溶媒依存性課題• 高効率連成計算
58

今後の展開• 理論的課題– 溶質分子の構造揺らぎ– 溶媒分子の内部構造揺らぎ• プログラムの課題– 連成計算の高度化– ポスト/プリプロセスの充実• インプット作成支援• アウトプット解析支援59

エクサスケールコンピューティングに向けて• エクサスケールコンピュータ– Many-cores?
– GPU?

• 3D-RISM for GPGPU– 丸山豊博士（大阪大）
Modified Anderson Method for Accelerating 3D-RISM
Calculations Using GPU, Y. Maruyama and F. Hirata, JCTC, 2012

• 大並列化にはmulti-GPU用3D-FFTが必要– 額田彰博士（東工大）
High Performance 3-D FFT using multiple CUDA GPUs, A.
Nukada, Y. Maruyama, S. Matsuoka, GPGPU-5 proceedings, 2012

• エクサに対応するなら・・・– 3D-RISMの理論と計算機科学の境界領域を押さえられる人材– 中〜長期的なプロジェクト推進が必要
60

以上です