AUTODOCK＞＞PYMOL
距離を測る・点線を描く
- Wizard → Measurement
- 距離を図りたい2原子を続けて左クリックすると、距離（単位はÅ）と結合（点線）が表示されます。 また、そのままクリックすることで続けて距離を測ることができます。
- 終わったら、オブジェクトパネルの下の "Done" をクリック
ラベルの調整
- ラベルを消す オブジェクトパネルの "measure01" のHide → labels
- ラベルのサイズや位置を調整する
- をご参照ください。
点線の調整
- 色の変更
- 赤字に色名を入力
- 太さの変更
- 赤字に任意の値を入力
- 点線の長さの変更
- 赤字に任意の値を入力
- 点線のギャップの変更
- 赤字に任意の値を入力
【PyMol+Chimera+AutoDock Vina】フリーソフトで雰囲気SBVSをやった時のメモ
1. SBVSってなに？
Structure-based virtual screeningの頭文字をとって(SBVS). Structure-based drug design (SBDD)の１分野で, 受容体タンパク質(薬の標的タンパク質)側の構造情報を元に親和性の高い薬をスクリーニングする手法です.
標的タンパク質を”鍵穴”, 医薬品候補化合物(リガンド)を”鍵”とした鍵と鍵穴問題として例えられる医薬品開発において, 鍵穴の情報であるタンパク質の立体構造を顕に扱うSBDDは創薬研究者にとって有用な数々の情報をもたらすそうです.
SBVSにおいて特に重要な役割を果たす手法は, 計算機によって標的タンパク質とリガンドの立体的な結合構造と結合親和性を予測するドッキング計算です. 例えば新型コロナウイルス(COVID-19)の治療薬候補を絞り込んだ時もこのドッキング計算が用いられました.
[タンパク質－リガンドドッキングの現状と課題[1]](https://www.jstage.jst.go.jp/article/cicsj/34/1/34_10/_article/-char/ja)によると, ドッキングは以下のような利点があるそうです. >ドッキング計算によってリード化合物の結合構造を正しく特定することができれば，標的タンパク質に対するリガンドの相互作用機序を*原子レベルで理解*することができる．さらにこれらの立体的な結合構造の情報を利用することによって， *医薬品候補化合物の結合親和性や選択性向上のための合理的な置換基の設計*も可能となる．
タンパク質は多数の分子で構成されているため, その立体構造をコンピュータにより計算しようとすると多くの計算コストを要します. そのため, ドッキング計算には以下のような簡略化のタイプがあるようです.
- Rigid Docking: 鍵と鍵穴モデル(Lock-and-Key Model)に従い, タンパク, リガンド共に剛体として扱う.
- Flexible Docking: リガンド側のみ柔軟性を考慮（タンパク質は剛体, Flexible Ligand Model）
- Flexible Side-Chains Docking: リガンドに加えて結合ポケットのアミノ酸側鎖の自由度も考慮（Paritally Flexible Protein Model）
- Flexible Protein Docking: タンパク、リガンド共に柔軟性を考慮（Induced-Fit Model）
2. スコア関数
ではどうやって最適なドッキング状態の計算を行うのか? ドッキング計算の手法には様々なアプローチが存在しますが, 共通する2つの目的は標的タンパク質に対するリガンドの結合構造と, 標的タンパク質に対するリガンドの結合親和性を高速に予測することです． この目的を達成するために, ドッキングプログラムには高精度なスコア関数と効率的な探査アルゴリズムが必要とされます.
力場に基づいたスコア関数では, タンパク質―リガンド原子間の相互作用はファンデルワールス相互作用 (vdW) や静電相互作用, これらに溶媒効果やリガンドの内部ポテンシャルの変化を記述したものによって表されるそうです. 以下のような物理化学的相互作用を考慮します.
Δ�=Δ����+Δ��������������+Δ��−����+Δ������������+Δ�����
- Δ�:自由エネルギーの変化量
- Δ����: ファンデルワールス相互作用
- Δ��������������: 静電相互作用
- Δ��−����: 水素結合による相互作用
- Δ������������: 脱溶媒効果
- Δ�����:リガンドのエントロピー変化の項
Δ����,Δ��������������,Δ��−����は直感的にも納得できるエネルギー項だと思いますが, Δ������������とΔ�����はなんでしょうか？ リガンドがタンパク質に結合するにあたって, 溶媒（水）に囲まれた環境から, 溶媒が外れ（脱溶媒）タンパク質の結合サイト（より疎水性の高い環境）に移行するという過程をたどると考えます. この時のエネルギー変化を考慮するためのものがΔ������������. また, 溶液中では自由な立体構造（配座、コンフォメーション）をとっているリガンドが, 結合に際してがポケットに合うように固定化されこの自由度が失われます. この分のエントロピー変化を考慮するものがΔ�����となっているようです.
これらのエネルギー変化量を考え, Δ�が負に大きいとき, それがより安定な構造であると考えられます.
それぞれの項のイメージはWhat is Docking?にわかりやすくまとまっております.
スコア関数は他にも色々あるらしいが詳しくは触れない.
3. Docking Study
3.1. 準備 (Macの場合)
- PyMol: 使えればなんでも問題ありませんが, 自分はGitHubからcloneしてビルド＆インストールをした, オープンソースのPyMolを使用しました.インストール方法はGitHubとここに従えばできると思います.
- Chimera: ドッキングの前処理に使います.ここからダウンロードしてインストール
- AutoDock Vina: ドッキングに用います. ここからダウンロード
- SARS-Cov-1のMain proteaseの共結晶構造: このタンパク質(PDB ID: 2GZ7)は2021年現在を世界を大きく騒がせているSARS-Cov-2(通称COVID-19)のMain proteaseと相同性が非常に高い. 右上のボタンからpmdファイルをダウンロード
3.2. ドッキングに向けての前処理 I 【PyMol】
3.2.1. 共結晶構造から標的タンパク質とリガンドを単離する
使うのはPyMolとダウンロードしてきた2gz7.pdbファイル.
前処理のベストプラクティスは正直よく知らない. 今回は共結晶構造から水分子を取り除き, 受容体(標的タンパク質)とリガンドに分割してそれぞれpdbファイルとして保存した.
- PyMolのコマンドラインに
remove solventを入力して水を除去 - 水除去後の構造
- リガンドをクリックしてremove→pdbファイルとして受容体を保存
- 「File」→「Reinitialize」→「Everything」してもう一度2gz7.pdbを読み込む. 再度
remove solventで水を除去し, リガンドをクリック→「action」→「invert」→「with object」で選択を反転 - remove
- リガンドをpdbファイルとして保存
3.2.2. ドッキングの中心座標を取得する
共結晶構造のリガンドがはまっていた座標をドッキングの中心座標にしたい.
リガンドのpdbファイルを開く.
リガンドをクリック
以下のコマンドをうつ
pseudoatom center_, sele
xyz = cmd.get_coords('center_', 1)
print xyz
↓実行結果↓
[[-24.638144 -41.232357 5.4246078]]
この座標を覚えておく.
3.3. ドッキングに向けての前処理 II 【Chimera】
単離した受容体とリガンドに水素を付加し, 電荷を割り当てる
このチャプターではChimeraを用います.
- Chimera起動
- 「File」→「Open」→先ほど作成した受容体のpdbファイルを開く
- 「Tools」→「Dock Prep」
- 「Delete solvent」, 「Add hydrogens」, 「Add charges」 にチェックが入っていることを確認してMol2ファイルに書き出します
- Add Hydrogens for Dock Prep とかいうウィンドウが出てくるのでOK
- Assign Charges for Dock Prep とかいうウィンドウが出てくるのでOKを押してmol2ファイルを生成します.
- リガンドのpdbファイルに関しても同様に水の除去と水素付加, 電荷の割り当てを行い, mol2ファイルを生成します.
3.4. ドッキングシミュレーション
準備は整いましたのでドッキングを行います.
一旦Chimeraの「File」→「Close Session」によりセッションを終了し, 先ほど前処理した受容体のmol2ファイルとリガンドのmol2ファイルの両方をChimeraで開きます.
以下のような画面になります.
「Tools」→「Surface/Binding Analysis」→「AutoDock Vina」を選択し, 以下の画面のように設定します.
要確認事項
- Receptorには前処理した受容体が選択されているか？
- Ligandには前処理したリガンドが選択されているか？
- Centerには先ほどPyMolで出力した座標が入っているか？
- Localのチェックボックスにチェックを入れ, AutoDock Vina由来の/bin/vinaまでのパスが指定されているか？
- Out fileに出力ファイルの名前が入っているか？
設定が完了したらOKを押すとドッキングが始まり, pdbqtファイルが生成されます
3.5. 結果の解析
3.5.1 結果の解析
画像のような感じで結果が出てきます.
Scoreはエネルギースコアを指しており, これが低い方が好ましい状態です. Root Mean Square Deviation(RMSD)はタンパク質構造の非類似性や誤りの指標として用いられるものです.
さまざまなドッキングポーズでのエネルギースコアを算出しましたが, スコアが小さい順に並んでいます. この状態で結果を保存します.保存したドッキング結果はChimeraの「Tools」→「Surface/Binding Analysis」→「View Dock」から確認できます.
| ID | Score | RMSD |
|---|---|---|
| 1 | -7.2 | 0.0 |
| 2 | -7.1 | 1.523 |
| 3 | -7.0 | 4.952 |
| 4 | -6.8 | 3.971 |
| 5 | -6.8 | 4.844 |
| 6 | -6.8 | 5.632 |
| 7 | -6.8 | 4.937 |
| 8 | -6.7 | 3.923 |
| 9 | -6.7 | 5.345 |
| という結果でした. |
3.5.2 ドッキングポーズをPyMolで詳しく見る
ドッキングした後に生成されたリガンド・受容体・結果のpdbqtファイルを全てPyMolで読み込み, 以下のように確認することができます. 赤がリガンド, 緑がドッキングで計算したドッキングポーズです. 結果のpdbqtファイル(緑)には上の９種のドッキングポーズが格納されています. 右側から「A」→「State」→「split」で全てのドッキングポーズを表示することができます. 
分割したポーズはスコアが低い順に並んでいました.
それぞれのドッキングポーズをgifにしてみました.
3.5. まとめ
SARS-Cov-1のMain proteaseの共結晶構造を用いてドッキングの動作テストを行いました. 使用したソフトは全て無料なので取り掛かりやすいと思います. 今回はドッキングの雰囲気は分かったかなという感触です. しかし実際の活用の面でこのような単純な解析で済むのかという面で複数の疑問が生じました.
- タンパク質の立体構造の最適化は行わないのか？
- 溶媒を除いたが, 生体内ではタンパク質は水に覆われている. 溶媒の影響は考慮しないのか？
- タンパク質のX線結晶構造解析の結果が得られていなかった場合ドッキングは不可能なのか？
この辺に関してまた新たな情報が入ったら記事を作ろうと思います.
計算の実施
計算の流れは大まかに以下のようになる。必要なツールを括弧内に記載いた。
- 計算対象のリガンドおよびタンパク質をPDB形式で書き出す (Pymol)
- 計算対象のタンパク質とリガンドに前処理を実施 (MGL tools)
- .confファイルを作成し計算オプションを指定 (Chimera, テキストエディタ)
- ドッキングを実施 (AutoDock Vina)
(Pymol) 計算対象のリガンドおよびタンパク質をPDB形式で書き出す
ドッキング計算では計算結果と結晶構造のリガンドを重ね合わせることで計算手法の妥当性を検討することがとても重要である。
今回はリガンドの結晶構造が既知のSARS-Cov-1のMain protease(PDB ID: 2GZ7)をターゲットタンパク質とし, 計算結果と結晶構造との比較を行う。このタンパク質は2020年現在を世界を大きく騒がせているSARS-Cov-2(通称COVID-19)のMain proteaseと相同性が非常に高く興味深い対象である。
以下に結合サイトをPymolで図示した。リガンドが結合サイト内のアミノ酸残基および水分子と水素結合を形成している。このことは水分子がリガンドの結合に重要な役割を果たしている可能性を示唆するが, 今回は計算を簡略化するために水分子は全て取り除いた。
PyMOL>remove solvent上のコマンドで水を取り除き, タンパク質およびリガンドをそれぞれligand.pdbおよびprotein.pdbというファイル名で保存した。詳しいPymolの操作方法は別のサイトを御覧いただきたい。
(MGL Tools) 計算対象のタンパク質とリガンドに前処理を実施
前述のリガンドおよびレセプターのPDBファイルをMGL Toolsで処理する。MGL ToolsはPython 2で書かれているためPython3だと動かないので注意。自分はanacondaでpython2.7の環境を作成してその中で動かした。
まず標的リガンドであるligand.pdbを以下のコマンドで処理する。”-A bonds_hydrogen”オプションを追加することでMGL toolsが自動的にリガンドの結合情報を変更し水素原子を付加してくれる。処理後はligand.pdbqtというファイルが生成され, これをAutodock Vinaのインプットとして用いる。
/(Path_to_MGLtools)/mgltools_i86Darwin9_1.5.6/MGLToolsPckgs/AutoDockTools/Utilities24/prepare_ligand4.py \
-l ligand.pdb -A bonds_hydrogen 次に標的タンパク質であるprotein.pdbも同様に処理する。得られたprotein.pdbqtをAutodock Vinaのインプットとして用いる。
/(Path_to_MGLtools)/mgltools_i86Darwin9_1.5.6/MGLToolsPckgs/AutoDockTools/Utilities24/prepare_receptor4.py \
-r protein.pdb -A bonds_hydrogen AutoDock Vinaは計算時に電荷情報を自動的に追加するため, インプットファイルで電荷情報を追加する必要はない(“Why don’t the results change when I change the partial charges?”参照)。
また, AutoDock VinaはUnited Atom modelを採用しているため, Non-polar hydrogenはAutoDock Vinaの計算実行時に全て削除される。一方でAutoDock Vinaは自動的に水素原子を付加してくれるわけではないので, polar hydrogenの位置を把握するためにも入力ファイルには水素原子の結合情報を与える必要がある。
(Chimera, テキストエディタ) .confファイルを作成し計算オプションを指定
AutoDock Vinaでは計算時にポケットの位置を指定する必要がある。ポケットの位置は.confファイル内に記入する。
まずChimeraでligand.pdbおよびprotein.pdbを開く。MGL Toolsで作成したpdbqtファイルは読み込めないので注意。次にTools -> AutoDock Vinaからポケットの位置を指定する。ligand.pdbファイル内の原子の座標を参考に”Receptor search volume options”からボックスの大きさおよび中心を指定し, ポケットが十分覆われるようにする。
前述のボックスの座標情報を元にAutoDock Vinaの計算に必要な.confファイルをテキストエディタで作成する。以下に作成したautodock.confの例を示す。
center_x = -28.00
center_y = 38.00
center_z = 4.00
size_x = 20.00
size_y = 20.00
size_z = 20.00
energy_range = 3
exhaustiveness = 8
num_modes = 9その他のオプションについて説明する。いずれの値もAutoDock Vinaのデフォルト値に設定した。
energy_rangeはベストスコアから何kcal/molまでの結果を出力するかを意味する。具体的にはenergy_range = 3ではベストスコアが-9.0 kcal/molの場合, -6.0 kcal/molまでの結果が出力される。
exhaustivenessはエネルギー地形をどれだけ徹底的に探索するかを示している(詳細)。大きくするほど計算時間がかかるがグローバルミニマムを逃す確率が下がる。
num_modesは最大の出力結果数を示す。9なら9つのモードが出力される。
(AutoDock Vina) ドッキングを実施
ついにドッキングを実施する。先ほど作成した２つの.pdbqtファイルおよびconfファイルをLinuxに転送。以下のコマンドでシミュレーションを実施する。
/(path_to_vina)/vina --receptor receptor.pdbqt --ligand ligand.pdbqt \
--config autodock.conf --out autodock.pdbqt --log autodock.log計算はAMD Ryzen7 3800Xの32並列で約10秒で終了。はやい。以下のRMSDはベストスコアのもの(1番目のmode）に対する相対値として計算されるもので, 結晶構造に対するRMSDではないので注意。
mode | affinity | dist from best mode
| (kcal/mol) | rmsd l.b.| rmsd u.b.
-----+------------+----------+----------
1 -7.1 0.000 0.000
2 -7.0 1.554 2.642
3 -7.0 4.893 6.867
4 -6.9 4.420 6.587
5 -6.9 1.526 3.694
6 -6.8 3.932 5.225
7 -6.8 4.876 7.154
8 -6.7 4.940 6.864
9 -6.7 5.309 7.497Code language: JavaScript (javascript)結果のautodock.pdbqtはPymolで可視化できる。
2番目と6番目のposeが結晶構造をよく再現した一方でそれ以外の構造は差異が大きい。6番目の結果にも関わらず結晶構造との再現度が高い辺り, スコアリング順が良い = ドッキング構造が良い, では必ずしもないことがわかる(参考）。
とはいえ正解のドッキングポーズが結果にきちんと含まれていることから, 今回の系に関してはin sillicoのアプローチが有効である可能性がある。今回の標的タンパク質と相同性が高いSARS-Cov-2のMpro阻害薬探索においてもドッキング計算は一定の有効性がありそうだ。
2024年12月10日 | カテゴリー：AUTODOCK VINA |
