2022 AMBER tutorial 3 with PDBID 1X70
"introduction"
In this tutorial, we will be doing molecular dynamics (MD) simulations using Amber 16. This will allow us to obtain much more detailed information about the binding affinities of a ligand with a receptor. Like the DOCK tutorial, we should start with setting up the directory. This can be done with the following command.
mkdir 001_structure 002_parameters 003_leap 004_equil 005_production
"Making the structures"
It is recommended to start fresh with the receptor and ligand files. Do the same process as in the DOCK tutorial. Start with importing pdb code 1X70, then isolate one of the proteins and its associated ligand. With an isolated protein, you can do the following command
Select -> residue -> all nonstandard
then do
actions -> atoms/bonds -> delete
This should isolate the protein and remove all of the nonstandard amino acids. With this, we can save this as a pdb and title it 1X70_fresh.pdb
Moving on to the ligand, we can repeat the same process as above but just isolate the ligand instead. So start with the fresh pdb code, then select the ligand and do
Select -> Invert (selected molecules)
and
actions -> atoms/bonds -> delete
This should isolate the ligand. Then we have to make sure it is charged properly and has hydrogens. We can do that with the following two commands.
Tools -> Structure Editing -> Add H Tools -> Structure Editing -> Add Charge -> (have Amber ff14SB and AM1-BCC selected) -> Ok
Once this is done, save the file as 1x70_lig_wH.mol2
After doing this with both files, upload them to seawulf and place them in the 001_structure directory.
"Simulation parameters"
In this section, we are generating the force field parameters that will be used in future steps. The first step is to swap to the parameters directory, 002_parameters, then use antechamber.
antechamber -i ../001_structure/1x70_lig_wH.mol2 -fi mol2 -o 1x70_ligand_antechamber.mol2 -fo mol2 -at gaff2 -c bcc -rn LIG -nc 1
The -nc tag at the end is the total charge of the ligand, so this may be different depending on the system. In this case, the net charge of the ligand is 1. After this step is completed, we need to parmch2 to generate another file.
parmchk2 -i 1x70_ligand_antechamber.mol2 -f mol2 -o 1x70_ligand.am1bcc.frcmod
"TLeap"
The next step is to create the files needed to simulate the system. Copy the following code into a file named leap.in
Last login: Mon Apr 25 06:08:28 on console (base) noahdean@Noahs-MacBook-Pro ~ % ssh nkdean@login.seawulf.stonybrook.edu nkdean@login.seawulf.stonybrook.edu's password: Reading $DUO_PASSCODE...
Pushed a login request to your device... Success. Logging you in... [nkdean@login2 ~]$ AMS_536 [nkdean@login2 noah]$ cd 1x70_AMBER/ [nkdean@login2 1x70_AMBER]$ cd 003_leap/ [nkdean@login2 003_leap]$ ls 1x70.complex.parm7 1x70.gas.ligand.rst7 1x70.wet.complex.pdb 1x70.gas.complex.pdb 1x70.gas.receptor.parm7 1x70.wet.complex.rst7 1x70.gas.complex.rst7 1x70.gas.receptor.rst7 leap.in 1x70.gas.ligand.parm7 1x70.wet.complex.parm7 leap.log [nkdean@login2 003_leap]$ vim leap.in [nkdean@login2 003_leap]$ cd ../001_structure/ [nkdean@login2 001_structure]$ ls 1x70_fresh.pdb 1x70_ligand_wH.mol2 1x70_receptor_fresh.pdb [nkdean@login2 001_structure]$ vim 1x70_fresh.pdb [nkdean@login2 001_structure]$ cd ../003_leap/ [nkdean@login2 003_leap]$ vim leap.in
#!/usr/bin/sh
###load protein force field
source leaprc.protein.ff14SB
###load GAFF force field (for our ligand)
source leaprc.gaff
###load TIP3P (water) force field
source leaprc.water.tip3p
###load ions frcmod for the tip3p model
loadamberparams frcmod.ionsjc_tip3p
###needed so we can use igb=8 model
set default PBradii mbondi3
###load protein pdb file
rec=loadpdb ../001_structure/1x70_fresh.pdb
#bond rec.649.SG rec.762.SG
#bond rec.454.SG rec.472.SG
#bond rec.444.SG rec.447.SG
#bond rec.328.SG rec.339.SG
#bond rec.385.SG rec.394.SG
###load ligand frcmod/mol2
loadamberparams ../002_parameters/1x70_ligand.am1bcc.frcmod
lig=loadmol2 ../002_parameters/1x70_ligand_antechamber.mol2
###create gase-phase complex
gascomplex= combine {rec lig}
###write gas-phase pdb
savepdb gascomplex 1x70.gas.complex.pdb
###write gase-phase toplogy and coord files for MMGBSA calc
saveamberparm gascomplex 1x70.complex.parm7 1x70.gas.complex.rst7
saveamberparm rec 1x70.gas.receptor.parm7 1x70.gas.receptor.rst7
saveamberparm lig 1x70.gas.ligand.parm7 1x70.gas.ligand.rst7
###create solvated complex (albeit redundant)
solvcomplex= combine {rec lig}
###solvate the system
solvateoct solvcomplex TIP3PBOX 12.0
###Neutralize system
addions solvcomplex Cl- 0
addions solvcomplex Na+ 0
#write solvated pdb file
savepdb solvcomplex 1x70.wet.complex.pdb
###check the system
charge solvcomplex
check solvcomplex
###write solvated toplogy and coordinate file
saveamberparm solvcomplex 1x70.wet.complex.parm7 1x70.wet.complex.rst7