クイックスタート: pdb2reaction all

pdb2reaction all は、反応物(R)と生成物(P)の 2 つの構造から、1 回の実行で反応経路を作ります。基質のまわりのクラスターモデルを切り出し、R と P の間の最小エネルギー経路(MEP)を探索します。--tsopt --thermo --dft を付けると、同じ実行のまま遷移状態(TS)の最適化・固有反応座標(IRC)の計算・振動数・DFT 一点計算まで進みます。反応座標の定義が手間なときに向いています。中間体や生成物の構造を PyMOL や GaussView で自分で作って入力できます。反応座標を自分で決めずに MEP を探索するので、新しい反応機構の候補が見つかることもあります。

以下のコマンドは、examples/ にある、ゲラニル二リン酸(GPP)の C6 位をメチル化する酵素 BezA の同梱例を使います。1.R.pdb が反応物、3.P.pdb が生成物です。同梱例は git clone https://github.com/t-0hmura/pdb2reaction && cd pdb2reaction/examples で取得できます。自分の反応では、全系の構造に置き換えてください。


主な用途

  • 全工程を初めて通す: 同梱例で、すべての段を 1 回実行

  • R と P の間の MEP を作る: 経路と、その最高エネルギーのイメージ(HEI、TS の候補)を取得

  • 同じ実行で TS・IRC・振動数・DFT まで進める: --tsopt --thermo --dft を付けて TS の候補を確認

最小コマンド

R と P を反応の順に渡し、切り出しの中心にする残基(-c)とリガンドの電荷(-l)を指定します。

pdb2reaction all -i 1.R.pdb 3.P.pdb -c 'SAM,GPP,MG' -l 'SAM:1,GPP:-3' \
 --out-dir ./result_all

端末の最後のほうの ====== Pipeline summary ====== の下に Scientific status: success と出れば成功で、summary.json の scientific_status にも同じ値が入ります。

(任意)同一実行で後処理まで行う

--tsopt で反応セグメント(ここでは seg_01)ごとの TS 最適化と IRC を、--thermo で振動数と熱化学を、--dft で R・TS・P の DFT 一点計算を追加します。--thermo と --dft には --tsopt が必要です。

pdb2reaction all -i 1.R.pdb 3.P.pdb -c 'SAM,GPP,MG' -l 'SAM:1,GPP:-3' \
 --tsopt --thermo --dft --out-dir ./result_all

実行の前に

構造にはすべての水素原子が要り、R と P は同じ原子を同じ順に並べている必要があります。詳しくは 複合体の構造を用意する を参照してください。

主な出力ファイル

最小コマンドは次のファイルを書き出します。

result_all/
├── summary.log                  # 実行の要約
├── summary.json                 # 結果(scientific_status を含む)
├── mep_trj.pdb                  # 全セグメントの MEP
├── energy_diagram_MEP.png       # 全セグメントの MEP のエネルギープロファイル
└── _work/                       # 途中のファイル(TS 候補の HEI を含む。実行後も残る)
    └── path_opt/                # MEP 探索(MEP を再帰的に詰める --refine-path のときは path_search/)
        ├── hei_seg_01.{xyz,pdb} # セグメント 1 の最高エネルギーのイメージ
        └── summary.json         # MEP 探索の結果

最小コマンドは MEP 探索で終わるため、segments/ は作られません。--tsopt を付けると、反応セグメントごとに segments/seg_NN/ ができ、R/TS/P の構造 reactant.pdb・ts.pdb・product.pdb と、ts/、irc/ が入ります。--thermo を付けると freq/ も加わります。

結果の確認

  1. 完了状況: scientific_status には、指定した段がすべて収束すると success、そうでなければ partial か failed が入り、理由は scientific_status_reasons に出ます。--tsopt のとき、虚振動のモードが、できる結合と切れる結合を動かしているか、端点が狙った R と P か、の 2 点は自分で確かめてください。

  2. TS の候補: 最初のセグメントの HEI _work/path_opt/hei_seg_01.pdb を開きます。--tsopt のときは、最適化した TS の segments/seg_01/ts.pdb も開きます。

  3. エネルギープロファイル: energy_diagram_MEP.png で、R と P の間にはっきりした障壁があるかを確かめます。

  4. TS(--tsopt のとき): TS 最適化が成功すると、反応モードの虚振動が 1 つ出ます。このとき端末に [tsopt] Converged (n_imag=1). と出て、summary.json の post_segments[].tsopt.n_imaginary_modes に本数が記録されます。segments/seg_01/ts/vib/imag_*_trj.xyz をビューアで開き、できる結合と切れる結合に沿って原子が動くかを確認してください。

  5. 端点(--tsopt のとき): segments/seg_01/irc/finished_irc_trj.xyz と、最適化した端点の segments/seg_01/reactant.pdb・product.pdb を開き、狙った R と P かを確かめます。IRC が収束しなくても、端点の最適化で狙った R と P に着けば、その結果は使えます。

all が各段をどう判定するかは 実行結果の判定 を参照してください。

使用上の注意点

  • DFT と GPU のメモリ: --dft には DFT 用の追加パッケージが要ります。その導入と GPU メモリについては MLIP の TS を DFT で確かめる の使用上の注意点を参照してください。

  • summary.json の障壁: segments[].barrier_kcal は TS 最適化の前の、MEP の上の障壁です。--tsopt を付けると、最適化した TS と端点から求めた障壁が post_segments[].mlip.barrier_kcal に入り、--thermo で post_segments[].gibbs_mlip.barrier_kcal、--dft で post_segments[].dft.barrier_kcal が加わります。

  • rate_limiting_step: rate_limiting_step.barrier_kcal は、すべてのセグメントにそろっている最も高いレベル(DFT//MLIP_Gibbs > DFT > MLIP_Gibbs > MLIP > MEP)で比べた、最も高い障壁です。使ったレベルは rate_limiting_step.method に入ります。

次のステップ