vsREMD-复合物体系案例

1.   lig（0415）pre

1.1  获取三维结构（ligprep/chem3D)

# 这个三维结构需要有氢（本地ligprep生成）

1.2  使用antechamber生成配体小分子参数（服务器：cpu256）

antechamber -i 0415.pdb -fi pdb -o 0415.gjf -fo gcrt
antechamber -i 0415.sdf -fi sdf -o 0415.gjf -fo gcrt


#修改0415.gjf文件
#HF改成B3LYP; 增加scrf=(solvent=water);  %mem=4GB； %nproc=16  SCF=NoVarAcc


nohup g16 0415.gjf &
[1] 250920


vi 0415.log

antechamber -i 0415.log -fi gout -o 0415.prep -fo prepi -c resp
# 这一步生成的NEWPDB.PDB 很重要
parmchk2 -i 0415.prep -f prepi -o 0415.frcmod -a y

1.3  通过对接调整小分子的位置信息

#蛋白质去氢
conda activate pdb2pqr预测
# 采用14个蛋白构象（13个apo聚类+晶体结构）


# 进行盲对接（drawgridbox)/ 根据口袋进行对接
# drawgridbox
import drawgridbox 
drawgridbox repc0, nx=5, ny=5, nz=5,  lw=1, g=0, b=0
Box dimensions (71.01, 60.11, 52.25)
Center of Mass: [ 81.581,  29.646,  12.226]


# 用高斯得到的 NEWPDB.PDB或转成sdf文件 去对接（获取小分子配体位置信息）
obabel NEWPDB.PDB -osdf -O 0415_g16.sdf


# 每个盲对接结果里选两个小分子位置 2*14
for i in {0..13};do smina --seed 1 --num_modes 10 -r ../protein_aligned/c${i}_H.pdb -l NEWPDB.PDB --center_x  81.581 --center_y 29.646 --center_z 12.226 --size_x 75 --size_y 65 --size_z 60 -o blind${i}.pdb;done
for i in {0..13};do smina --seed 1 --num_modes 10 -r ../protein_aligned/c${i}_H.pdb -l 0415_g16.sdf --center_x  81.581 --center_y 29.646 --center_z 12.226 --size_x 75 --size_y 65 --size_z 60 -o blind${i}.pdb;done

• 这里的 c13_H.pdb 就是 6lk0fix_H.pdb 

#将dock.pdb的atom name修改成gaussian参数生成的类型，注意通过pymol可视化检查
#修改 dock.pdb 
for i in {0..27};do cat dock${i}.pdb | cut -b 29-100 > location${i};done
grep -v "H" NEWPDB.PDB > num_test
# 检查原子序列是否一致，若不一致，手动进行调整
cat num_test| cut -b 1-28 > num


for i in {0..27};do paste -d "" num location${i} > 6274_${i}.pdb;done

2.  28个副本初始结构准备

for i in {0..13};do mkdir p${i};done
for i in {0..13};do cp parameters/pre.sh ./p${i}/pre${i}.sh;done
for i in {0..13};do cp protein/c${i}.pdb ./p${i};done
for i in {0..13};do sed -i "s/repcM/c${i}/g" ./p${i}/pre${i}.sh;done




vi tleap16.in
vi tleap20.in
vi tleap21.in
solvatebox com TIP3PBOX 11


vi tleap0.in
vi tleap1.in
vi tleap18.in
vi tleap19.in
solvatebox com TIP3PBOX 12


vi tleap26.in
vi tleap27.in
solvatebox com TIP3PBOX 8


tleap -s -f tleap${i}.in
tleap -s -f tleap1.in
tleap -s -f tleap2.in


sh change.sh
vi gromacs.gro 
#以p8/17 为基准
原子数54365
溶剂数16324
#topol文件；posre1.itp;posre2.itp;index.ndx 全部采用 p8/17


gmx_mpi make_ndx -f gromacs.gro -o index.ndx 
1 /| 17|


echo 1|gmx_mpi genrestr -f gromacs.gro -o posre1.itp
echo 17|gmx_mpi genrestr -f gromacs.gro -o posre2.itp #限制小分子


source ../../../lywu.sh
nohup sh adjust.sh &


#gro
for i in {0..27};do mv npt$i.gro pre$i.gro;done


# mdtest
cp ../vst/md.mdp ./
cp ../vst/npt.mdp ./
cp ../vst/temp ./
for i in {0..27};do cp md.mdp md$i.mdp;done
for i in {0..27};do cp npt.mdp npt$i.mdp;done
for a in {0..27};do b=`expr ${a} + 1`; c=`sed -n "${b}p" temp`; sed "s/300.00/${c}/g" md.mdp > md${a}.mdp; done
for a in {0..27};do b=`expr ${a} + 1`; c=`sed -n "${b}p" temp`; sed "s/300/${c}/g" npt.mdp > npt${a}.mdp; done

• 检查mdp文件

which gmx_mpi
source ~/lywu/lywu.sh
which gmx_mpi


for i in {0..27};do gmx_mpi grompp -f npt${i}.mdp -c ../gro/pre${i}.gro -p ../topol/topol.top -n ../topol/index.ndx -o npt${i}.tpr -maxwarn 1; done


nohup sh -c 'for i in {0..27};do mpirun -np 50 gmx_mpi mdrun -v -deffnm npt${i}; done' &

• npt 18后面的副本均崩溃，处理方法另外记录

# 0.1ns test
for i in {0..27};do sed -i 's/50000000/50000/g' md$i.mdp md$i.mdp;done


mkdir 608_0.1ns_test
for i in {0..27};do gmx_mpi grompp -f ../md${i}.mdp -c ../npt${i}.gro -p ../../topol/topol.top -n ../../topol/index.ndx -o md${i}.tpr; done


nohup mpirun -np 84 gmx_mpi mdrun -v -deffnm md -multi 28 -replex 1000 &
[1] 109663




# 2ns test
for i in {0..27};do sed -i 's/50000/1000000/g' md$i.mdp md$i.mdp;done


mkdir 608_2ns_test
for i in {0..27};do gmx_mpi grompp -f ../md${i}.mdp -c ../npt${i}.gro -p ../../topol/topol.top -n ../../topol/index.ndx -o md${i}.tpr; done


nohup mpirun -np 84 gmx_mpi mdrun -v -deffnm md -multi 28 -replex 1000 &
#查看任务
1.tail -f nohup.out  #这个-f 能让我们看到正在进行的任务
2.jobs
3.vi nohup.out
4.vi md0.log

HPC_Pre

#HPC
# 500ns
for i in {0..27};do sed -i 's/2250000000000/250000000/g' md$i.mdp md$i.mdp;done
for i in {0..27};do sed -i 's/1000000/250000000/g' md$i.mdp md$i.mdp;done
for i in {0..27};do gmx_mpi grompp -f ../md${i}.mdp -c ../npt${i}.gro -p ../../topol/topol.top -n ../../topol/index.ndx -o md${i}.tpr; done
# test
for i in {0..27};do sed -i 's/250000000/50000/g' md$i.mdp md$i.mdp;done
for i in {0..27};do gmx_mpi grompp -f ../md${i}.mdp -c ../npt${i}.gro -p ../../topol/topol.top -n ../../topol/index.ndx -o md${i}.tpr; done
nohup mpirun -np 140 --oversubscribe gmx_mpi mdrun -v -deffnm md -multi 28 -replex 1000 &
nohup mpirun -np 84 gmx_mpi mdrun -v -deffnm md -multi 28 -replex 1000 &
test  [1] 236409
500ns 

• 先不测试，cpu交满
• 以及，一定要 cat md0.log 检查！