Prepare a reference structure
The INPHARMA program uses the PDB atom numbering for the calculations of the theoretical INPHARMA NOEs. It is therefore necessary to properly prepare the PDB files prior to the calculations. These steps models missing residues or atoms, and re-maps the atom numbering such that it is consistent between difference PDB files.
Contents
Process the structure file
As a reference, we will download a PDB file (PDB id 3DND) containing ligand (with id LL1).
pdb="3DNE"
ligand="LL1"
# Fetch the PDB file
wget -Oxray.pdb http://www.pdb.org/pdb/files/${pdb}.pdb
To fetch the ligand you can use the following PyMOL script. Note: This script will also align (superimpose) the x-ray structure onto a reference PDB file (all structures your prepare for INPHARMA should be aligned on the binding site).
# Extract the ligand
cat<<EOF>align.pml
# Load the xray PDB file
load xray.pdb
# A reference to align to
load reference.amber.pdb, reference
align xray, reference
save xray_fit.pdb, xray
select ligand, resname LL1
save ligand.pdb, ligand
EOF
# Run the PyMOL script
pymol align.pml
# Extract ATOM records only
grep "^ATOM" xray_fit.pdb > xray_ATOMs.pdb
At this point your should have the following files in your directory:
- xray_ATOMs.pdb
- raw PDB file containing the receptor structure ONLY. Should be aligned to a reference PDB file.
- ligand.pdb
- raw PDB file of the ligand ONLY. Should be in the binding pocket of the receptor.
Modeller
The next step is to ensure that our receptor structure is complete (no missing residues/atoms, and the correct sequence). For this we will use two Modeller-script (Python).
The first script aligns the sequence of the structure (xray_ATOMs.pdb) and a sequence file. The sequence file should contain the protein sequence in Modeller format (named target_sequence.ali).
Requirement: Modeller should be installed and working.
>P1;target
sequence:target:::::::0.00: 0.00
HLDQFERIKTLGTGSFGRVMLVKHKETGNHYAMKILDKQKVVKLKQIEHTLNEKRILQAV
NFPFLVKLEFSFKDNSNLYMVMEYVPGGEMFSHLRRIGRFSEPHARFYAAQIVLTFEYLH
SLDLIYRDLKPENLLIDQQGYIQVTDFGFAKRVKGRTWTLCGTPEYLAPEIILSKGYNKA
VDWWALGVLIYEMAAGYPPFFADQPIQIYEKIVSGKVRFPSHFSSDLKDLLRNLLQVDLT
KRFGNLKNGVNDIKNHKWFATTDWIAIYQRKVEAPFIPKFKGPGDTSNFDDYEEEEIR*
from modeller import *
env = environ()
aln = alignment(env)
mdl = model(env, file='xray_ATOMs.pdb', model_segment=('FIRST:A','LAST:A'))
aln.append_model(mdl, align_codes='XRAY', atom_files='xray_ATOMs.pdb')
aln.append(file='target_sequence.ali', align_codes='target')
aln.align2d()
aln.write(file='alignment.ali', alignment_format='PIR')
aln.write(file='alignment.pap', alignment_format='PAP')
Run the sequence alignment with the following command:
python 01_align.py
Next script makes a new receptor PDB file based on the sequence alignment generated above:
from modeller import *
from modeller.automodel import *
log.verbose()
env = environ()
# directories for input atom files
env.io.atom_files_directory = ['.']
a = automodel(env, alnfile='alignment.ali',
knowns='XRAY', sequence='target',
assess_methods=(assess.DOPE, assess.GA341))
a.md_level = None
a.rand_method = None
a.max_var_iteration = 1
a.starting_model = 1
a.ending_model = 1
a.make()
Run the script with the following commands, and re-name the output PDB file:
python 02_build_model.py
mv target.B99990001.pdb receptor_model.pdb
This step should provide you with following file(s):
- receptor_model.pdb
- a PDB file of the receptor built with Modeller without missing residues, and with the correct sequence.
Amber
Finally, we will use Amber to add hydrogen atoms to our model structures. Using Amber will also allow for energy minimization, and molecular dynamics simulations, which can come in handy.
The first step is to generate parameters for the ligand. To do this, the ligand must be correctly protonated. There are several options to this. Babel is one of them:
babel -h -ipdb ligand.pdb -opdb ligand_H.pdb
Alternatively, you can download the ligand (already protonated) with this command:
wget -Oligand_H.pdb http://ligand-expo.rcsb.org/reports/${ligand:0:1}/${ligand}/${ligand}_${model}.pdb
Note: Protonating ligands might not be trivial, and there is no guarantee that these programs yields reasonable protonation.
Once hydrogens are (correctly) added to the ligand, we will use Antechamber to build ligand parameters:
antechamber -i ligand_H.pdb -fi pdb -o prepi -fo prepi -c bcc -s 2 -nc 0
parmchk -i prepi -f prepi -o frcmod
Amber paramter/topolgy, and coordinate files can be generated with the following Leap script:
cat<<EOF> build
source leaprc.ff10
source leaprc.gaff
loadAmberPrep prepi
loadamberparams frcmod
set default PBradii mbondi2
receptor = loadpdb receptor_model.pdb
ligand = loadpdb ligand.pdb
saveamberparm receptor receptor.prmtop receptor.inpcrd
saveamberparm ligand ligand.prmtop ligand.inpcrd
complex = combine { receptor ligand }
saveamberparm complex complex.prmtop complex.inpcrd
EOF
tleap -f build
This step should provide you with following file(s):
- receptor.prmtop receptor.inpcrd
- Amber topology and coordinate files of the receptor.
- complex.prmtop complex.inpcrd
- Amber topology and coordinate files of the receptor-ligand complex.
- ligand.prmtop ligand.inpcrd
- Amber topology and coordinate files of the ligand.
Convert back to PDB files
# Convert to PDB-files
ambpdb -p complex.prmtop < complex.inpcrd > complex.amber.pdb
ambpdb -p receptor.prmtop < receptor.inpcrd > receptor.amber.pdb
# Convert to SYBYL mol2-files
top2mol2 -p receptor.prmtop -c receptor.inpcrd -o receptor.sybyl.mol2
top2mol2 -p ligand.prmtop -c ligand.inpcrd -o ligand.sybyl.mol2
# Fetch the ligand, and save to a seperate PDB file
grep ${ligand} complex.amber.pdb > ligand.amber.pdb
Energy minimization
To get rid of unreasonable hydrogen-hydrogen distances, we will perform a energy minimization using Amber:
CPUS=4
ligand="LL1"
ligand_resi="299"
cat<<EOF>min.in
short minimization
&cntrl
imin=1, maxcyc=1000, ncyc=700,
cut=15., rgbmax=12., ntb=0,
igb=2, saltcon=0.2,
ntr=1, restraintmask=':${ligand_resi}',
restraint_wt=0.5,
ntpr=100
/
EOF
$MPI_HOME/bin/mpirun -np $CPUS $AMBERHOME/bin/sander.MPI -O -i min.in -o min.out \
-p complex.prmtop -c complex.inpcrd -r min.rst -ref complex.inpcrd
# Convert to PDB files
ambpdb -p complex.prmtop < min.rst > min.amb.pdb
grep ${ligand} min.amb.pdb > ${ligand}.amb.pdb
Convert complexed structure to a receptor structure:
cat<<EOF>stripLigand.ptraj
trajin min.rst
strip :${ligand_resi}
trajout receptor.min.amb.rst restart
EOF
ptraj complex.prmtop stripLigand.ptraj
mv receptor.min.amb.rst.1 receptor.min.amb.rst
ambpdb -p receptor.prmtop < receptor.min.amb.rst > receptor.min.amb.pdb
# Convert to SYBYL mol2-files
top2mol2 -p receptor.prmtop -c receptor.min.amb.rst -o receptor.sybyl.mol2
Summary
This description should provide the steps needed to build the reference PDB files needed for INPHARMA experiments. The INPHARMA program will depend on the atom numbering of the files complex.amber.pdb, receptor.amber.pdb, and ligand.amber.pdb.