此前,何博士為大家在北鯤云的直播間分享了Amber熱力學積分計算相對自由能變化(直(zhi)播回放(fang)可在視頻號:北鯤云-直播回放中查看)
直播結束(shu)后有很多小伙伴(ban)來向我們要PPT資料(liao),這(zhe)里何博(bo)士也為大家準備了(le)文(wen)字版本(ben)的教程(cheng)。將為大家介紹(shao)如何在北鯤云計(ji)算(suan)平臺(tai)上利用(yong)Amber熱力學積分(fen)計(ji)算(suan)相對自由能變化,體系包(bao)括小分(fen)子-蛋(dan)白(小分(fen)子改變),小分(fen)子-蛋(dan)白(蛋(dan)白突變),蛋(dan)白-蛋(dan)白相互作用(yong)。
本教程(cheng)要求(qiu)使用者一定程(cheng)度了解Amber動力學模擬程(cheng)序。
Amber是美國(guo)加州(zhou)大學(xue)Peter Kollman等(deng)開發的(de)一款著名的(de)分子(zi)動力學(xue)模(mo)擬軟件(jian)包。Amber主要(yao)適(shi)用(yong)于蛋白質,小分子(zi)和多糖等(deng)生(sheng)物分子(zi)體系的(de)模(mo)擬。
本文所需的(de)所有文件(jian)請在//github.com/Xinheng-He/ti_toturial上(shang)下載。

該應用場(chang)景解決將蛋白口袋內的小分子(zi)A變為小分子(zi)B所產生(sheng)的相對自由能變。

將蛋白口袋(dai)內的苯(ben)轉化為苯(ben)酚(fen)。
首先,使用pymol將分(fen)子打開(kai),并選中(zhong)小分(fen)子,保存為mol2文件(jian),如(ru)下圖所示,我們使用
save *my_case\ben_ligand.mol2
save *\my_case\benfen_ligand.mol2
這(zhe)兩個(ge)命令保存了(le)變(bian)化前后的配(pei)體。

(保(bao)存pymol中的sele對象)
開啟一(yi)個北鯤(kun)云(yun)管理(li)節點加載環境。
module add Anaconda3/2020.02
source/public/software/.local/easybuild/software/amber/aber20/amber.sh
ulimit-s unlimited
ulimit-l unlimited
對苯生成具有電荷(he)(he)的可用mol2文(wen)件(jian),總(zong)電荷(he)(he)為0,殘基(ji)名為BEN。
antechamber-i ben_ligand.mol2 -fi mol2 -o ben_real.gaff2.mol2 -fo mol2
生成frcmod力場(chang)參數文(wen)件(jian)。
parmchk2-i ben_real.gaff2.mol2 -f mol2 -o BEN.gaff2.frcmod-s gaff2-a yes
上述操作對苯酚再來一次。
antechamber -i benfen_ligand.mol2 -fi mol2 -o benfen_real.gaff2.mol2 -fo mol2 -rn FEN-atgaff2-anyes-drno-pfyes-cbcc-nc0
parmchk2 -i benfen_real.gaff2.mol2 -f mol2 -o FEN.gaff2.frcmod -s gaff2 -a yes
根據frcmod文(wen)件,生成兩(liang)個分(fen)子(zi)的(de)文(wen)庫文(wen)件,該文(wen)件描述了分(fen)子(zi)內部的(de)原子(zi)類型和鍵連信息。
tleap-f-<<_EOF
source leaprc.gaff2
loadamberparams FEN.gaff2.frcmod
FEN = loadmol2 benfen_real.gaff2.mol2
saveoff FEN FEN.lib
savepdb FEN FEN.pdb
quit
_EOF
tleap -f - <<_EOF
source leaprc.gaff2
loadamberparams BEN.gaff2.frcmod
BEN = loadmol2 ben_real.gaff2.mol2
saveoff BEN BEN.lib
savepdb BEN BEN.pdb
quit
_EOF
注意前后的力場(chang)要保持一(yi)致。
將(jiang)兩(liang)個(ge)pdb文件(FEN.pdb和BEN.pdb)中同樣(yang)位(wei)置的(de)(de)全(quan)部重原子,保存成同樣(yang)的(de)(de)坐標,注意名字要和lib中的(de)(de)一(yi)樣(yang),放成一(yi)個(ge)lig.pdb,在下面的(de)(de)tleap過程中,tleap會自動根據(ju)lib文件將(jiang)complex中的(de)(de)原子變成真(zhen)實的(de)(de)樣(yang)子,這樣(yang)做是為了保證一(yi)樣(yang)原子的(de)(de)位(wei)置完全(quan)一(yi)致,減(jian)少不必要的(de)(de)變量。

(pdb文件的內容)
使用pdb4amber,檢查蛋白是否(fou)有二(er)硫鍵(jian),或需要編輯(ji)的殘基(ji)。
pdb4amber pure_protein.pdb -o pure_check.pdb
cat pure_check_sslink
沒有二硫鍵,之后使用pure_check.pdb。
接下來(lai)在tleap中加載配體(ti)和受體(ti)。
tleap -f - <<_EOF
source leaprc.protein.ff14SB
source leaprc.gaff2
source leaprc.water.tip3p
loadAmberParams frcmod.ionsjc_tip3p
loadoff BEN.lib
loadoff FEN.lib
loadamberparams BEN.gaff2.frcmod
loadamberparams FEN.gaff2.frcmod
ligands = loadpdb lig.pdb
complex = loadpdb pure_check.pdb
complex = combine {ligands complex}
check complex
solvatebox ligands TIP3PBOX 15
addions ligands Na+ 0
savepdb ligands ligands_vdw_bonded.pdb
saveamberparm ligands ligands_vdw_bonded.parm7 ligands_vdw_bonded.rst7
solvatebox complex TIP3PBOX 15
addions complex Cl- 0
savepdb complex complex_vdw_bonded.pdb
saveamberparm complex complex_vdw_bonded.parm7 complex_vdw_bonded.rst7
quit
_EOF
注意根據實際電(dian)(dian)荷情況調整addions,如果ligand/complex帶負電(dian)(dian),加(jia)Na+,反之加(jia)Cl-,離子(zi)類型可以自己選擇。
complex_vdw_bonded.pdb并檢查其中的配(pei)(pei)體(ti)分子區別,是否只是不(bu)一樣(yang)的配(pei)(pei)體(ti)部分不(bu)一樣(yang),確定(ding)配(pei)(pei)體(ti)不(bu)一樣(yang)部分的原子。設(she)置(zhi)出發和結束原子。

紅框是有區(qu)別的(de)原子(zi),其他原子(zi)需要手動編(bian)輯(ji)使其保持一致的(de)坐標(biao),紅線(xian)是編(bian)輯(ji)后的(de)原子(zi)。
修改initial中的in文(wen)件下述部分。
timask1=':BEN',timask2=':FEN',
scmask1=':BEN@H6',scmask2=':FEN@O1,H6'
使6個in文件符合實(shi)際更改的原子(zi)(zi)情況,只改變不同的原子(zi)(zi),該步(bu)驟的目標是優化vdw變化過程中氫原子(zi)(zi)的位置。
使用sbatch run_v100.slurm,將任(ren)務提交到北鯤云(yun)的單卡V100集(ji)群上,該(gai)任(ren)務大(da)概耗(hao)時(shi)1分鐘,注(zhu)意該(gai)任(ren)務如(ru)果(guo)有(you)報錯,可能是以下問題:
如果上述過程報關于ambmask的(de)錯(比(bi)如mask中不能用_符號),可(ke)以(yi)改寫ambmask來獲得正確可(ke)識別的(de)mask。
ambmask-p complex_vdw_bonded.parm7 -c complex_vdw_bonded.rst7 -find :FEN
是ambmask的(de)輸入方(fang)式,在parm7和rst7中尋找這些原子,并(bing)用pdb的(de)形式輸出(chu)。
如(ru)果報(bao)錯Error : Atom 12 does not have match in V1 ,說明(ming)這個原子(zi)(zi)在(zai)兩個小分子(zi)(zi)配(pei)(pei)體中間的位置區別太大,TI不(bu)能識別這兩個分子(zi)(zi)作為一樣的背景(jing),因此(ci)將(jiang)這兩個原子(zi)(zi)(初(chu)始和結束配(pei)(pei)體的)加(jia)入坐(zuo)標(biao)中,就可(ke)以解決這個問題。
運(yun)行(xing)結束后,提取優化過(guo)后的分子結構。
p ligands_vdw_bonded.rst7 ligands_vdw_bonded.rst7.leap
cp press_lig.rst7 ligands_vdw_bonded.rst7
cp complex_vdw_bonded.rst7 complex_vdw_bonded.rst7.leap
cp press_com.rst7 complex_vdw_bonded.rst7
cpptraj -p ligands_vdw_bonded.parm7 <<_EOF
trajin ligands_vdw_bonded.rst7
strip ":1,2"
outtraj ligands_solvated.pdb onlyframes 1
unstrip
strip ":2-999999"
outtraj ligands_BEN.pdb onlyframes 1
unstrip
strip ":1,3-999999"
outtraj ligands_FEN.pdb onlyframes 1
_EOF
cpptraj -p complex_vdw_bonded.parm7 <<_EOF
trajin complex_vdw_bonded.rst7
strip ":1,2"
outtraj complex_solvated.pdb onlyframes 1
unstrip
strip ":2-999999"
outtraj complex_BEN.pdb onlyframes 1
unstrip
strip ":1,3-999999"
outtraj complex_FEN.pdb onlyframes 1
_EOF
再次使用tleap生成decharge,vdw和recharge的文件,decharge是配(pei)體(ti)1,recharge是配(pei)體(ti)2,修改閱讀的pdb的名(ming)字即可。
tleap -f- <<_EOF
# load the AMBER force fields
source leaprc.protein.ff19SB
source leaprc.gaff2
source leaprc.water.tip3p
loadAmberParams frcmod.ionsjc_tip3p
loadOff BEN.lib
loadOff FEN.lib
loadamberparams BEN.gaff2.frcmod
loadamberparams FEN.gaff2.frcmod
# coordinates for solvated ligands as created previously by MD
lsolv = loadpdb ligands_solvated.pdb
lbnz = loadpdb ligands_BEN.pdb
lphn = loadpdb ligands_FEN.pdb
# coordinates for complex as created previously by MD
csolv = loadpdb complex_solvated.pdb
cbnz = loadpdb complex_BEN.pdb
cphn = loadpdb complex_FEN.pdb
# decharge transformation
decharge = combine {lbnz lbnz lsolv}
setbox decharge vdw
savepdb decharge ligands_decharge.pdb
saveamberparm decharge ligands_decharge.parm7 ligands_decharge.rst7
decharge = combine {cbnz cbnz csolv}
setbox decharge vdw
savepdb decharge complex_decharge.pdb
saveamberparm decharge complex_decharge.parm7 complex_decharge.rst7
# recharge transformation
recharge = combine {lphn lphn lsolv}
setbox recharge vdw
savepdb recharge ligands_recharge.pdb
saveamberparm recharge ligands_recharge.parm7 ligands_recharge.rst7
recharge = combine {cphn cphn csolv}
setbox recharge vdw
savepdb recharge complex_recharge.pdb
saveamberparm recharge complex_recharge.parm7 complex_recharge.rst7
quit
_EOF
生成好這些(xie)過程的文件(jian)(jian)后,使用setup_run.sh來產生三個步驟的輸入(ru)文件(jian)(jian),在修改setup_run.sh時,注(zhu)意以下(xia)部分。
decharge_crg=":2@H6"
vdw_crg=":1@H6|:2@O1,H6"
recharge_crg=":1@O1,H6"
decharge="ifsc=0,crgmask='$decharge_crg',"
vdw_bonded="ifsc=1,scmask1=':1@H6',scmask2=':2@O1,H6',crgmask='$vdw_crg'"
recharge="ifsc=0,crgmask='$recharge_crg',"
適配修改(gai),注意H6的(de)都改(gai)成初(chu)始的(de)ambmask,O1,H6的(de)都改(gai)成目(mu)標的(de)ambmask,:前面的(de)1或(huo)者(zhe)2不要改(gai)。
如有(you)必要(yao)修改(gai)λ,改(gai)變prod.tmpl和(he)setup.sh中的(de)值(0.00922開始的(de)那(nei)一串)。
該文件將(jiang)直接生成(cheng)所需的(de)(de)slurm文件,并提交(jiao)到對(dui)應的(de)(de)機(ji)器上,默認使用g-v100-1,運(yun)(yun)行pmemd.cuda,大致運(yun)(yun)行時(shi)間(jian)1小時(shi)左右。
有時候pmemd.cuda會運行失敗(bai),此時轉(zhuan)用(yong)來運行,使(shi)用(yong)run_mpi.py,提交(jiao)命令(ling):
python run_mpi.py ligands 500000 mpi
將會(hui)檢(jian)查(cha)所有ligands下(xia)的(de)文件,對于5分(fen)鐘內沒更新,且info中運(yun)行步驟在500000 以下(xia)的(de),會(hui)提交(jiao)到32核(he)CPU機器上運(yun)行后續的(de)模擬,直到結束。
這個運(yun)(yun)行(xing)步(bu)驟(zou)(zou)非常(chang)緩慢(man),使用32核CPU算(suan)完1納秒的步(bu)驟(zou)(zou)可能需要(yao)7-8天,可以嘗(chang)試(shi)先運(yun)(yun)行(xing)一(yi)段(duan)CPU,再運(yun)(yun)行(xing)一(yi)段(duan),改為python run_mpi.py ligands 500000 cuda即可。
運行結束(shu)后,使用alchemical-analysis/alchemical_analysis/alchemical_analysis.py來分析結果,注意(yi)該(gai)文(wen)件(jian)在python2下運行。
pip2 install matplotlib
pip2 install scipy
pip2 install numpy
pip2 install pymbar==3.0.3
運行:
mkdir -p ana_recharge && cd ana_recharge
../../alchemical-analysis/alchemical_analysis/alchemical_analysis.py -a AMBER -d . -p ../recharge/[01]*/ti00[1-9] -q out -o . -t 300 -v -r 5 -u kcal -f 50 -g -w
mkdir -p ../ana_decharge && cd ../ana_decharge
../../alchemical-analysis/alchemical_analysis/alchemical_analysis.py -a AMBER -d . -p ../decharge/[01]*/ti00[1-9] -q out -o . -t 300 -v -r 5 -u kcal -f 50 -g -w
mkdir -p ../ana_vdw && cd ../ana_vdw
../../alchemical-analysis/alchemical_analysis/alchemical_analysis.py -a AMBER -d . -p ../vdw_bonded/[01]*/ti00[1-9] -q out -o . -t 300 -v -r 5 -u kcal -f 50 -g -w
此(ci)時只能輸出三種(zhong)變化的結果(guo),將其TI 一(yi)列加和得到最終的結果(guo)。
如果見到(dao)pymbar的warning,只要注釋掉對應的assert就(jiu)可(ke)以了。
vim /home/cloudam/.local/lib/python2.7/site-packages/pymbar/timeseries.py
去第162行。


該應用場景解(jie)決將蛋白口(kou)袋內的殘(can)基A變(bian)為(wei)殘(can)基B所(suo)產生的相(xiang)對自(zi)由能(neng)變(bian)。

將蛋白口袋內(nei)的(de)LEU轉化(hua)GLN。
首先,使用(yong)pymol將分(fen)子打開,并(bing)選中殘基,使用(yong)wizard-mutagenesis-protein完(wan)成突(tu)變,或者使用(yong)命令行(xing)完(wan)成突(tu)變:
load*.pdb
cmd.wizard("mutagenesis")
cmd.do("refresh_wizard")
cmd.get_wizard().set_mode("GLN")
cmd.get_wizard().do_select("86/")
cmd.get_wizard().apply()
cmd.set_wizard("done")
save*out_name.pdb,enabled
將突變后(hou)的蛋白文件(jian)(jian)保存為pdb文件(jian)(jian)。
將突(tu)變(bian)(bian)前(qian)的蛋白(WT.pdb)和突(tu)變(bian)(bian)后的蛋白(L86Q.pdb)的pdb文件上傳。使用tleap,讀取(qu)野生型結構。
tleap
source leaprc.protein.ff14SB
source leaprc.gaff2
source leaprc.water.tip3p
loadamberparams frcmod.ionsjc_tip3p
zn = loadpdb WT.pdb
check zn
solvateBox zn TIP3PBOX 10
addIons znCl- 0
savepdb zn box_check.pdb
quit
記(ji)錄盒子的范德華半(ban)徑,并將結構中(zhong)的水(shui)提取出(chu)來備用(yong)。

python dry_for_TI.py box_check.pdb wat.pdb WT_receptor.pdb
同樣的方式讀取突變(bian)型結構,使其保持(chi)Amber的原子順序。
tleap
source leaprc.protein.ff14SB
source leaprc.gaff2
source leaprc.water.tip3p
loadamberparams frcmod.ionsjc_tip3p
zn = loadpdb L86Q.pdb
check zn
solvateBox zn TIP3PBOX 10
addIons zn Cl - 0
savepdb zn L86Q_leap.pdb
quit
python dry_for_TI.py L86Q_leap.pdb wat1.pdbL86Q_dry.pdb
對(dui)比兩個(ge)去水后的文件,發(fa)現由于Amber重編號,突變(bian)的殘基(ji)變(bian)為84位(wei),此時使用check_diff_online.py,來保證不是突變(bian)的殘基(ji)的位(wei)置都絕對(dui)一致,以允許進(jin)行TI過(guo)程。
python check_diff_online.py L84Q L86Q_dry.pdb WT_receptor.pdb 84 L86Q_check.pdbWT_
print的是空字典(dian),說明所有的原子都是匹配的。
再次用tleap讀取(qu)受體和(he)(he)配(pei)體(之前第一部分保存的(de)mol2文件),并讀取(qu)水盒(he)子,其中(zhong)ligand是剛才保存的(de)配(pei)體(不(bu)變(bian)(bian)部分),m1和(he)(he)m2分別是突(tu)變(bian)(bian)前后的(de)部分(注(zhu)意突(tu)變(bian)(bian)只改變(bian)(bian)側鏈,不(bu)改變(bian)(bian)主鏈),注(zhu)意這里使(shi)用了剛才記錄的(de)范德華半(ban)徑。
tleap
source leaprc.protein.ff14SB
source leaprc.gaff2
loadOff FEN.lib
loadamberparams FEN.gaff2.frcmod
source leaprc.water.tip3p
loadamberparams frcmod.ionsjc_tip3p
ligand = loadmol2 FEN.gaff2.mol2
m1 = loadpdb L86Q_check.pdb
m2 = loadpdb WT_check.pdb
w = loadpdbwat1.pdb
protein = combine {m1 m2 w}
complex = combine {m1 m2 ligand w}
set default nocenter on
setBox protein vdw{39.415 43.577 52.292}
savepdb protein protein.pdb
saveamber parmprotein protein.parm7 protein.rst7
setBox complex vdw {39.415 43.577 52.292}
savepdb complex complex.pdb
saveamberparm complex complex.parm7 complex.rst7
quit
使用(yong)parmed處理protein.parm7和complex.parm7,以保證正確(que)的(de)所(suo)改變(bian)的(de)位置(zhi)提供給TI運算,此(ci)處的(de)162為WT或者(zhe)突(tu)變(bian)體(ti)的(de)殘基數,@后面的(de)內容是python get_mutation.py L86Q得(de)到的(de)mapping結果,是在那個突(tu)變(bian)殘基上但也沒有變(bian)化的(de)殘基,84是突(tu)變(bian)位置(zhi),246是84+162。
parmed protein.parm7 <<_EOF
loadRestrt protein.rst7
setOverwrite True
tiMerge :1-162 :163-324 :84&!@CA,C,O,N,H,HA,CB :246&!@CA,C,O,N,H,HA,CB
outparm merged_L86Q_protein.parm7 merged_L86Q_protein.rst7
quit
_EOF
parmed complex.parm7 <<_EOF
loadRestrt complex.rst7
setOverwrite True
tiMerge :1-162 :163-324 :84&!@CA,C,O,N,H,HA,CB :246&!@CA,C,O,N,H,HA,CB
outparm merged_L86Q_complex.parm7 merged_L86Q_complex.rst7
quit
_EOF
正(zheng)確獲得這些文(wen)件(jian)后,使(shi)用:
python auto_gene_inp_run.py 84 162 CA,C,O,N,H,HA,CB L86Q
來生成tmpl文件(jian)(jian)(run.tmpl文件(jian)(jian)需要上傳),并(bing)將tmpl文件(jian)(jian)轉成真正的(de)ti文件(jian)(jian)放進文件(jian)(jian)夾(jia),同時使用slurm提交最小化,加熱和運行(xing)步驟。
注意(yi),這(zhe)里直接(jie)使用(yong)cuda很容易(yi)斷,可以適當自己修改(gai)之(zhi)前的run_mpi.py來(lai)使用(yong)cpu續跑中(zhong)斷的模擬。
運(yun)行(xing)結束(shu)后的分析略。
本案(an)例(li)將實際運行一(yi)個蛋(dan)(dan)白(bai)(bai)-蛋(dan)(dan)白(bai)(bai)相互作用上的(de)突變(bian)。我們計算新冠(guan)病毒受(shou)體(ti)結合區域(rbd)到人(ren)ACE2受(shou)體(ti)(ace2)復合物(wu)上發生(sheng)Omicron的(de)突變(bian)之一(yi)的(de)返回突變(bian)A484E后的(de)結合自由能(neng)變(bian)化。
生成連接了二硫(liu)鍵(jian)的大復合體(ti)水(shui)(shui)盒,記錄vdw盒子大小,額外(wai)加(jia)入(ru)0.15M/L NaCL (總水(shui)(shui)數量*0.002772)。
tleap
source leaprc.protein.ff14SB
source leaprc.gaff
source leaprc.water.tip3p
loadamberparams frcmod.ionsjc_tip3p
zn = loadpdb omi_SS.pdb
bond zn.333.SG zn.358.SG
bond zn.376.SG zn.429.SG
bond zn.388.SG zn.522.SG
bond zn.477.SG zn.485.SG
bond zn.637.SG zn.645.SG
bond zn.848.SG zn.865.SG
bond zn.1034.SG zn.1046.SG
checkzn
solvateBox zn TIP3PBOX 10
addIons zn Na+ 0
addIons znNa+8 0
addIons znCl- 0
savepdb zn box_check.pdb
quit
python dry_for_TI.py box_check.pdb wat1.pdb omi_rbd.pdb,omi_ace2.pdb
同(tong)時(shi)生成一個小的(de)水(shui)盒,用于跑蛋白部(bu)分的(de)TI。
tleap
source leaprc.protein.ff14SB
source leaprc.gaff
source leaprc.water.tip3p
loadamberparams frcmod.ionsjc_tip3p
m1 = loadpdb omi_rbd.pdb
bond m1.4.SG m1.29.SG
bond m1.47.SG m1.100.SG
bond m1.59.SG m1.193.SG
bond m1.148.SG m1.156.SG
check m1
solvateBox m1 TIP3PBOX 10
addIons m1 Na+ 0
addIons m1Na+ 28
addIons m1Cl- 0
savepdb m1 ligands_recharge.pdb
quit
python dry_for_TI.py ligands_recharge.pdb rbd_wat.pdb test_rbd.pdb
檢查(cha)輸(shu)入的(de)是(shi)否在TI區域之外沒(mei)有區別,注意輸(shu)入的(de)順序,前面是(shi)突變后的(de)蛋白,后面是(shi)原始(shi)的(de)蛋白。
python check_diff_online.py A152EA484E_rbd.pdb omi_rbd.pdb 152 A484E_check.pdb omi_check.pdb
設定正確的二硫鍵(jian),分別加(jia)載兩種水(shui)盒,輸出(chu)protein和complex的拓撲學文件和坐(zuo)標文件。
tleap
source leaprc.protein.ff14SB
source leaprc.gaff
source leaprc.water.tip3p
loadamberparams frcmod.ionsjc_tip3p
ligand = loadpdbomi_ace2.pdb
bond ligand.308.SG ligand.316.SG
bond ligand.519.SG ligand.536.SG
bond ligand.705.SG ligand.717.SG
m1 = loadpdb omi_check.pdb
bond m1.4.SG m1.29.SG
bond m1.47.SG m1.100.SG
bond m1.59.SG m1.193.SG
bond m1.148.SG m1.156.SG
m2 = loadpdb A484E_check.pdb
bond m2.4.SG m2.29.SG
bond m2.47.SG m2.100.SG
bond m2.59.SG m2.193.SG
bond m2.148.SG m2.156.SG
w1 = loadpdb rbd_wat.pdb
w2 = loadpdb wat1.pdb
protein = combine {m1m2w1}
complex = combine {m1m2ligandw2}
set default nocenter on
setBox protein vdw{43.21553.42159.922}
savepdb protein protein.pdb
saveamber parmprotein protein.parm7 protein.rst7
setBox complex vdw{64.171 75.490 114.587}
savepdb complex complex.pdb
saveamberparm complex complex.parm7 complex.rst7
quit
使用(yong)parmed處理(li)protein.parm7和complex.parm7,以保證正確的(de)位置(zhi)。
parmedprotein.parm7 <<_EOF
loadRestrt protein.rst7
setOverwrite True
tiMerge :1-193 :194-386 :152&!@CA,C,O,N,H,HA,CB :345&!@CA,C,O,N,H,HA,CB
outparm merged_A484E_protein.parm7 merged_A484E_protein.rst7
quit
_EOF
parmed complex.parm7 <<_EOF
loadRestrt complex.rst7
setOverwrite True
tiMerge :1-193 :194-386 :152&!@CA,C,O,N,H,HA,CB :345&!@CA,C,O,N,H,HA,CB
outparm merged_A484E_complex.parm7 merged_A484E_complex.rst7
quit
_EOF
正確獲(huo)得這些文件后,使(shi)用(yong):
python auto_gene_inp_run.py 152 193 CA,C,O,N,H,HA,CB A484E
來生成(cheng)tmpl文件(jian)(run.tmpl文件(jian)需(xu)要上傳(chuan)),并將tmpl文件(jian)轉成(cheng)真正的ti文件(jian)放進文件(jian)夾,同時(shi)使用(yong)slurm提交最小化,加熱和運行步驟(5ns)。