准备好对接后的蛋白pro.pdb和配体环肽文件test.mo

第一步：

python 00_make_single_lig.py test.mol2 test_LIG_raw.mol2

00_make_single_lig.py文件如下：

#!/usr/bin/env python3 import sys if len(sys.argv) != 3: print("Usage: python 00_make_single_lig.py input.mol2 output.mol2") sys.exit(1) infile = sys.argv[1] outfile = sys.argv[2] section = None with open(infile, "r") as fin, open(outfile, "w") as fout: for line in fin: if line.startswith("@<TRIPOS>"): section = line.strip() fout.write(line) if section == "@<TRIPOS>SUBSTRUCTURE": fout.write(" 1 LIG 1 GROUP 0 **** **** 0 ROOT\n") continue if section == "@<TRIPOS>ATOM": parts = line.split() if len(parts) >= 9: atom_id = int(parts[0]) atom_name = parts[1] x = float(parts[2]) y = float(parts[3]) z = float(parts[4]) atom_type = parts[5] charge = float(parts[8]) fout.write( f"{atom_id:7d} {atom_name:<8s} " f"{x:10.4f} {y:10.4f} {z:10.4f} " f"{atom_type:<8s} {1:5d} {'LIG':<8s} {charge:10.6f}\n" ) else: fout.write(line) elif section == "@<TRIPOS>SUBSTRUCTURE": continue else: fout.write(line)

第二步：antechamber，由于环肽原子量很多，这一步要持续几个小时，移植到能量收敛到0.0005.

antechamber -i test_LIG_raw.mol2 -fi mol2 \ -o LIG.mol2 -fo mol2 \ -at gaff2 -c bcc -nc 0 \ -rn LIG -s 2

第三步：

parmchk2 -i LIG.mol2 -f mol2 -o LIG.frcmod -s gaff2

第四步：处理蛋白

pdb4amber -i pro.pdb -o pro_amber.pdb --reduce --dry

第五步：继续处理单边加氢等：

awk ' /^ATOM/ || /^HETATM/ { atom=substr($0,13,4) gsub(/ /,"",atom) if (atom ~ /^H/) next } {print} ' pro_amber.pdb > pro_amber_noH.pdb

第六步：
将配体受体合并：

tleap -f 01_tleap_complex.in | tee 01_tleap_complex.log

01_tleap_complex.in文件如下：

source leaprc.protein.ff19SB source leaprc.gaff2 source leaprc.water.opc loadamberparams LIG.frcmod pro = loadpdb pro_amber_noH.pdb lig = loadmol2 LIG.mol2 complex = combine {pro lig} check complex charge complex solvateoct complex OPCBOX 10.0 addionsrand complex Na+ 0 addionsrand complex Cl- 0 saveamberparm complex complex.prmtop complex.inpcrd savepdb complex complex_solvated.pdb quit

第八步：就是正常MD流程：
1. 限制性最小化：02_min1.in

Minimization 1: solvent and ions, solute restrained &cntrl imin=1, maxcyc=10000, ncyc=5000, cut=10.0, ntb=1, ntr=1, restraint_wt=10.0, restraintmask='!:WAT,Na+,Cl-', / END

2. 全体系最小化 03_min2.in

Minimization 2: whole system &cntrl imin=1, maxcyc=10000, ncyc=5000, cut=10.0, ntb=1, / END

3. 升温04_heat.in

Heating from 0 to 300 K &cntrl imin=0, irest=0, ntx=1, ntb=1, cut=10.0, nstlim=250000, dt=0.002, tempi=0.0, temp0=300.0, ntc=2, ntf=2, ntt=3, gamma_ln=2.0, ntpr=1000, ntwx=5000, ntwr=5000, ntr=1, restraint_wt=5.0, restraintmask='!:WAT,Na+,Cl-', nmropt=1, / &wt TYPE='TEMP0', istep1=0, istep2=250000, value1=0.0, value2=300.0 / &wt TYPE='END' / END

5. NPT 平衡，有弱限制 05_eq1.in

NPT equilibration with weak restraint &cntrl imin=0, irest=1, ntx=5, ntb=2, ntp=1, taup=2.0, cut=10.0, nstlim=500000, dt=0.002, temp0=300.0, ntc=2, ntf=2, ntt=3, gamma_ln=2.0, ntpr=1000, ntwx=5000, ntwr=5000, ntr=1, restraint_wt=1.0, restraintmask='!:WAT,Na+,Cl-', / END

6. NPT 平衡，无限制 06_eq2.in

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

8. 生产模拟，单段 10 ns 07_prod.in

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

9.执行运行命令run_md.sh

#!/bin/bash set -e PMEMD=pmemd.cuda echo "Step 1: restrained minimization" $PMEMD -O \ -i 02_min1.in \ -o 02_min1.out \ -p complex.prmtop \ -c complex.inpcrd \ -r 02_min1.rst7 \ -ref complex.inpcrd echo "Step 2: full minimization" $PMEMD -O \ -i 03_min2.in \ -o 03_min2.out \ -p complex.prmtop \ -c 02_min1.rst7 \ -r 03_min2.rst7 echo "Step 3: heating" $PMEMD -O \ -i 04_heat.in \ -o 04_heat.out \ -p complex.prmtop \ -c 03_min2.rst7 \ -r 04_heat.rst7 \ -x 04_heat.nc \ -ref 03_min2.rst7 echo "Step 4: NPT equilibration with restraint" $PMEMD -O \ -i 05_eq1.in \ -o 05_eq1.out \ -p complex.prmtop \ -c 04_heat.rst7 \ -r 05_eq1.rst7 \ -x 05_eq1.nc \ -ref 04_heat.rst7 echo "Step 5: NPT equilibration without restraint" $PMEMD -O \ -i 06_eq2.in \ -o 06_eq2.out \ -p complex.prmtop \ -c 05_eq1.rst7 \ -r 06_eq2.rst7 \ -x 06_eq2.nc echo "Step 6: production MD" prev=06_eq2.rst7 for i in $(seq -w 1 10) do echo "Production segment ${i}" $PMEMD -O \ -i 07_prod.in \ -o prod_${i}.out \ -p complex.prmtop \ -c ${prev} \ -r prod_${i}.rst7 \ -x prod_${i}.nc \ -inf prod_${i}.info prev=prod_${i}.rst7 done echo "MD finished."

最后运行：

chmod +x run_md.sh ./run_md.sh