用語

圧力

分子動力学における圧力は、運動エネルギーとバイラルから計算できます。

変動

シミュレーション内で圧力を利用する場合でも、シミュレーションボックスの圧力値は大幅に変動します。瞬間的な圧力は意味がなく、定義されていません。数ピコ秒の時間スケールでは、通常、真の圧力を正確に示すものではありません。この変動は、圧力はマクロな性質であり、微視的なスケールで測定および調整される際に、時間平均としてのみ正確に測定できるため、完全に正常なものです。変動の大きさや速度は、システムの原子数、使用される圧力カップリングの種類、およびカップリング定数の値によって異なります。数百バールのオーダーの変動は一般的です。216個の水分子の箱の場合、500〜600バールの変動が標準です。変動は粒子の数の平方根に比例するため、21600個の水分子(100倍大きい)のシステムでも、50〜60バールの圧力変動が残ります。

周期境界条件

周期境界条件 (PBC) は、分子動力学シミュレーションにおいて、有限サイズによる境界効果の問題を回避し、システムを無限に近い状態にするために使用されます。ただし、周期的な効果が生じる可能性があります。

初心者が軌跡を視覚化する際、しばしば問題を観察しているかのように誤解することがあります。

  • 分子(分子)は箱の中心に留まりません。

  • 分子(の一部)が箱から拡散しているように見えます。

  • 穴が作成される、または

  • 分子が破壊されたり、

  • その単位格子は菱面体十二面体または立方体八面体でしたが、シミュレーション後には傾いた立方体のように見えます。

  • シミュレーションセル全体に、異常な結合が多数発生します。

これは問題やエラーではなく、期待されることです。

PBCが存在すると、ある面(例えば右側)からシミュレーションボックスを通過し、別の面(例えば左側)から再びボックスに入ってくる原子が存在することを示します。例えば、大きなタンパク質の場合、タンパク質が突出している面とは反対側のシミュレーションボックスの面を見ると、溶媒の中には穴が開いているように見えることがあります。これは、分子がシミュレーション中に、初期位置から自由に拡散できるためです。分子は、特定の場所に固定されているわけではありません。シミュレーション中は、ボックスが特定の場所に固定されているわけではありません。また、分子は、シミュレーション中に自然に全体として存在するというわけではありません。さらに、任意の形状の周期性セルは、直方体(別名トリキイン体)として表現でき、|Gromacs|は、ボックスの初期形状に関係なく、内部でこれを実現します。

これらの視覚的な問題は、シミュレーションの終了後に、gmx trjconv を使用してオプションの入力を使用することで修正できます。これにより、軌跡ファイルが処理されます。同様に、原子の位置の RMSD などの分析は、周期性効果を考慮して参照構造と比較する場合に誤る可能性があります。この問題を解決するには、gmx trjconv を使用します。いくつかの複雑なケースでは、複数の操作が必要となるため、gmx trjconv を複数回実行する必要があります。

詳細については、リファレンス マニュアル の対応するセクションを参照してください。

推奨されるワークフロー

gmx trjconv を使用して周期的な効果を修正し、可視化や分析に適した結果を得ることは、複雑な場合があります。 複数の呼び出しが必要になることがあります。 例えば、タンパク質とリガンドを一緒に保持するために、カスタムのインデックスグループを作成する必要があるかもしれません。 以下の手順を順番に実行(必要な手順のみ)することで、望ましい結果を得られるはずです。 各手順の詳細については、gmx trjconv -h を参照する必要があります。 これは意図的なものであり、「私が望むように」という魔法のようなレシピはありません。 まず、自分が何をしたいのかを決定する必要があります。 :-)

  1. まず、分子を完全に組み立てることで、完全にすることができます。

  2. 分子/粒子をグループ化したい場合は、それらをグループ化してください。

  3. もし、ジャンプを削除したい場合は、軌跡から最初のフレームを基準として抽出し、その後、その最初のフレームを基準として -pbc nojump を使用してください。

  4. システムを基準に基づいて中央に配置します。これにより、システムが移動するため、この手順の後に -pbc nojump オプションを使用しないでください。

  5. もしかしたら、他の -pbc または -ur オプションと一緒に、すべてを同じ箱に入れるのが良いかもしれません。

  6. 必要に応じて、結果の軌跡を別の参照構造に適合させ、その後、PBCに関連するオプションは使用しないでください。

ポイント3に関連して、問題は、gmx trjconv`が-sオプションで提供された参照構造を使用して、最初のフレームからのジャンプを削除することです。参照構造(入力ファイル)がクラスタリングされていない場合、または全体として使用されない場合、`-pbc nojump``を使用すると、ステップ1と2を元に戻します。

サーモスタット

サーモスタットは、シミュレーションサンプルを正しいエンブレム(つまり、NVTまたはNPT)に配置するのに役立つように設計されています。具体的には、システム内の温度をある程度調整することでです。まず、私たちが「温度」という言葉をどのように理解しているのかを明確にする必要があります。シミュレーションでは、「瞬間(運動)温度」は、通常、システム全体の運動エネルギーを使用して、均等分割の定理から計算されます。つまり、温度はシステム全体の運動エネルギーから計算されます。

したがって、熱源の目的とは何か? 実際には、目的は一定の温度を維持することではなく、これは総運動エネルギーを固定することであり、これはナンセンスであり、NVTまたはNPTの目的にも合わないものです。 むしろ、システムの平均温度が正しいことを保証することです。

この理由を理解するには、グラスに置かれた水のイメージを考えてみてください。もし、グラス内の小さな領域にあるいくつかの分子を非常に詳細に観察し、それらの運動エネルギーを測定できるとします。この少数の粒子における運動エネルギーが常に正確に一定であるとは予想しません。むしろ、粒子数が少ないため、運動エネルギーに変動が生じることを予想します。粒子数をどんどん大きくしていくと、平均値の変動は小さくなり、最終的にはグラス全体を観察したときに、「一定の温度」であると判断します。

分子動力学シミュレーションは、水の一杯と比較して一般的に非常に小さいため、より大きな変動が生じます。したがって、熱源の役割を、温度を一定に保ち、

  1. 適切な平均気温、および

  2. 適切なサイズの変動。

詳細については、リファレンス マニュアル の関連セクションをご覧ください。このセクションには、温度結合の適用方法と、現在利用可能な種類に関する情報が記載されています。

どうすればいいか

一般的な推奨事項に関するヒント:

  • 理想的には、適切な温度分布を測定するサーモスタットを使用し、同時に正しい平均温度も表示するようにしてください(詳細については、対応するマニュアルの該当箇所を参照)。

  • 少なくとも: 正しい平均温度を示すサーモスタットを使用し、それが妥当な場合に、システム内のコンポーネントに適用してください(「やってはいけないこと」の最初の項目を参照)。 一部の場合は、``tc-grps = System``を使用すると、「高温の溶媒/低温の溶質」の問題が発生する可能性があります(「参考文献」の3番目の項目を参照)。

やってはいけないこと

一般的に推奨されない、いくつかのプラクティスのヒント:

  • システム内の各コンポーネントに対して、個別のサーモスタットを使用しないでください。一部の分子動力学サーモスタットは、熱力学の限界内でのみ効果的です。サーモスタットを使用するグループは、十分な大きさである必要があります。たとえば、小さな分子、タンパク質、水に対してそれぞれ別のサーモスタットを使用する場合、予測が困難なエラーやアーティファクトを導入する可能性があります。特に、水溶液中のイオンを、別のグループとは別に制御しないでください。タンパク質シミュレーションの場合、「tc-grps = Protein Non-Protein」を使用するのが一般的です。

  • 少数の自由度を持つシステムでは、多数の自由度でしか機能しないサーモスタットの使用は避けてください。たとえば、エンド状態において、システム内に非常に少数の自由度を持つコンポーネント(つまり、相互作用のない小さな分子)が存在する場合の、Nosé-HooverまたはBerendsenのサーモスタットの使用は避けてください。

Further reading

  1. Cheng, A. & Merz, K. M. 蛋白質のダイナミクス研究におけるNosé-Hoover 鎖アルゴリズムの応用。 J. Phys. Chem. 100 (5), 1927–1937 (1996).

  2. Mor, A., Ziv, G. & Levy, Y. 不均一な自由度のタンパク質のシミュレーション:熱板の効果。 J. Comput. Chem. 29 (12), 1992–1998 (2008).

  3. Lingenheil, M., Denschlag, R., Reichold, R. & Tavan, P. 「熱溶媒/低温溶媒」問題の再検討。 J. Chem. Theory Comput. 4 (8), 1293–1306 (2008).

エネルギーの節約

原則として、分子動力学シミュレーションは、全エネルギー、全運動量、および(非周期系の場合)全角動量を保存する必要があります。しかし、いくつかのアルゴリズムおよび数値上の問題により、必ずしもそうとは限りません。

  • カットオフ処理と/または長距離静電相互作用処理(詳細は、Van Der Spoel, D. & van Maaren, P. J.「シミュレーションにおける液体の層構造アーチファクトの起源」*J. Chem. Theor. Comp.* 2, 1–11 (2006)。)

  • ペアリストの扱い、

  • 制約アルゴリズム(例:Hess, B. P-LINCS: 分子シミュレーションのための並列線形制約ソルバー。 J. Chem. Theor. Comp. 4, 116–122 (2008)。)。

  • 統合ステップサイズ。

  • 温度連動圧力連動.

  • 丸め誤差(特に単精度の場合)、たとえば大きな数を引く場合(Lippert, R. A. et al. J. Chem. Phys. 126, 046101 (2007)。)。

  • 統合アルゴリズムの選択(|Gromacs|では通常、leap-frog)

  • 重心の動きの除去: 複数のグループでこれを実行すると、エネルギー保存則が破られます。

平均的な構造

さまざまな GROMACS ユーティリティは、平均構造を計算できます。おそらく、これは、例えば、エンザイム平均NMR構造から得られたアイデアに基づいています。場合によっては、平均構造を計算することが有益です(例えば、平均二乗変動(RMSF)を計算するためのステップとして、すべての原子の位置の平均が必要になります)。

ただし、平均的な構造が必ずしも意味のあるものではないことに注意することが重要です。例を挙げると、もし私がボールを左手で、次に右手で交互に保持するとしたら、ボールの平均的な位置はどこでしょうか?それは、常にボールが左手か右手かという状態の中間地点です。同様に、構造の場合、別々の安定な構造状態が存在する場合、平均は意味をなさなくなる傾向があります。これは、側鎖レベル、または一部のバックボーン領域、またはさらには、完全ならせんや二次構造の構成要素で発生する可能性があります。

したがって、分子動力学シミュレーションから平均構造を導き出し、物理的に不自然な結合長や奇妙な構造などのアーティファクトが見つかったとしても、必ずしも問題があるとは限りません。これは、シミュレーションから得られた平均構造が、必ずしも物理的に意味のある構造とは限らないことを示しています。

爆発

シミュレーションの失敗 は、非常に専門的な用語で、一般的なシミュレーションの失敗を説明するために使用されます。具体的には、通常、受け入れられないほど大きな力が作用し、その結果、インテグレーターが故障するという状況を指します。

もう少し背景を説明するために、分子動力学はニュートンの運動方程式を数値的に積分するために、小さな離散的な時間ステップを使用し、これらの時間ステップを使用して、前の時間ステップでの速度、位置、および力の値から、新しい速度と位置を決定します。もしある時間ステップで力が大きくなりすぎると、次の時間ステップに進む際に、速度/位置に極めて大きな変化が生じる可能性があります。通常、これはエラーの連鎖を引き起こします。ある原子が1つの時間ステップで非常に大きな力を受けるため、次の時間ステップで制御不能にシステム全体を移動し、希望する位置を超えて移動したり、別の原子やそれに類似した物体に衝突したりします。これにより、さらに大きな力が作用し、制御不能な動きが繰り返され、最終的にはシミュレーションパッケージがクラッシュします。制約付きのシミュレーションの場合、この現象の最初の兆候は、通常、LINCSまたはSHAKEに関する警告またはエラーであることになります。これは、制約が問題の原因ではないだけで、単に最初にクラッシュするものであるためです。また、表形式または1-4の相互作用が、表でサポートされている距離を超えていることに関する警告も同様です。これは、あるコンピューターシステムでは問題が発生する一方で、別のシステムでは、異なるコンピューターシステムでのこれらの計算の数値的な再現が不可能であるため、安定したシミュレーションが実現する可能性があります。

考えられる原因としては:

  • 十分に小さく表示しなかった。

  • あなたの初期構造は不適切であり、おそらく立体的な衝突が起きています。

  • タイムステップが大きすぎ(特に制約条件を選択した場合)、

  • あなたは、ソフトコアを使用せずに、自由エネルギー計算において粒子挿入を行っています。

  • 現在、不適切な圧力を使用しています(たとえば、平衡状態でない場合にBerendsenを使用することは有効ですが、その場合は後でより正確な圧力を制御するアルゴリズムに切り替える必要があります)。

  • 不適切な温度結合を使用している可能性があります。たとえば、不適切なグループを使用している場合や、

  • あなたの位置制限が、システムに存在する座標と比べて異なっており、

  • システム内に、他の水分子とは隔離された単一の水分子が存在します。

  • あなたは gmx mdrun に存在するバグを経験しています。

「爆発」が発生する主な原因は、特定のタイムステップサイズに対して大きすぎる力が作用することであるため、いくつかの基本的な解決策があります。

  • 「力を過度に大きくしないようにしてください。」

  • より小さいタイムステップを使用してください。

より良いシステム準備は、シミュレーションの初期段階で問題が発生した場合に、力が大きくなるのを防ぐための方法です。

システムが不安定な場合の診断

システムが不安定になっている問題を解決することは、特に初心者にとっては困難な場合があります。以下は、そのような状況に対処する際に役立つ一般的なヒントです。

  1. もしクラッシュが比較的早い段階(数ステップ程度)で発生している場合は、「nstxout」(または「nstxout-compressed」)を1に設定し、すべての可能なフレームをキャプチャします。これにより得られる軌跡を観察し、最初に不安定になる原子/残基/分子を特定できます。

  2. 問題を簡略化して、原因を特定することを試みる:

    • もし新しい溶媒の容器をお持ちの場合、不安定さが分子の構造的な問題によるものか、それとも初期設定における衝突によるものかを調べるために、分子を最小化し、シミュレーションを実行してみてください。

    • もし、タンパク質とリガンドのシステムをお持ちの場合は、目的の溶媒中でタンパク質だけをシミュレーションしてみてください。もし安定であれば、真空中でリガンドをシミュレーションし、その構造が安定な配置、エネルギーなどを与えるかどうかを確認してください。

    • 複雑なアルゴリズムの使用を避け、特に、事前に十分に調整を行っていない場合は、その使用を控えてください。

  3. システム内のさまざまなコンポーネントのエネルギーを監視するには、gmx エネルギー を使用します。分子内的な項が急激に変動する場合は、これが不適切な結合パラメータを示している可能性があります。

  4. 必ず、エラーメッセージ(例:`gmx pdb2gmx`を実行する際にアトムが不足している、`gmx grompp`を実行する際に名前が一致しないなど)を無視したり、本来避けるべき回避策(例:`gmx grompp -maxwarn`を使用するなど)を使用したりして、あなたのトポロジーが正しく構成され、正しく解釈されていることを確認してください。

  5. 必ず、選択した力場とシステムの種類に適した設定が mdp`ファイルに設定されていることを確認してください。特に重要な設定は、カットオフの処理、適切な近傍探索間隔(``nstlist`)、および温度カップリングです。不適切な設定は、システムの初期設定が妥当であっても、モデルの物理的な計算に問題を引き起こす可能性があります。

明示的な溶媒を使用しない場合、生成ランよりも短い時間ステップでイキベーションを開始することで、エネルギーの均一化をより安定的に行うことができます。

以下のような一般的な状況で、不安定性が頻繁に発生することがあります。通常、これはシステムに新しい種(リガンドやその他の分子)を導入した場合に起こります。問題の原因を特定するために、以下の方法でシステムを簡略化してください(例:タンパク質とリガンドの複合体の場合)。

  1. タンパク質(水溶液中)は、それ自体で十分に折りたたまれていますか? これは、座標とシステムの整合性、および準備のテストです。 これが失敗した場合、おそらく gmx pdb2gmx を実行する際に問題が発生したか、または gmx genion がタンパク質に非常に近いイオンを配置した可能性があります(これはランダムなことです)。

  2. 以下のテキストを翻訳してください

  3. (もし前の項目が成功した場合) 結合子は水中で最小化されますか、または、水中で結合子の短いシミュレーションは成功しますか?

その他の問題の原因となりうるものは、生体分子の構造自体にあります。

  1. 「-missing」オプションを使用して`gmx pdb2gmx`を実行した場合は、実行しないでください。代わりに、欠落した座標を再構築してください。

  2. 長/短結合警告を結合長の変更によって上書きしましたか? もしそうなら、それは避けてください。おそらく、原子が欠けているか、非常に不適切な入力ジオメトリを使用している可能性があります。

分子動力学

分子動力学 (MD) は、原子や分子が物理学の基本的な法則に従って相互作用するコンピュータシミュレーションです。 GROMACS の |参照マニュアル <md>| は、この分野に関する包括的な入門書として、また GROMACS を使用するための具体的な情報を提供しています。 GROMACS を使用したい場合は、まず数章を読むことが必須であり、時間を無駄にすることなく、効率的に学習を進めることができます。

  • 導入:分子モデリング(スライド、ビデオ) - 理論的枠組み、モデリングレベル、限界と可能性、システムと方法(エリック・リンダール)

書籍

周囲にはいくつかの教科書があります。

優れた入門書としては、以下のようなものがあります。

  • A. リーチ (2001) 分子モデリング:原理と応用。

  • T. シュリック (2002) 分子モデリングとシミュレーション

プログラミングの知識がある場合:

  • D. ラパポルト (1996) 分子動力学シミュレーションの技法

  • D. フレネル、B. スミス (2001) 分子シミュレーションの理解

より物理学者の視点から:

  • M. アレン、D. ティルデレー (1989) 液体コンピュータシミュレーション

  • H.J.C. Berendsen (2007) 量子力学から流体ダイナミクスまでの階層的モデリング:物理世界のシミュレーション

種類 / 集合

  • NVE - 粒子の数 (N)、系の体積 (V) とエネルギー (E) は一定 / 保全される。

  • NVT(粒子数(N)、システム体積(V)および温度(T))は一定/保存されます。 (「恒温」に関する詳細は、thermostats を参照してください。)

  • NPT - 粒子数 (N)、システム圧力 (P) および温度 (T) は一定 / 保恒されます。 (「一定の圧力」に関する詳細は、pressure coupling を参照してください)。

強制場

強制場は、物理システムを研究するために使用できる、潜在関数とパラメータ化された相互作用の集合です。これらに関する一般的な歴史、機能、および使用方法に関する情報は、このガイドの範囲外であり、ユーザーは関連する文献を参照するか、関連する Wikipedia ページから始めることを推奨します。