Difference between revisions of "2022 AMBER tutorial 3 with PDBID 1X70"

Revision as of 20:08, 3 May 2022

"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