天天看點

【Amber】帶小分子配體分子動力學模拟

目錄

      • 适用于Amber18\20版本
      • 一、檔案準備
        • 1.檔案檢查與拆分
        • 2.小分子檔案處理
        • 3.蛋白質檔案處理
        • 4.複合體檔案處理
      • 二、分子動力學模拟
        • 1.能量最小化
        • 2.體系加熱
        • 3.恒壓平衡
        • 4.全局平衡
        • 5.正式模拟

适用于Amber18\20版本

使用Amber做帶小分子配體的MD時,由于Amber無法識别小分子,需要對小分子做處理,相比于單純蛋白分子的MD步驟會繁瑣一些,在此記錄一下。

一、檔案準備

1.檔案檢查與拆分

含配體的複合體pdb檔案可以來自晶體或對接結果,其中小分子與蛋白分子的相對位置将會是MD的初始位置,在準備過程中需要檢查确認。

【1】首先在可視化軟體如薛定谔等檢查複合體中小分子,確定其鍵連準确,沒有原子,特别是氫原子的缺失。

【2】檢查完畢後儲存複合體檔案,使用文本編輯器打開,推薦使用

notepad++

,找到對應小分子的行,将其剪切到新的文本中,命名為

Y.pdb

,剩餘的蛋白分子命名為

X.pdb

2.小分子檔案處理

【3】使用Amber的工具antechamber将小分子pdb檔案轉為mol2檔案

如果小分子配體較為複雜,此處會耗費較長的時間。

注意:此處使用bcc計算分子的電荷分布,其準确性低于RESP,如果懂得使用gaussian或Multiwfn,建議使用RESP計算的電荷

【4】得到

Y.mol2

檔案後,使用parmchk2産生小分子參數檔案

parmchk2 -i Y.mol2 -f mol2 -o Y.frcmod
           

【5】建立檔案

Y.tleap

,寫入以下内容

source leaprc.gaff
Y = loadmol2 Y.mol2 
check Y
loadamberparams Y.frcmod
check Y
saveoff Y Y.lib 
saveamberparm Y Y.prmtop Y.inpcrd
quit
           

【6】使用

tleap

程式,生成小分子

Y.lib

檔案

tleap -f Y.tleap
           

至此得到

Y.frcmod

Y.lib

Y.prmtop

Y.inpcrd

檔案。

這裡碰到一個問題,在處理帶有磷酸基團的小分子時,如果為磷酸基團添加氫原子,在MD過程中,磷酸基團的氫原子會快速無序擺動,導緻整個模拟體系崩潰,在體系加熱過程中的展現為溫度大範圍波動,并遠遠超出設定的上限,報錯顯示浮點超出上限或非法記憶體通路,去除磷酸基團的氫原子後一切恢複正常,報錯資訊如下。

Error: an illegal memory access was encountered launching kernel kNLSkinTest
           

3.蛋白質檔案處理

【7】使用文本編輯器打開

X.pdb

檔案,删掉除原子資訊外的所有行,隻保留如下的行

...
...
ATOM  16100  N   LEU B 537      18.387 -77.040 -41.003  1.00  0.80           N  
ATOM  16101  CA  LEU B 537      19.579 -76.228 -41.346  1.00  0.80           C  
ATOM  16102  C   LEU B 537      20.327 -75.718 -40.033  1.00  0.80           C  
ATOM  16103  O   LEU B 537      20.000 -76.288 -38.963  1.00  0.80           O  
ATOM  16104  CB  LEU B 537      19.148 -74.980 -42.190  1.00  0.80           C  
ATOM  16105  CG  LEU B 537      18.636 -75.198 -43.626  1.00  0.80           C  
ATOM  16106  CD1 LEU B 537      18.176 -73.876 -44.243  1.00  0.80           C  
ATOM  16107  CD2 LEU B 537      19.672 -75.869 -44.530  1.00  0.80           C  
ATOM  16108  OXT LEU B 537      21.005 -74.658 -40.138  1.00  0.80           O1-
ATOM  16109  H   LEU B 537      18.021 -76.812 -40.085  1.00  0.00           H  
ATOM  16110  HA  LEU B 537      20.334 -76.831 -41.849  1.00  0.00           H  
ATOM  16111  HB3 LEU B 537      20.001 -74.318 -42.301  1.00  0.00           H  
ATOM  16112  HB2 LEU B 537      18.389 -74.435 -41.628  1.00  0.00           H  
ATOM  16113  HG  LEU B 537      17.759 -75.840 -43.564  1.00  0.00           H  
ATOM  16114 HD11 LEU B 537      17.370 -74.053 -44.953  1.00  0.00           H  
ATOM  16115 HD12 LEU B 537      17.814 -73.195 -43.476  1.00  0.00           H  
ATOM  16116 HD13 LEU B 537      18.989 -73.377 -44.769  1.00  0.00           H  
ATOM  16117 HD21 LEU B 537      19.253 -76.103 -45.506  1.00  0.00           H  
ATOM  16118 HD22 LEU B 537      20.558 -75.250 -44.669  1.00  0.00           H  
ATOM  16119 HD23 LEU B 537      19.998 -76.804 -44.088  1.00  0.00           H  
TER   16120      LEU B 537
           

【8】使用Amber的工具

pdb for amber

删除模型内氫原子

pdb4amber -i X.pdb -o X_noH.pdb -y --dry
           

【9】以Amber的氫命名方式重新添加氫原子

reduce X_noH.pdb > X_H.pdb
           

【10】使用pdb4amber對加氫後的pdb檔案進行處理

pdb4amber -i X_H.pdb -o X.pdb
           

至此得到Amber标準的蛋白質分子檔案

X.pdb

4.複合體檔案處理

【11】将

Y.pdb

使用文本編輯器打開,将原子資訊複制到

X.pdb

檔案後,重命名為

com.pdb

【12】建立新檔案

com.tleap

寫入以下内容,儲存

source leaprc.protein.ff14SB `此處可以使用ff19SB` 
source leaprc.water.tip3p
source leaprc.gaff
loadamberparams Y.frcmod
loadoff Y.lib
complex = loadpdb com.pdb
check complex
solvatebox complex TIP3PBOX 12.0
addions2 complex Na+ 0
addions2 complex Cl- 0
saveamberparm complex X_box.prmtop X_box.inpcrd
savepdb complex X_box.pdb
quit
           

此處可以使用ff19SB

【13】使用

tleap

程式,生成複合體參數檔案

X_box.prmtop

X_box.inpcrd

檔案以及水盒檔案

X_box.pdb

tleap -f com.tleap
           

至此,得到

X_box.prmtop

X_box.inpcrd

X_box.pdb

二、分子動力學模拟

1.能量最小化

【1】限制主鍊最小化

min1.in

&cntrl  
  imin=1, 
  maxcyc=10000, 
  ncyc=5000, 
  ntb=1, 
  ntr=1,
  restraintmask=':1-1104', 
  restraint_wt=2,
  cut=8.0 
 / 
 END
           

執行指令

pmemd.cuda -O -i min1.in -o X_box_min1.out -p X_box.prmtop -c X_box.inpcrd -r X_box_min1.rst -ref X_box.inpcrd
           

【2】無限制最小化

min2.in

&cntrl  
  imin=1, 
  maxcyc=100000, 
  ncyc=5000, 
  ntb=1, 
  ntc=1,
  ntf=1, 
  ntpr=10,
  cut=8.0 
 / 
 END
           

執行指令

pmemd.cuda -O -i min2.in -o X_box_min2.out -p X_box.prmtop -c X_box_min1.rst -r X_box_min2.rst
           

2.體系加熱

【3】限制主鍊恒容50ps加熱

heat.in

explicit solvent initial heating: 50ps
 &cntrl
  imin=0,
  irest=0,
  nstlim=25000, dt=0.002,
  ntc=2, ntf=2, ntx=1,
  cut=8.0, ntb=1,
  ntpr=500, ntwx=500,
  ntt=3, gamma_ln=2.0,
  tempi=0.0, temp0=300.0, ig=-1,
  ntr=1,
  restraintmask=':1-1104',
  restraint_wt=2.0,
  iwrap=1
  nmropt=1
  /
  &wt TYPE='TEMP0', ISTEP1=0, ISTEP2=25000,
  VALUE1=0.0, VALUE2=300.0, /
  &wt TYPE = 'END' / 
 END
           

執行指令

pmemd.cuda -O -i heat.in -o X_box_heat.out -p X_box.prmtop -c X_box_min2.rst -r X_box_heat.rst -x X_box_heat.nc -ref X_box_min2.rst
           

3.恒壓平衡

【4】50ps平衡

press.in

explicit solvent density: 50ps
 &cntrl
  imin=0,
  irest=1,
  ntx=5,
  nstlim=25000, dt=0.002,
  ntc=2, ntf=2,
  cut=8.0, ntb=2, ntp=1, taup=2.0,
  ntpr=500, ntwx=500,
  ntt=3, gamma_ln=2.0,
  temp0=300.0, ig=-1,
  ntr=0,
  / 
 END
           

執行指令

pmemd.cuda -O -i press.in -o X_box_press.out -p X_box.prmtop -c X_box_heat.rst -r X_box_press.rst -x X_box_press.nc -ref X_box_heat.rst
           

4.全局平衡

【5】10ns平衡

eq.in

&cntrl
  imin=0, irest=1, ntx=5,
  nstlim=5000000, dt=0.002,
  ntc=2, ntf=2,
  cut=10.0, ntb=2, ntp=1, taup=2.0,
  ntpr=500, ntwx=500, ntwr=5000,
  ntt=3, gamma_ln=2.0,
  temp0=300.0,
 / 
 END
           

執行指令

pmemd.cuda -O -i eq.in -o X_box_eq.out -p X_box.prmtop -c X_box_press.rst -r X_box_eq.rst -x X_box_eq.nc
           

5.正式模拟

【6】100ns模拟

md.in

explicit solvent production run: 100ns
 &cntrl
  imin=0,
  irest=1,
  ntx=5,
  nstlim=50000000, dt=0.002,
  ntc=2, ntf=2,
  cut=8.0, ntb=2, ntp=1, taup=2.0,
  ntpr=5000, ntwx=5000, ntwr=50000,
  ntt=3, gamma_ln=2.0,
  temp0=300.0, ig=-1,
  iwrap=1
  / 
 END
           

執行指令

pmemd.cuda -O -i md.in -o X_box_md.out -p X_box.prmtop -c X_box_eq.rst -r X_box_md.rst -x X_box_md.nc
           

繼續閱讀