Free energy landscape (FEL), is an important outcomes of molecular dynamics (MD) simulations, providing essential information for understanding the relations between conformational changes and functions. However, the FEL of biomolecular systems has multiple local energy minimum states, which makes it difficult to explore broad conformational space using conventional MD (CMD) simulations. Biomolecular systems cannot escape from one of multiple local energy minimum states separated by high energy barriers in a limited computational time. This is called multiple energy minima problem. Trapping into local energy minimum states typically leads to an inaccurate picture of the FEL. Hence, the sampling efficiency of the broader conformational space is essential for calculating the FEL.

In this thesis, I propose a novel multi-scale molecular dynamics simulation method for calculating free-energy landscapes of biomolecular systems, namely Multi-Scale Free Energy Landscape calculation method (MSFEL) to solve the multiple energy minima problem. This method is designed to simultaneously achieve high accuracy and sampling efficiency using the combination of coarse-grained (CG-) and all-atom (AA-) models. Firstly, the CG-model is used for calculating a rough energy landscape over broad conformational space. Subsequently, the AA-model is used to increase the sampling accuracy of local energy landscapes. This multi-scale approach attempts to achieve high sampling efficiency in the CG-model and high-resolution sampling in the AA-model. MSFEL has also an advantage of compatibility with massively parallel computing as time consuming AA-MD simulations are performed completely in parallel. Therefore, MSFEL will be one of efficient free energy calculation methods using huge computer resources in the future.

MSFEL consists of the following four steps. STEP1: Wide conformational space is sampled by the CG-MD simulation. STEP2: Multiple CG-structures are selected and AA-model structures are generated by using databases. STEP3: Non-interactive massively parallel AA-MD simulations are performed to investigate the local energy landscapes. STEP4: The multiple trajectories obtained from the independent AA-MD simulations are combined with the weighted histogram analysis method (WHAM) and the FEL is calculated. Figure 1 shows the flowchart of MSFEL.

For a benchmark test, we firstly applied MSFEL to short peptide, Met-enkephalin (5-residue) in vacuum. Figure 2(a) shows the 100-ns trajectories of Met-enkephalin in vacuum simulated with all-atom CMD simulation started from two local energy minimum states, LM1 and LM2 projected on to the first principal component axis. As expected from the preceding works, Met-enkephalin is shown to be trapped in the initial states at least for 100-ns. Figure 2(b) shows the CG-MD trajectory projected on to the same first principal axis as in Figure 2(a), indicating frequent conformational transitions between LM1 and LM2 whereas conventional AA-MD failed to attain the transition within 100-ns. Thus more efficient conformational sampling is achieved in the CG-MD simulation.

After CG-MD simulation, we modeled 100 AA-structures from the CG-structures and then performed 1-ns×100 independent AA-MD simulations in parallel. In the free energy calculation, the total probability distribution is calculated from 100 trajectories obtained in the independent AA-MD simulations with the WHAM. The FELs obtained from MSFEL and the CMD starting from LM1 projected onto the first and second principal axes are shown in Figure 3(a) and (b), respectively.

Due to trapping in LM1, the CMD failed to calculate accurate FEL. In contrast, MSFEL is successful in identifying LM1 and LM2 as well as the transition state TS between LM1 and LM2. In addition to these known states, a novel local energy minimum state LM3 that has not been reported was found by MSFEL.

For other applications, we applied MSFEL to the following typical mini-proteins, tryptophan zipper (12-residue), chignolin (10-residue) and villin headpiece subdomain (35-residue) to calculate folding FELs. The results of each application also showed high computational efficiencies than the CMD simulations. Furthermore, MSFEL could also detect folding pathways and intermediate states during folding processes on each application. Figure4 (a) and (b) show one of the applications, the folding FEL of chignolin obtained from MSFEL and the CMD (100-ns) simulation, respectively. Compared to the CMD simulation, the misfolded and the intermediate conformations were also sampled by applying MSFEL.

In conclusion, the benchmark test and the applications of MSFEL demonstrated that this method is an efficient and useful method in the free energy landscape calculations of biomolecules the FEL obtained from MSFEL is useful information for understanding the relevant biological functions.

Figure 1. The flowchart of MSFEL

Figure 2. The time series of Met-enkephalin dynamics in vacuum projected on to the first principal component axis

Figure 3. The FEL of Met-enkephalin in vacuum projected onto the first and the second principal axes.

Figure 4. The folding FEL of Chignolin projected onto the donor and the acceptor distance of main-chain atoms.

本論文は4章からなる。第1章は序論であり,研究の背景と目的が述べられている。第2章では本研究で新しく開発したシミュレーション手法Multi-Scale Free Energy Landscape calculation method (MSFEL)について説明している。第3章ではMSFELをペプチドやタンパク質の自由エネルギー地形の計算に適用した結果について、記述している。先行研究や通常の分子動力学計算と比較し、その有効性を明らかにしている。第4章では、本研究の結論が述べられている。

タンパク質の機能を解明する上で、その天然構造と機能発現における構造変化を明らかにすることは非常に重要である。構造を解明する上で数値計算は重要な研究手法である。特に、計算で求められる自由エネルギー地形はタンパク質の天然構造への折り畳みや構造変化を理解する上で重要である。最近、コンピューターの性能向上において、並列計算可能なCPUコア数が著しく、増加しており、建設中の次世代スーパコンピューター「京」などでは数十万コアが利用可能になる。これらの計算資源を有効活用するには並列計算に適した計算手法の開発が重要となる。本研究では従来の方法に比べて並列計算により適した新しい自由エネルギー地形の計算法MSFELを提案している。

本論文では開発した計算法MSFELは4ステップから構成される。1)まず、粗視化タンパク質模型の分子シミュレーションを行い、自由エネルギー地形に広く分布した、複数(例では100個)の立体構造をサンプリングする。2)次に、その粗視化模型の立体構造から、全原子模型の立体構造を作成する。3)各々の立体構造から、初期構造に調和振動子型の拘束をかけた全原子シミュレーションを行う。このステップがもっとも計算量を必要とするが、初期構造ごとに完全に独立に計算できるので、並列効率100%の並列計算が可能である。4)重みをかけ直してカノニカル集団での確率分布を求める計算手法であるweighted histogram analysis method (WHAM)実際の使用で有用な結果をえるためには、粗視化模型や全原子模型の立体構造の構築手法の選択が重要である。ここでは、粗視化模型として溶媒、タンパク質の側鎖を陽に含まず、主鎖の炭素原子を残した模型(Cα一based CG model)を用いている。全原子模型の立体構造の構築に、既成のプログラムMaxSproutとBBQ、SCWR上を用いている。

計算法MSFELを2つのペプチドMet-enkephalin、Tryptophan zipper(Trpzip2)、タンパク質Chignolin、タンパク質VilinのサブドメインHP35に適用して、自由エネルギー地形を求めている。Met-enkephalinとTrpzip2では粗視化模型のシミュレーションで、カノニカル集団の分子動力学計算を用いているが、ChignolinとHP35ではレプリカ交換分子動力学計算を用いている。大きな分子の計算では後者がより有効である。どの系においても、通常のカノニカル集団の分子動力学計算では天然構造の周辺しか探索できていないが、MSFELでは広い範囲の自由エネルギー地形が得られている。また、Met-enkephalinにおいては、先行研究の結果を再現するだけでなく、新しい準安定構造を発見している。今回用いた粗視化模型には全原子模型の立体構造構築の際、側鎖の作成配置に任意性ができるが、側鎖の回転緩和が十分速く行われており、結果に影響を及ぼさないことを確認している。このように結果の妥当性と手法の有効性が4つの適用例全てにおいて、確認されている。

ChignolinとHP35において先行研究では溶媒を連続体近似した模型でしか、自由エネルギー地形を求めていなかった。本研究では、溶媒を含んだ全原子計算でより詳細なエネルギー地形を明らかにしている。Chignolinでは知られていた準安定な折り畳み構造と天然構造の中間に新しい準安定構造を見つけている。また、HP35では従来言われていなかった折り畳み経路の存在を示唆する結果が得られている。このように本論文では、新しい計算法を提案するだけでなく、実際にタンパク質計算に使用することによって、新しい構造を解明している。

以上のように本論文では,タンパク質の自由エネルギー地形を並列計算機を用いて効率的に計算できる新計算手法を開発している。4つのペプチドもしくはタンパク質に適用し、その有効性を示すとともに、新しい中間体、折り畳み経路の発見という新しい知見も得ている。提案された新手法はタンパク質だけではなく、他の高分子などにも、今後拡張可能な計算手法であり、生物物理学、計算科学の発展が期待できる。また、本論文の内容の一部はChemical Physics Lettersで出版予定である。なお本論文は指導教員である北尾彰朗氏との共同研究であるが,論文提出者が主体となって研究を行ったもので,論文提出者の寄与が十分であると判断する。

したがって,審査員全員の一致により,博士(理学)の学位を授与できると認める。