复合物MD

配体准备

# 22224 
# 172.21.85.24  
# 使用薛定谔优化结构（图形化界面操作）
# pH = 8.0，手性从3D结构定（晶体结构），力场为2017版的默认力场OPLS3
# 高斯优化
# 原始文件 /home/databank/lywu/trypsin/inhibitiors/gaussian_job 


# 注意检查体系的电荷和自选多重度
antechamber -i lig.sdf -fi sdf -o lig.gjf -fo gcrt -at gaff2 -gn "%nproc=4" -gm "%mem=4GB" -gk "#B3LYP/6-31G* em=gd3bj pop=MK iop(6/33=2,6/42=6) opt" -rn MOL  -nc 0  


###################################################################################
# 85.24 /home/databank/zzy/project/MD/koff/trypsin/ligprep-opt
# 基于高斯优化日志文件，拟合RESP电荷，生成AMBER需要的参数文件（prep)
conda activate ambertools
for i in `ls`;do cd $i;antechamber -i lig.log -fi gout -o lig.prep -fo prepi -c resp;parmchk2 -i lig.prep -f prepi -o lig.frcmod -a y;cd ../;done


# 这一步本质上是获得坐标，尝试直接从晶体结构中抠出小分子，提取坐标
# gv里分别打开 NEWPDB.PDB 与 5eg4_ligand_native_277.mol2（晶体结构配体文件）
# 打开GV tools-atom list-鼠标点击symbol列-rows-sort selected -拖动鼠标选择所有H-edit-delete-selected atoms），接着仍在atom list中，edit-z matrix-standardize，分别保存。
# 得到 g16-gv.pdb crystal-gv.pdb
# GV里 检查 g16-gv.pdb crystal-gv.pdb，更改crystal-gv.pdb （atom list - Tag)，使其Tag 与g16-gv.pdb 完全一致。


dos2unix *
grep -v "CONECT" g16-gv.pdb |cut -b 1-26 > symbol
cat crystal-gv.pdb | cut -b 29-78|sed '1,2d'|sed '1i Combineing the coordinates from the crystal structure with the atomic numbers from the Gaussian optimized structure'> coordinate
paste -d " " symbol coordinate|sed '1i TITLE      lig.pdb'> lig.pdb
### 可视化检查

受体准备

# PDBBIND
# /home/databank/zzy/project/MD/koff/trypsin/pro-pre/pdbbind
conda activate pdb2pqr
pdb2pqr30 5eg4_ligand_native_277_protein.pdb 5eg4_pred.pqr --ff AMBER --ffout AMBER --with-ph 8.0 --pdb-output 5eg4_pred.pdb


# 检查蛋白是否有二硫键 （pymol-show-disulfide, select resn CYS 检查位置)
# 如果蛋白中存在二硫键，需要进行特殊的处理（教程待完善）

复合物准备

# /home/databank/zzy/project/MD/koff/trypsin/complex-pre/trypsin-31
# tleap


# 蛋白质力场
source leaprc.protein.ff19SB
# 配体力场
source leaprc.gaff2
# 溶剂力场
source leaprc.water.opc
# 导入配体参数
loadamberparams /home/databank/zzy/project/MD/koff/trypsin/ligprep-opt/g16/31/lig.frcmod
loadamberprep /home/databank/zzy/project/MD/koff/trypsin/ligprep-opt/g16/31/lig.prep
# 加载配体文件
mol=loadpdb lig.pdb 
# 加载蛋白文件
pro=loadpdb pro-tleap.pdb
# 组装蛋白配体
com=combine {pro mol}


# 保存amber参数（不加溶剂的）
saveamberparm pro pro.prmtop pro.inpcrd
saveamberparm mol lig.prmtop lig.inpcrd
saveamberparm com native.prmtop native.inpcrd
savepdb com com-dry.pdb


# 检查电荷
charge com
# 加溶剂
solvatebox com OPCBOX 10
# 平衡电荷 （体系原本带正电补氯离子，带负电补钠离子）
addionsrand com Cl- 6
charge com
saveamberparm com complex.prmtop complex.inpcrd
savepdb com com.pdb

# AMBER 参数转 GMX 参数
#进入python 
import parmed as pmd
amber = pmd.load_file('complex.prmtop','complex.inpcrd')
amber.save('gmx.top')
amber.save('gmx.gro')

# 在蛋白质结束和小分子开始前添加位置限制（moleculetype）
sed -i '36213 a ; Include Position restraint file\n#ifdef  POSRES\n#include "posre1.itp\n#endif\n' gmx.top
# 在小分子结束和离子开始前添加位置限制（moleculetype）
sed -i '37443 a ; Include Position restraint file\n#ifdef  POSRES\n#include "posre2.itp\n#endif\n' gmx.top


# 生成index (GMX 里没有 protein_MOL 这种现成的 group，我们可以通过index文件生成）
# 例如我们想生成Protein_MOL 和 Water_Cl- 这两个组，可如下操作
gmx_mpi make_ndx -f gmx.gro -o index.ndx
> 1|15
> 17|16
# 注意，group序号以你自己的为准，上面的只是示例


# 生成位置限制文件，规定哪些原子被限制
# 这个用于限制蛋白
gmx_mpi genrestr -f gmx.gro -o posre1.itp
# 这个用于限制配体 （下面的62 替换为你的配体原子数+4）
sed -n '1,62p' posre1.itp > posre2.itp

体系平衡

mkdir em nvt npt mdp


# 85.24
cp /home/databank/lywu/vsremd_plus/protein-ligand/trypsin-benzamidine/pre-equilibrium/mdp-file/* ./mdp/


# 把 gmx.gro gmx.top 复制过来


### 能量最小化
gmx_mpi grompp -f ../mdp/em.mdp -c ../gmx.gro -r ../gmx.gro -p ../gmx.top -o em.tpr 
nohup mpirun -np 100 gmx_mpi mdrun -v -deffnm em &


### 等温等容预平衡
# 修改为1fs
gmx_mpi grompp -f ../mdp/nvt.mdp -c ../em/em.gro -r ../em/em.gro -p ../gmx.top -o nvt.tpr
nohup mpirun -np 100 gmx_mpi mdrun -v -deffnm nvt &
# 修改为2fs
gmx_mpi grompp -f ../mdp/nvt.mdp -c nvt.gro -r nvt.gro -p ../gmx.top -o nvt2.tpr
nohup mpirun -np 100 gmx_mpi mdrun -v -deffnm nvt2 &


### 等温等压预平衡
# 修改为1fs
gmx_mpi grompp -f ../mdp/npt.mdp -c ../nvt/nvt2.gro -r ../nvt/nvt2.gro -p ../gmx.top -o npt.tpr
nohup mpirun -np 100 gmx_mpi mdrun -v -deffnm npt &
# 修改为2fs
gmx_mpi grompp -f ../mdp/npt.mdp -c npt.gro -r npt.gro -p ../gmx.top -o npt2.tpr
nohup mpirun -np 100 gmx_mpi mdrun -v -deffnm npt2 &


# 检查npt2.gro的构象