Re: [AMBER] Setting MM-PBSA script

From: William Lees <>
Date: Wed, 26 Oct 2016 15:14:01 +0100


I started by editing PDB files too, but I found it much easier in the
end to use The only 'trick' with that is that you have
to know the residue numbers that tleap uses (numbered sequentially from
1 throughout the file, without regard to chain breaks).

Here are the commands I used to create an unsolvated complex, and then
three sets of receptor+ligand, using those tleap residue numbers: -p 3gbm_wat.prmtop -c 3gbm_com.prmtop -s :2791-7660,:WAT -p 3gbm_com.prmtop -r 3gbm_HA_1.prmtop -l
3gbm_AB_1.prmtop -n :497-930 -p 3gbm_com.prmtop -r 3gbm_HA_2.prmtop -l
3gbm_AB_2.prmtop -n :1427-1860 -p 3gbm_com.prmtop -r 3gbm_HA_3.prmtop -l
3gbm_AB_3.prmtop -n :2357-2790

If you want to check what residue is at a given position in tleap, you
can use desc. For example to see what is residue 99:

> foo=loadpdb "xyz_refined.pdb"


> desc foo.99
RESIDUE sequence number: 99
RESIDUE PDB sequence number: 99

Alternatively I wrote a rather trivial python script which lists out the
residue numbers in a convenient form. It's called AmberNum and you can
find it at . There are some
other scripts there which you might find useful for MM/PBSA, including
some that are intended to help with the analysis of multiple identical

All the best


On 25/10/2016 14:56, Jag Silwal wrote:
> HI William,
> thank you so much for your input.
> Ideally what I would like to do is with the same MD run I would like to
> calculate ddG for each of the ligands (as you have suggested) and also ddG
> for dimer interface as well. Basically when I got that high ddG for my
> first trial I attempted to use other protein as my ligand. But I ran into
> few problems.
> When I prepared the individual PDB file, I simply edited in a text editor
> and saved each as a separate PDB. For eg. in this case if I have four
> chains, ABCD, then I made PDBs as: A, B, C, D,ABC, ABD, ACD and BCD. But
> when I do this how would I handle the residue and atoms number issue? I got
> error message right at the beginning of the calculation that the residues
> are not sequential. How do I go about renumbering the residues?
> Also, is that the correct way to split the PDBs? It would be great if you
> could share your experience.
> the last simulation I ran was with ABC as a ligand and D as receptor and at
> least the MM-PBSA calculation was finished because residue number were
> sequential in this case. I checked all the residues and atoms for each
> individual PDBs and they all look ok to me. So I don't understand why the
> ddG came up that high.
> If you could share some details on how you handled your trimer I would
> really appreciate your help.
> Thanks again for your help
> Sincerely,
> Jag
> On Tue, Oct 25, 2016 at 6:12 AM, William Lees <> wrote:
>> Just for info, I've been working with a trimer to which three ligands
>> bind, and I can confirm that MM/PBSA and MM/GBSA analysis works as
>> expected with this kind of setup. In the longer term, you'll probably
>> want to run the calculations against both of the ligands so that you get
>> a 'two for one' benefit from your trajectiories.
>> Is it possible that you have defined the receptor/ligand boundaries
>> slightly incorrectly, so that one or two ligand residues are in the
>> receptor, or vice versa? I have found this easy to do and it does give
>> rise to large delta delta Gs.
>> All the best
>> William
>> On 24/10/2016 19:03, Jag Silwal wrote:
>>> Dear all,
>>> I ran a simulation for protein complexes which includes four chains as
>>> shown below. Protein "C" is a dimer and protein A binds to protein C and
>> I
>>> am interested to run MM-PBSA to estimate binding free energy between
>>> protein A and protein C in its dimeric form.
>>> Fo MM-PBSA, I created a separate PDB file for A and another PDB file with
>>> the rest of the proteins (C-C dimer +another A). Then for the MMPBSA
>>> command I treated A as the ligand protein and C-C dimer + A as receptor
>>> protein. I ran with the following script:
>>> &general
>>> |startframe=4501, endframe=5550, interval=5,
>>> |verbose=2, keep_files=2, strip_mask=':WAT,CL',
>>> |/
>>> |&gb
>>> |igb=8, saltcon=0.150,
>>> |/
>>> |&pb
>>> |inp=1, radiopt=0, istrng=0.15, fillratio=4.0
>>> |/
>>> When the run was done I checked the Delta delta G and it came to be about
>>> 1000 kcal/mol which I think is way too high for the interface.
>>> Previously, I have separately ran and analyzed the Delta delta G without
>>> the dimer where there is only one A and one C the Delta delta G was -60
>>> kcal/mol.
>>> So I am sure I am not either setting these different PDBs right or may
>> be I
>>> am missing something important while setting up for MM-PB/GBSA analysis.
>>> How do I approach this problem?
>>> Any insight would be really appreciated.
>>> [image: Inline image 2]
>>> Sincerely,
>>> Jag
>>> Graduate Student,
>>> Michigan State University
>>> Department of Chemistry
>>> _______________________________________________
>>> AMBER mailing list
>> _______________________________________________
>> AMBER mailing list
> _______________________________________________
> AMBER mailing list

AMBER mailing list
Received on Wed Oct 26 2016 - 07:30:03 PDT
Custom Search