カスケード選択型分子動力学法によるタンパク質機能の動的秩序解析

Abstract

生体機能はタンパク質の構造遷移と深く関わっており，その動的秩序と機能の相関を明らかにする必要がある．しかしながら，従来の分子動力学シミュレーションで到達可能な時間スケールは生体機能の発現する時間スケールよりも短く，生体機能に関わる構造変化をサンプリングする研究手法が求められている．この総説では我々が提唱するカスケード選択型分子動力学法の概要を示すとともに，生物学的に重要な現象に対する応用例を紹介する．

1 はじめに

タンパク質は極めて複雑な構造と精緻な機能をもつ高分子化合物であり, 分子認識, 情報伝達, 酵素反応など, 生体内でおこる様々な生命現象の根幹をなす. タンパク質は化学的性質の異なる20種類のαアミノ酸から構成されているが, アミノ酸残基間に生じる種々の弱い相互作用 (水素結合, ファンデルワールス相互作用, 静電相互作用) のため, 1つのタンパク質であったとしても複数の (準) 安定構造が存在する. 生体内においては, 各々のタンパク質はある特定の構造 (Native状態) をとり, そのタンパク質特異的な機能を発現しているものと考えられている. このように立体構造と機能の間には大きな相関 (構造-活性相関) があることが期待されていることから, これまでX線回折実験や核磁気共鳴法(NMR) などの実験的手法により, 数多くのタンパク質の立体構造が明らかにされてきた.現在, 10万件にも及ぶ構造が Protein Data Bank (PDB) に登録されている. しかし, 実験的手法は膨大な時間と労力を要し, また, 構造が決定できないタンパク質も多く存在することから, 理論計算に寄せられる期待は高い. 一方で, タンパク質の機能発現に重要な構造変化は, 通常の分子動力学シミュレーション (MD) で追跡可能な時間より長時間の確率過程において観測される「レアイベント」である. 従って, タンパク質フォールディングやドメイン運動などの大規模な構造変化を追跡するためには, 長時間ダイナミクス [1,2], もしくは, マルチカノニカル法 [3,4] やレプリカ交換法 [5] などの効率的なサンプリング手法が必要であった. 前者はAntonなどのスーパーコンピュータ [6] でのみ達成されており, 後者は構造サンプリングを達成するために手法ごとにある種のノウハウの修得が必須であり,誰しもが簡単に計算を実行できる訳ではない.

我々は, より簡便かつ高速なサンプリング手法として, (1) 構造変化を誘起する可能性が高い初期構造選択と, (2) 短時間MDに基づく初期構造の構造リサンプリングから成る効率的構造サンプリング手法を開発した. 本手法で鍵となるのは, 構造変化を誘起する可能性が高い初期構造を適切に選択することである. 遷移経路が得られた後は, その経路に沿って精密な自由エネルギー解析を行えば良いことから, 従来の自由エネルギー解析法よりも効率的である. 本総説では, タンパク質の構造変化を抽出する手法として, 4つの計算スキームとその具体的な適用例を概観する. それぞれの開発手法の詳細については, レビュー論文 [7] を参照されたい.

2 理論的背景

2.1 タンパク質の分子動力学シミュレーション

分子動力学シミュレーション(MD)とは, 各原子を質点と見なしてNewtonの運動方程式   

$m{\ddot{R}}_i=F_i=-{\nabla }_{i}V\left(\left\{R_j\right\}\right)$

の数値解を求めることで, 原子位置やその運動について解析するための方法である. ここで, R_iはi番目の原子の座標, m_iはその質量, F_iはその原子に掛かる力であり, ポテンシャルV({R_j})の座標微分から計算される. 特に, タンパク質の古典MDシミュレーションでは, ポテンシャルはタンパク質の安定構造を表現するためにモデル化されている. 例えば, Amber分子力場は   

$\begin{array}{c}{V}^{\text{MM}}=\sum _{\text{bonds}}{K}_{ij}^R{\left(R_{ij}-R_{ij}^{eq}\right)}^2+\sum _{\text{angle}}{K}_{ijk}^{\theta }{\left({\theta }_{ijk}-{\theta }_{ijk}^{eq}\right)}^2\\ +\frac{1}{2}\sum _{\text{torsion}}{K}_{ijll}^{\phi }\left[1+cos\left(n{\phi }_{ijk}-{\phi }_{ijk}^{eq}\right)\right]\\ +\sum _{i>j}\left[\frac{A_{ij}}{R_{ij}^{12}}-\frac{B_{ij}}{R_{ij}^6}\right]+\sum _{i>j}\frac{q_i{q}_j}{\epsilon R_{ij}}\end{array}$

のように与えられ, 第1項, 第2項, 第3項はそれぞれ, 分子内における2原子間の結合距離R_ij, 3原子間の結合角度θ_ijk, 4原子からなる二面角ϕ_ijklの関数として表される (Figure 1を参照). ここで, $K_A^X$ならびに$X_A^{\text{eq}}\left(X=R,\theta \text{, }\phi ,A=i,j,\cdots \right)$は, アミノ酸残基ごとに決められるパラメータである. 一方, 第4項, 第5項はそれぞれ, van der Waals (vdW)相互作用, 原子間の静電相互作用の寄与を表し, いずれも2原子間の距離の関数である. また, A_ijやB_ijはvdW半径から求められる原子対ごとのパラメータであり, q_iは原子の電荷を表す.

Figure 1

 . Definitions of bond distance, bond angle, and dihedral angle in the AMBER force field

Newtonの運動方程式は2階微分方程式なので, 数値的に解くためには, 時刻t = 0における位置r_i(0), および速度v_i(0)などの初期条件を与える必要がある. 初期位置については, PDBに登録されているX線結晶構造解析などの実験手法で得られた構造からモデリングすることが多い. 一方, 初期速度v_i(0)を与えるには, 温度(T)と運動エネルギー(E_kin)の関係式$E_{\text{kin}}=3N{k}_{\text{B}}T/2$を利用し, 設定したい温度TでのMaxwell-Boltzmann (MB) 分布   

$f_{\text{MB}}\left(\left\{v_i\right\}\right)=\prod _{i=1}^N\left(\frac{m_i}{2\pi k_{\text{B}}T}\right)exp\left(-\frac{m_i{v}_i^2}{2k_{\text{B}}T}\right)$

に従うように, 乱数を用いて各原子に初期速度v_i(0)を割り振る. ここで, k_Bはボルツマン定数である. しかし, これだけでは系の温度は設定したい温度にならないため, 系の温度をある温度Tに制御するための様々な計算手法(Berendsen熱浴 [8], Langevin熱浴 [9], Nose-Hoover熱浴 [10,11] など)が考案されている.

平衡状態に到達し, 物理量を計算(測定)するためには, 無限時間となるくらいの十分長い時間のMDが必要である. 例えば, 前述のタンパク質フォールディングに関しては, ミリ秒から秒の時間スケールで折りたたむことが知られていることから, 単一の構造から出発して平衡状態に達成するためには, 更に長時間のMDが必要となる. 一方で, 異なる初期座標•速度から出発した複数のMDの結果を合算して平衡状態の分布を計算することも, 原理的には可能である. しかしながら, もし初期状態が類似していると, ほぼ同じ統計アンサンブルを重複して合算することも考えられる. また, ある温度で許される同程度のエネルギー状態の安定構造がいくつ存在しているかは系に依存しているため, 簡単な例を除いて事前に知ることができない. 以上のことから, 初期条件の選び方には細心の注意が必要であり, 正しい統計アンサンブルを得るためには, ある温度で許される構造を全て網羅することが肝要である.

2.2 効率的な構造サンプリング手法

我々の提案する方法は, 極めてシンプルである. Figure 2 にそのアルゴリズムのフローチャートを示す. 以下に, 構造サンプリング法の計算手順を簡単に解説する.

Figure 2

 . A flowchart of our conformational sampling scheme

1. 選ばれた複数の初期状態から, 比較的短時間 (100 psなど) のカノニカルアンサンブルでのMDを実行.

2. 予め決められた選択則 (selection rule) に基づき, 各々のMDトラジェクトリの構造をランキング.

3. 高いランクのスナップショットが優先的に選択されるよう, 確率的に次のサイクルの初期構造 (seed)を選択.

4. MB分布に従い, ターゲット温度の初期速度を再生成し, 選択した初期構造から短時間MDをリスタート.

5. (1) ~ (4) のサイクルを繰り返し, ある反応座標に粗視化した分布がこれ以上変化しない, あるいは, 望みの構造に十分近い構造が得られた時点でサイクルを終了.

6. 全サイクルをつなぎ合わせた反応経路 (reactive trajectory) に対して, Umbrella Sampling method (U.S.法) [12,13], およびWeighted Histogram Analysis Method (WHAM) [14,15] を適用し, 望みの反応座標に射影した自由エネルギー面を描画.

前述のように, 本手法では「どのようにして遷移確率が高い初期構造を選択するのか?」ということが, 効率的な構造探索を実現する上で重要な鍵となる. 我々はこれまで, Table 1 に示すような複数の構造サンプリング法を提案してきた. これらの方法論は, 対象となる系に合わせて, 計算手順 (2) の選択則のみが異なり, 骨格となるアルゴリズムは同一である. 計算手順 (6) 以外は, 構造と初期速度を何度も再生成させて短時間MDを繰り返すだけなので, あらゆる市販およびフリーソースのMDコードを修正せずに適用できる. いずれの方法においても, 同時に複数の短時間MDを実行すれば, 多くの配置サンプリングが達成できることから, 分散処理が容易である. また, ひとたびreactive trajectoryが得られれば, 計算手順 (6) のU.S.法において分散処理可能であることから, 本手法は並列計算機環境に適した極めて実用的な計算手法である. 以下では, 本総説でとりあげる各計算手法を概説する.

Table 1. Categorization of conformational sampling method

(A) Fluctuation Flooding Method (FFM)

機能発現に関係しているタンパク質のダイナミクスは異方的である場合が多く, 非調和性の高いある特定の振動モードが支配的である. 我々は, 主成分分析により抽出される異方性の高い構造揺らぎをもつタンパク質構造が高い確率で構造遷移するものと仮定し, これらを構造リサンプリングしていくことで効率的に構造遷移を誘起するFluctuation Flooding Method (FFM) を開発した [16]. Figure 3 (a) にFFMの基本概念を示す. 概念図から分かるように, 主成分分析により抽出したj番目の主成分モード方向に射影した主成分座標 (PC_j) の最大および最小値は, 現在探索している中で最も準安定構造から近い, もしくは離れていることが予想されるので, $\left({\text{PC}}_j^{min},{\text{PC}}_j^{max}\right)$近傍に存在する構造を次のサイクルの初期構造として選択し, 短時間MDをリスタートするサイクルを繰り返すことにより, 構造遷移を効率的に誘起することが可能となる.

Figure 3

(a) Concept of FFM

(b) Concept of OFLOOD method

(c) Concept of PaCS-MD

(B) Outlier FLOODing (OFLOOD) Method (OFLOOD法)

タンパク質は多自由度複雑系であるが, 高次元の反応座標空間上に複数の密度分布の高い準安定状態(クラスタ)が存在し, 大きな構造変化の際に異なるクラスタ間を遷移する. 密度分布の低い疎(スパース)な構造は, 複数のクラスタの淵(時には2つのクラスターの間)に存在することから, 密度分布の低い状態を集中的に選択することで構造変化を誘起する. 情報科学的に, クラスタに属さないスパースな分布は「はずれ値」あるいは "Outlier" と呼ばれ, 階層的クラスタリング手法であるFlexDice [17] などを用いて検出できる. 我々は, 反応座標空間における「はずれ値」に対応するタンパク質構造を集中的に構造リサンプリングすることで, 効率的に構造遷移を促進する Outlier FLOODing Method (OFLOOD法) [18] を開発した. OFLOOD法の概念図を Figure 3 (b) に示す. 概念図において, 反応座標空間に分布するタンパク質の状態分布 (黒) は階層的にクラスタリングされ, それぞれの階層においてクラスタに属さないはずれ値 (赤) が検出される. クラスタリング終了後, はずれ値に対応するタンパク質構造を次のサイクルの初期構造として選択し, 構造リサンプリングを繰り返すことで, サイクルを増すごとに構造分布が広がっていく. OFLOOD法の拡張形として, はずれ度合いの強さを考慮し, より出現頻度の小さい状態を優先的に初期構造に選択し, 構造探索していくことで, 更なる探索効率の向上に成功している. 詳細に関しては, 本論文 [19,20,21] を参照されたい.

(C) Parallel Cascade Selection Molecular Dynamics(PaCS-MD)

PaCS-MD [22] は, 短時間MDより得られたトラジェクトリとターゲット構造との類似度 (例えば, 構造間の平均自乗距離, root-mean square deviation (RMSD)) に基づいてランキングし, その上位を初期構造として短時間MDをリスタートする (このプロセスを1サイクルとする) . PaCS-MDの概念図を Figure 3 (c) に示す. PaCS-MDは, この一連のサイクルを繰り返すことにより, 外部摂動を用いることなく構造遷移を誘起することが可能である. FFMとOFLOOD法とは異なり, PaCS-MDはターゲットが必要であり, 始構造と終構造が予め既知である場合に, その構造間遷移を抽出するために用いられる. ターゲットを必要としないPaCS-MDの拡張版であるnon-targeted PaCS-MD (nt-PaCS-MD)に関しては, 本論文 [23] を参照されたい. また, PaCS-MD実行時のパラメータ (例えば, サイクル毎に選択する初期構造数, MDの時間スケール, 温度など) が探索効率にどのように影響してくるかに関しても議論した. 詳細は, 本論文 [24,25,26,27] を参照されたい.

(D) Taboo SeArch algorithm (TBSA)

TBSA [28] は, タンパク質状態分布の逆状態分布を利用して出現確率の小さいタンパク質状態を特徴付け, 構造リサンプリングの初期構造に選択することで, 構造遷移を促進する. 逆状態分布を定義するための反応座標は様々考えられるが, 一番単純な例として, エネルギー空間に射影された状態分布を考え, その逆状態分布を定義し, これを重みとして, 構造リサンプリングの初期構造を選択する. 計算手順として, 構造リサンプリングにより得られるスナップショットを用いて, 逆分布関数をアップデートしながら初期構造選択を繰り返し, 短時間MDにより構造リサンプリングを繰り返す. 収束条件として, タンパク質状態分布関数が変化しなくなるまで 初期構造選択と構造リサンプリングを繰り返す. エネルギーを反応座標としたTBSAの適用例として, Chignolinのフォールディング経路を抽出した [29]. また, よりシステムサイズの大きなタンパク質に関して, 水中のグルタミン結合タンパク質(223残基)にTBSAを適用し, Open-Closed構造間の大規模構造変化を抽出した [30].

TBSAと類似したアルゴリズムで, 探索構造の中から出来るだけ構造相関が小さくなるように初期構造を選択しながら構造リサンプリングを繰り返し, 効率的な構造サンプリングを実現するStructural Dissimilarity Sampling (SDS) も, ターゲットが未知である場合に非常に有効な手法である. 計算手順の詳細に関しては, 本論文 [31] を参照されたい.

2.3 自由エネルギー解析

反応経路の近似解であるreactive trajectoryが求まれば, U.S.法とWHAMを組み合わせることで, 自由エネルギーを計算することができる. 具体的には, 注目している系に対してある特定の参照構造 (距離や角度など) を安定化 (逆に言うと, 参照構造以外を不安定化) させるバイアスポテンシャル$w_i^{\text{biased}}\left(\xi \right)=k_{0,i}^{\text{ref}}{\left(\xi -{\xi }_{0,i}^{\text{ref}}\right)}^2/2$を加えて得られた (偽の) 分布関数${\rho }_i^{\text{biased}}\left(\xi \right)$の情報から, 元々のバイアスのない分布関数${\rho }^{\text{unbiased}}\left(\xi \right)$を再構築し, 自由エネルギーを計算する. ${\rho }^{\text{unbiased}}\left(\xi \right)$を${\rho }_i^{\text{biased}}\left(\xi \right)$の線形結合で表すと,   

${\rho }^{\text{unbiased}}\left(\xi \right)=\frac{\sum _{i=0}^L{n}_i{\rho }_i^{\text{biased}}\left(\xi \right)}{\sum _{i=0}^L{n}_{i}exp\left(-\beta \left(w_i^{\text{biased}}\left(\xi \right)-g_i\right)\right)}$

のようになる. ここで, n_iは${\rho }_i^{\text{biased}}\left(\xi \right)$を計算するために使用したMDステップ数であり, 分母にある$g_i$は   

$g_i=-k_{\text{B}}Tln\left(\int dR{e}^{-\beta w_i^{\text{biased}}\left(\xi \right)}{\rho }^{\text{unbiased}}\left(\xi \right)\right)$

から求められる. g_iも${\rho }^{\text{unbiased}}\left(\xi \right)$に依存していることから, これら2式は繰り返し計算によって収束するまで評価することになる. 実際の計算で注意しなければならないのが, U.S.法により得られたトラジェクトリの反応座標空間における分布の重なりである. WHAMでは隣のサンプル窓の分布関数との重なりが小さいと収束が遅い, あるいは収束しない恐れがあるため, サンプル窓を逐次増やしながら収束するまで計算を実行する必要がある.

3 結果と考察

3.1 FFMによるT4リゾチームのドメイン運動解析

FFMの適用例として, 加水分解酵素T4リゾチームのOpen-Closed構造変化を抽出した計算結果を紹介する. 出発構造として, T4リゾチーム (野生型) のOpen構造を選択した. 構造遷移を誘起するための振動モードとして, Open構造を出発した MD (10 ns) を実行しておき, 得られたトラジェクトリを用いて主成分分析を行い, 第1–2主成分座標を反応座標として選択した. 各主成分モードの寄与率を調べたところ, 上位2つの主成分モードにより全体のタンパク質構造揺らぎの80%以上を占めていたことから, Open-Closed構造変化を十分記述できると判断した.

ここで, FFMを用いて抽出したT4リゾチームのOpen構造からClosed構造へ至る構造遷移プロセスを示す [16]. Figure 4 は, サイクルごとに得られたトラジェクトリの第1–2主成分座標空間への射影である (青). 比較のため, Open構造 (Closed構造) を出発構造とした通常の長時間MDのトラジェクトリの射影も示す (赤) (緑). サイクルを繰り返すごとに, Open構造からその周辺へと射影点が拡大していき, 15サイクル (数十ns) 程でClosed構造 (マゼンタ) へ達している. 比較計算として, 通常のMDシミュレーションを, Open構造をスタートし1μs実行したが, Closed構造へ至る構造遷移を抽出できなかったことからも, FFMの構造探索効率の高さがうかがえる.

Figure 4

 : Conformational transition pathways from the open to the closed states of T4 lysozyme reproduced by FFM

3.2 OFLOOD法によるタンパク質Villinのフォールディング解析

OFLOOD法の適用例として, 小タンパク質Villin (35残基) のフォールディング経路抽出過程を示す. 反応座標として, 先行研究 [32] で定義されたRMSDを用いた. 具体的には, 2次元RMSDとして, Helix1-Helix2 (segment A)およびHelix2-Helix3 (segment B)からなる部分構造に対して天然構造から測定したRMSDを計算し, 状態を2次元空間に射影した. 初期構造として完全に伸びきったアミノ酸鎖をモデリングし, 溶媒をGB/SAでモデル化した後, OFLOOD法による構造サンプリングを20サイクル繰り返した.

Figure 5 に「はずれ値」と短時間MDによる構造リサンプリングから得られたトラジェクトリの射影を示す. サイクルを繰り返すごとに, 「はずれ値」 (黒) が分布の端を広げていき, 10サイクル程でVillinの天然構造 ( × ) をサンプルしている. 20サイクル終了後のX戦結晶構造から測定した最小C_α RMSDは0.60 Åであり, 精度よく天然構造をサンプルできている. 更に, 先行研究では抽出できていなかった副フォールディング経路 (Figure 5を参照) も抽出することができた. 計算効率としては, 天然構造をサンプルする (C_α RMSD < 1.0 Å) 累積計算時間が135.6 nsであり, 先行研究 [3] のレプリカ交換MD (8 μs) と比較して, 非常に効率的にフォールディング経路を抽出できた. 我々の計算結果は, nsオーダーの計算コストで, μsオーダーのタンパク質フォールディング過程を抽出可能であることを示しており, OFLOOD法を基軸としたアプリケーションが期待される.

Figure 5

 . Major and minor folding pathways of Villin reproduced by OFLOOD method

3.3 OFLOOD法によるその他のタンパク質のフォールディング解析

OFLOOD法のVillinへの適用例から示唆されるように, タンパク質フォールディングでは, 初期段階において2次構造の形成が重要になる. そこで, 2次構造の形成と全体のフォールディング機構の相関を調べるため, Villinとは別の小タンパク質1FMEと1E0LにOFLOOD法を適用した [33]. これら1 FMEと1E0Lは, αヘリックスとβシートを構成要素とする典型的なファーストフォールディングのタンパク質である. (Figure 6を参照)

Figure 6

 : The native structures of 1FME and 1E0L

Figure 6に示すように, 1FMEは, βシート(赤), ループ(緑), αヘリックス(青)で構成される. また1E0Lは, 3つのβシート(赤, 緑, 青)で構成される. ここで, 反応座標を定義するため, 部分的なRMSDを定義した. 具体的には, 1FMEに関して中央のループ領域, 1E0Lに関してβシート領域(1E0L)をオーバーラップする2つのセグメントAとBを定義し, 反応座標とした. 最終的に, 両セグメントに対する部分的なRMSD (RMSD_segmentAとRMSD_segmentB)が共に小さな値になった時, 天然構造へのフォールディンが完了したことを意味する. これらの反応座標を用いて, 完全に伸びきった変性構造をスタートし, Villinと同様の計算スキームでOFLOOD法による構造サンプリングを50サイクル繰り返した. Figure 7とFigure 8 にその計算結果を示す. 各サイクルにおいて100個のはずれ値 (緑) を初期構造として選択し, 100psの短時間MDをリスタートさせた. 両構造共に, 50サイクル終了後, 十分広い構造空間を探索できており, 天然構造近傍も高精度にサンプルできている. また, OFLOOD法により計算したポテンシャルエネルギー面から, 両構造に関して2つの支配的なフォールディング経路が存在することが示唆された. 興味深いことに, Villinと同様に2次構造の形成が全体のフォールディングの律速となっている可能性が高いことが明らかになった. 特に1FMEの場合には, フォールディング初期段階においてαヘリックスとβシートの形成がトリガーとなって, 全体のフォールディングを誘導していることが示唆された. (Figure 9を参照) これに対し1E0Lでは, N端とC端のβシートの形成に伴い, 全体のフォールディングが開始されることが示唆された. (Figure 10を参照)

Figure 7

 . (left) Potential energy surface of 1FME estimated by OFLOOD method. (right) Conformational subspace of 1FME sampled by OFLOOD method

Figure 8

 . (left) Potential energy surface of 1E0L estimated by OFLOOD method. (right) Conformational subspace of 1E0L sampled by OFLOOD method

Figure 9

 : Multiple folding pathways of 1FME predicted by OFLOOD method

Figure 10

 : Multiple folding pathways of 1E0L predicted by OFLOOD method

いずれのケースも, 2つの支配的なフォールディング経路が検出され, 2次構造の形成に伴いフォールディングが開始される普遍的なメカニズムが明らかになった. これらの計算結果は, 従来の2状態遷移のような単純なフォールディングメカニズムではなく, フォールディング途中に存在する中間体構造 (準安定状態) を経由して天然構造へ至る, 複雑なメカニズムが存在することを示唆している.

3.4 PaCS-MDによる基質-酵素複合体形成における誘導適合解析

PaCS-MD法の適用例 [34] として, ナイロンオリゴマー分解酵 (NylB) と基質であるALD (Ahx-liner dimer) の結合後の誘導適合経路抽出過程を示す. 野生株 (WT) では, 基質結合後にループ領域 (Gln166-Val177) にあるTyr170がALDと水素結合を形成し, 大規模な構造変化(誘導適合)を起こす (Figure 11 (a)を参照). このTyr170をPheに置換したY170F変異体では, 水素結合ができなくなったことで基質結合の安定性が低下し, 酵素活性が落ちるものと考えられていた.

Figure 11

 . Structure of active site

(a) Free energy analysis upon the induced fit (b) WT (c) Y170F mutant

Figure 11 (b)と(c)はWT, およびY170F変異体に対してPaCS-MDを20サイクル実行し, その後U.S.法-WHAMによって求められた自由エネルギー面を示している. 解析の結果, WTではループ移動前 (Open構造) からループ移動後 (Closed構造) への構造変化過程において, 4カ所の準安定点が存在することが明らかとなった. また, ループ部位の詳細な移動経路としては, 今回判明した準安定点をうまく段階的に通過することによって構造変化していることが判明した. さらに, Open構造からClosed構造への構造変化に伴い, 約1.4 kcal/molの安定化が得られることも明らかになった. 一方で, Y170F変異体では基質のアミド結合と水素結合できない代わりに, フェニル基がアルキル鎖と相互作用することで安定化していた. 準安定状態は2つに減少し, 最安定構造のベイスンの形状からわかるようにWTと比較して揺らぎ幅が大きい. 一方で, 全体的に自由エネルギー地形は大きく変化するものの, Open構造からClosed構造に至る安定化エネルギーはほぼ同程度であることが明らかになった. 以上のようにレアイベントである誘導適合過程に関しても, PaCS-MDが簡便かつ高速に反応経路解析ができることが分かる.

3.5 PaCS-MDによる細胞分裂ダイナミクスに関する動的秩序解析

細菌の細胞分裂に関係しているタンパク質FtsZは, Figure 12に示すように細胞膜の内側にリング状のフィラメント (Zリング) を形成し, ダイナミックに離合集散を繰り返すことで細胞膜に陥入を生じさせる. この細胞膜の陥入はZリングの収縮により起きると考えられるが, その分子メカニズムには未解明な部分が多い. 近年, 立命館大学•松村教授の研究グループにより黄色ブドウ球菌由来のGDP結合型FtsZのX線結晶構造が決定された. この結晶構造には, 同一結晶中に立体構造が大きく異なるFtsZが2状態(T-状態とR-状態)が同定されている. 同一種で状態の異なる2構造が得られた例はこれまでに報告されておらず, 大変貴重な実験結果である. そこで, T-R状態間の構造遷移経路の探索を行い細胞分裂タンパク質の動的秩序過程の解析を行った. 具体的にはPaCS-MDを適用し, 2状態構造間の構造遷移経路を探索した [35].

Figure 12

 . T-R state transition mechanism of FtsZ monomer from Staphylococcus aureus

PaCS-MDにより得られた遷移経路を解析したところ, 状態遷移における重要なアミノ酸残基のメカニズムを突き止めることができた. 具体的には, 29番目のアルギニン残基(Arg29) の側鎖のフリップがスイッチとなり, 状態遷移を制御していることが明らかになった. Figure 12 に示すように, Arg29の側鎖がフリップすることで188番目のアルギニン(Asn188)の側鎖と水素結合を形成し, FtsZ中央に存在しているヘリックスが捩れた構造から, 捩れが解消された直線的な形状に構造遷移するメカニズム解明した. また, PaCS-MDにより抽出したArg29-Asn188の水素結合距離の時系列データから, T-R状態遷移に伴いArg29-Asn188間に水素結合が形成されていることを明らかにした (Figure 12 を参照). 更に, FtsZモノマーは2構造間の構造遷移において基質であるGDPを段階的に認識•解除していることを明らかにし, 中間体構造を経て多段階に状態遷移していることを突き止めた. 本研究において, 実験と計算化学が密に連携することにより, FtsZモノマーの構造揺らぎとT-R状態間構造遷移を解明し, FtsZポリマーの離合集散の関係を解明するための足がかりを築くことができた. 今後はFtsZポリマーのシミュレーションも検討し, より現実の環境に近いモデルを構築して細胞分裂過程の動的秩序解明を進めて行く予定である.

4 まとめ

タンパク質などの構造が柔らかで複雑な系では実効的な自由度の数が多く, エネルギー最小経路を結ぶ唯一の反応経路ばかりではなく, 反応経路の周りの様々な配置の影響も無視できなくなってくる. 特に, 周りの水分子の配置まで考えると, 莫大な計算量が必要になることが考えられる. このような系を理論的に取り扱うためには, さらに効率的な経路探索手法, およびその周りの揺らぎの影響を含めた自由エネルギー計算手法が必須であることから, 本解説で提案した手法などは, 今後益々その重要性を増してくるものと期待される. 一方で本手法においては, 望みの構造変化を誘発する選択則をどのように選ぶのかという問題が残されている. 言い換えると, 反応座標をどのように設定するのかで, 得られる反応経路が変わる可能性もある. 本手法と従来のサンプリング手法やダイナミクス解析手法と組み合わせることで, さらなる効率的な方法が開発されるものと思われる.

Acknowledgment

本総説は, 文部科学省新学術領域研究「動的秩序と機能」の京都大学佐藤教授より依頼があり, 我々のグループの最新の結果に関する記事を執筆させて頂きました. 本研究を遂行するにあたり, 研究費の支援(課題No. JP26102525, No. JP 26105012, No. JP 26107004, JP15J03797, JP16H06164)を受けました.

参考文献

• [1]   D. E. Shaw, P. Maragakis, K. Lindorff-Larsen, S. Piana, R. O. Dror, M. P. Eastwood, J. A. Bank, J. M. Jumper, J. K. Salmon, Y. Shan, W. Wriggers, Science, 330, 341 (2010). doi:10.1126/science.118740920947758
• [2]   K. Lindorff-Larsen, S. Piana, R. O. Dror, D. E. Shaw, Science, 334, 517 (2011). doi:10.1126/science.120835122034434
• [3]   N. Nakajima, H. Nakamura, A. Kidera, J. Phys. Chem. B, 101, 817 (1997). doi:10.1021/jp962142e
• [4]   U.H.E. Hansmann, Y. Okamoto, F. Eisenmenger, Chem. Phys. Lett. 259, 321, (1996). DOE: doi:10.1016/0009-2614(96)00761-0
• [5]   Y. Sugita, Y. Okamoto, Chem. Phys. Lett., 314, 141 (1999). doi:10.1016/S0009-2614(99)01123-9
• [6]   D. E. Shaw, J. C. Chao, M. P. Eastwood, J. Gagliardo, J. P. Grossman, C. R. Ho, D. J. Lerardi, I. Kolossváry, J. L. Klepeis, T. Layman, C. McLeavey, M. M. Deneroff, M. A. Moraes, R. Mueller, E. C. Priest, Y. Shan, J. Spengler, M. Theobald, B. Towles, S. C. Wang, R. O. Dror, J. S. Kuskin, R. H. Larson, J. K. Salmon, C. Young, B. Batson, K. J. Bowers, Commun. ACM, 51, 91 (2008). doi:10.1145/1364782.1364802
• [7]   R. Harada, Y. Takano, T. Baba, Y. Shigeta, Phys. Chem. Chem. Phys., 17, 6155 (2015). doi:10.1039/C4CP05262E25659594
• [8]   H. J. C. Berendsen, J. P. M. Postma, W. F. van Gunsteren, A. DiNola, J. R. Haak, J. Chem. Phys., 81, 3684 (1984). doi:10.1063/1.448118
• [9]   M. P. Allen, G. P. Tildesley, Computer Simulation of Liquid, Oxford University Press, New York, (1987).
• [10]   S. Nosé, J. Chem. Phys., 81, 511 (1984). doi:10.1063/1.447334
• [11]   W. G. Hoover, Phys. Rev. A, 31, 1695 (1985). doi:10.1103/PhysRevA.31.1695
• [12]   G. M. Torrie, J. P. Valleau, J. Comput. Phys., 23, 187 (1977). doi:10.1016/0021-9991(77)90121-8
• [13]   G. M. Torrie, J. P. Valleau, J. Chem. Phys., 66, 1402 (1977). doi:10.1063/1.434125
• [14]   A. M. Ferrenberg, R. H. Swendsen, Phys. Rev. Lett., 63, 1195 (1989). doi:10.1103/PhysRevLett.63.119510040500
• [15]   M. Souaille, B. Roux, Comput. Phys. Commun., 135, 40 (2001). doi:10.1016/S0010-4655(00)00215-0
• [16]   R. Harada, Y. Takano, Y. Shigeta, J. Chem. Phys., 140, 125103 (2014). doi:10.1063/1.4869594
• [17]   T. Nkamura, Y. Kamidoi, N. Wakabayashi, N. Yoshida, IPSJ J. Database, 46, 40 (2005), DOI.
• [18]   R. Harada, T. Nakamura, Y. Takano, Y. Shigeta, J. Comput. Chem., 36, 97 (2015). doi:10.1002/jcc.2377325363340
• [19]   R. Harada, T. Nakamura, Y. Shigeta, Chem. Phys. Lett., 639, 269 (2015). doi:10.1016/j.cplett.2015.09.031
• [20]   R. Harada, T. Nakamura, Y. Shigeta, J. Comput. Chem., 37, 724 (2016). doi:10.1002/jcc.2425526611770
• [21]   R. Harada, T. Nakamura, Y. Shigeta, Bull. Chem. Soc. Jpn., 89, 1361 (2016). doi:10.1246/bcsj.20160244
• [22]   R. Harada, A. Kitao, J. Chem. Phys., 139, 035103 (2013). doi:10.1063/1.481302323883057
• [23]   R. Harada, A. Kitao, J. Chem. Theory Comput., 11, 5493 (2015). doi:10.1021/acs.jctc.5b0072326574337
• [24]   R. Harada, Y. Shigeta, J. Comput. Chem., 38, 2671 (2017). doi:10.1002/jcc.2506028861895
• [25]   R. Harada, Y. Shigeta, Bull. Chem. Soc. Jpn., 90, 1236 (2017). doi:10.1246/bcsj.20170177
• [26]   R. Harada, Y. Shigeta, Chem. Lett., 46, 862 (2017). doi:10.1246/cl.170207
• [27]   R. Harada, Y. Shigeta, Mol. Simul., 44, 206 (2018). DOI.org/10.1080/08927022.2017.1362696 doi:10.1080/08927022.2017.1362696
• [28]   R. Harada, Y. Takano, Y. Shigeta, J. Comput. Chem., 36, 763 (2015). doi:10.1002/jcc.2385425691321
• [29]   R. Harada, Y. Takano, Y. Shigeta, Chem. Phys. Lett., 630, 68 (2015). doi:10.1016/j.cplett.2015.04.039
• [30]   R. Harada, Y. Takano, Y. Shigeta, J. Chem. Theory Comput., 12, 2436 (2016). doi:10.1021/acs.jctc.6b0008227070761
• [31]   R. Harada, Y. Shigeta, J. Comput. Chem., 38, 1921 (2017). doi:10.1002/jcc.2483728558119
• [32]   H. Lei, C. Wu, H. Liu, Y. Duan, Proc. Natl. Acad. Sci. USA, 104, 4925 (2007). doi:10.1073/pnas.060843210417360390
• [33]   R. Harada, Y. Takano, Y. Shigeta, J. Comput. Chem., 38, 790 (2017). doi:10.1002/jcc.2474828192630
• [34]   T. Baba, R. Harada, M. Nakano, Y. Shigeta, J. Comput. Chem., 35, 1240 (2014). doi:10.1002/jcc.2361424752461
• [35]   J. Fujita, R. Harada, Y. Maeda, Y. Saito, E. Mizohata, T. Inoue, Y. Shigeta, H. Matsumura, J. Struct. Biol., 198, 65 (2017). doi:10.1016/j.jsb.2017.04.00828456664