Hi Olivia,
I just realize that my solution will not work the way you plan to run your TI simulations. You obviously want to run "full" dual topology i.e, each residue is an independent copy (own set of coordinates) where one needs to be completely annihiliated (transformed to "nothing") while simultaneously the other is created from "nothing". So you are "absolutely" right to compare this to "absolute" FE calculation because you are doing effectively two of those at the same time.
If you wish to procede with separating electrostatics from vdW you will have to completely discharge both PHE and ALA (reverse the sign as you actually want recharge or run it in the opposite direction) in a linear transformation and then transform the discharged residued into each other. I would think of this as a 3-step protocol: 1) PHE -> PHE(0) 2) PHE(0) -> ALA(0) 3) ALA(0) -> ALA.
The softcore implementation in sander/pmemd is for _both_ electrostatic and vdW softcore functions and on for both by default, so you could try and use this and transform in a single step. I had some issue with this for certain relative calculations but for some tests with absolute FE it seemed to be fine. If you wanted to test this you can compare it to the 3-step protocol above, but only compare "double deltas" i.e. the ddGs from completely closed therodynamic cycles.
To explain what I am talking about in the tutorial and in my previous emails: You can "map" equivalent atoms in a transformation pair to each other which means that those atoms share coordinates and the only transformations that take place are between atoms types which can be done linearly. In the example of ALA->VAL the alanine can be fitted completely into VAL allowing for the 2-step protocol explained earlier. In PHE->ALA you would not want to break the ring and so could only map the backbone atoms. But thanks for pointing me to this. I have to clarify this in the tutorial to avoid further confusion.
Cheers,
Hannes.
________________________________________
From: Olivia Pierce [olivijuly.gmail.com]
Sent: 18 July 2015 00:36
To: AMBER Mailing List
Subject: Re: [AMBER] Mutation of a protein residue in TI computation with pmemd.MPI
Hi Hannes,
Thanks so much for clarification. I have one last question - just to make
sure that I understood everything correctly and got my input files are what
they should be.
My understanding for your response on 2) in the previous email, "So you
have to do the PHE->INT transformation and the INT->ALA transformation
separately" is that each of these two simulations, PHE -> INT and INT ->
ALA, is JUST SOFT CORE transformation (i.e., vdW transformation) because 1)
crgmask doesn't "transform" the charge from one to another and 2) INT
already has the atomic charges of ALA for the atoms corresponding to ALA. Am
I right?
So the TI part of my input files look like:
icfe = 1, clambda = $i (this will be my lambda),
ifsc = 1,
timask1 = ':65', timask2 = ':178',
scmask1 = ':65', scmask2 = ':178',
where the residue 65 is PHE and 178 is INT residue. However, what I am
still not clear about is whether this soft core simulation will transform
the atomic charges of PHE to those of INT, together with the vdW radii
(even though I called this soft core transformation above). Could you
clarify me?
I could think of doing two step transformations (charge and vdW separately)
for both PHE -> INT and INT -> ALA legs but that will be essentially same
as PHE -> nothing -> ALA. My mind is still viewing this type of computation
from the perspective of absolute FE calculation so please correct/clarify
me if I am wrong in any place.
Thank you very much, really appreciate your help.
Best,
Olivia
On Fri, Jul 17, 2015 at 4:35 PM, <hannes.loeffler.stfc.ac.uk> wrote:
> Hi Olivia,
>
> good that you found the tutorial. It has been online only since today...
> :)
>
> 1) In relative calculations you need to do the two steps according to what
> state you have the disappearing/appearing atoms in. In fact, this can only
> work if you have those in one end-state but not in both. (In the latter
> case you would have to work with partial or complete de/recharge as
> faciliated with cgrmask and exemplified in the ligand part of tutorial
> A9). Yes, it is about the abnormalities you refer two. You just need to
> realize the ALA->VAL is the setup in the opposite direction i,e, while in
> PHE->ALA you have disappearing atoms in PHE, ALA->VAL has _appearing_
> atoms. Look at the figure and turn it around and this should become clear.
>
> 2) When you set up the INT state you explicitly set the charges. So there
> is no need to use crgmask. You could only use it to set the charges zero
> for the disappearing atoms anyway. There is currently no way in the code
> to only transform certain force field parameters (Gromacs calls this
> "lambda-paths") as you describe. This can only be emulated through careful
> setup of topology files. I have asked for this feature to the developers
> (I think Dan would be doing this) but I do not know if they have the time
> to pick up on this. From a user perspective this would certainly make this
> kind of transformations easier. So you have to do the PHE->INT
> transformation and the INT->ALA transformation separately.
>
> Cheers,
> Hannes,
>
> ________________________________________
> From: Olivia Pierce [olivijuly.gmail.com]
> Sent: 17 July 2015 20:57
> To: AMBER Mailing List
> Subject: Re: [AMBER] Mutation of a protein residue in TI computation with
> pmemd.MPI
>
> Hi Hannes,
>
> Thank you so much! I actually found your tutorial A9 for Amber as well
> (while the transformation therein is Ala -> Val). I have a couple more
> questions though....
>
> 1) In the "Multi-stepping" part of your "Side chain mini" tutorial (
> http://ambermd.org/tutorials/advanced/tutorial9/#sidechain_mini), you
> suggested doing vdW transformation prior to the charge transformation. I've
> always done (absolute binding free energy) TI calculation with the charge
> transformation first and then vdW transformation - mainly in order to
> remove any possibility of any abnormality when remove the vdW terms while
> atomic charges are still there. However, should it be the other way (charge
> first, then vdW) in this case? Is there any reason?
>
> 2) So I created the "INT" residue that has the atomic charges of ALA but
> bond terms for PHE (therefore the disappearing dummy atoms derived from
> PHE). However, the 'crgmask' command removes the entire charges of the
> specified atoms. And there is nothing equivalent to scmask1 & scmask2 for
> crgmask. Is there any way to transform charges from PHE to ALA? Or do I not
> need this charge transformation step at all because I created the INT
> residue that has the atomic charges of ALA already? Or still need to
> transform the charges of PHE to INT?
>
> I appreciate your help!
>
> Best,
> Olivia
>
> On Fri, Jul 17, 2015 at 1:22 PM, <hannes.loeffler.stfc.ac.uk> wrote:
>
> > Hi Olivia,
> >
> > so, you are going for the 2-step protocol?
> >
> > In that case you can create the intermediate from the regular PHE
> > (depending on your force field) and simply replace the charges with those
> > of ALA. The disappearing atoms ("dummies") would have zero charges but
> > still have the same bonded and vdW terms as PHE so they are
> "electrostatic
> > dummies" only. One practical approach is to create a template of your
> > own. The other option that comes to mind is to do some post-processing
> > with parmed but I haven't looked into that yet.
> >
> > I have attached two files which are part of the upcoming A9 tutorial on
> > side-chain mutations. You would need to run the .cmd (it's a leap
> script)
> > first to create INT.lib. Next you loadoff INT.lib before you read a
> > modified PDB file where your PHE (or ALA) is replaced with the new INT
> > residue. So for pmemd you would create prmtops for PHE/INT and INT/ALA.
> > First step is a linear transformation (only charges change) and for the
> > second step you would define the "dummies" as softcore atoms.
> >
> > Please note that this example is for VAL->ALA and the atoms in VAL have
> > been reordered because I am also playing with mappings. In case you
> can't
> > open the attachments I have also pasted their contents below.
> >
> > Cheers,
> > Hannes.
> >
> > -----
> >
> > # create a template for VAL maximally mapped to ALA
> > # NOTE: force field dependent, here FF14SB
> >
> > source leaprc.ff14SB
> > INT = loadmol2 INT.mol2
> >
> > set INT restype protein
> > set INT head INT.1.N
> > set INT tail INT.1.C
> > set INT.1 connect0 INT.1.N
> > set INT.1 connect1 INT.1.C
> >
> > saveoff INT INT.lib
> >
> > quit
> >
> > -----
> >
> >
> > .<TRIPOS>MOLECULE
> > INT
> > 16 15 0 0 0
> > PROTEIN
> > USER_CHARGES
> > ****
> > VAL template from ff14SB, charges as in ALA with zero for dummies
> > .<TRIPOS>ATOM
> > 1 N 15.1767 16.7230 17.6889 N 1 INT -0.4157
> > 2 H 15.0727 15.7440 17.9149 H 1 INT 0.2719
> > 3 CA 16.1687 17.0960 16.7279 CX 1 INT 0.0337
> > 4 HA 16.8367 17.8110 17.2079 H1 1 INT 0.0823
> > 5 CB 16.9957 15.8960 16.2529 3C 1 INT -0.1825
> > 6 HB 16.3417 15.1760 15.7619 HC 1 INT 0.0603
> > 7 CG1 17.7567 16.2340 15.5499 CT 1 INT 0.0603
> > 8 CG2 17.4767 15.4240 17.1089 CT 1 INT 0.0603
> > 9 C 15.5647 17.7420 15.5159 C 1 INT 0.5973
> > 10 O 16.1667 18.6510 14.9469 O 1 INT -0.5679
> > 11 HG11 18.4111 16.9550 16.0409 HC 1 INT 0.0000
> > 12 HG12 18.3449 15.3804 15.2125 HC 1 INT 0.0000
> > 13 HG13 17.2759 16.7064 14.6931 HC 1 INT 0.0000
> > 14 HG21 16.7160 15.0856 17.8124 HC 1 INT 0.0000
> > 15 HG22 18.0642 14.5708 16.7713 HC 1 INT 0.0000
> > 16 HG23 18.1307 16.1444 17.6014 HC 1 INT 0.0000
> > .<TRIPOS>BOND
> > 0 1 2 1
> > 1 1 3 1
> > 2 3 4 1
> > 3 3 5 1
> > 4 3 9 1
> > 5 5 6 1
> > 6 5 7 1
> > 7 5 8 1
> > 8 7 11 1
> > 9 7 12 1
> > 10 7 13 1
> > 11 8 14 1
> > 12 8 15 1
> > 13 8 16 1
> > 14 9 10 1
> >
> > ________________________________________
> > From: Olivia Pierce [olivijuly.gmail.com]
> > Sent: 17 July 2015 17:52
> > To: AMBER Mailing List
> > Subject: Re: [AMBER] Mutation of a protein residue in TI computation with
> > pmemd.MPI
> >
> > Hi,
> >
> > Thanks Dan and Hannes for your comments! I also thought about placing the
> > dummy atoms into the places of the disappearing atoms of Phe, but the
> > problem I thought was that tleap might complain - the mutating residue
> name
> > is Ala but obviously it has more atoms than what Ala should have. Do you
> > think this would not be the problem? I would appreciate any further
> > suggestions.
> >
> > Thanks!
> > Olivia
> >
> > On Fri, Jul 17, 2015 at 2:30 AM, <hannes.loeffler.stfc.ac.uk> wrote:
> >
> > > Dan has already told you how to principally solve this. But the
> > > description of what you are planning to do does not seem to match your
> > > input file. I understand that you want to do an elestrostatic
> > > transformation only first but what crgmask does is to only set the atom
> > > charges to zero. Apparently you switch off the charges of ALA
> > completely.
> > > If that is really your plan you would do a linear transformation of
> > > ALA(charge) to ALA(uncharged). Although if you mutate PHE->ALA the
> steps
> > > would seem to be: 1) PHE(charged)->PHE(uncharged); 2)
> > > PHE(uncharged)->ALA(uncharged); 3) ALA(uncharged)->ALA(charged).
> > >
> > > An electrostatic only transformation would (currently) require you to
> > > create and intermediate that has all vdW and bonded parameters of PHE
> but
> > > the charges of ALA. As you would do that in a linear transformation
> you
> > > would need to fill up all disappearing atoms of PHE with dummies in
> ALA.
> > > These dummies are then defined as softcore atoms in the second step.
> > > Overall, this protocol seems to be more attractive to me but requires
> > more
> > > work to set up (I am writing a tool to do exactly this but it is still
> > only
> > > dealing with non-covalently bound molecules yet). Alternatively, pmemd
> > > would implement lambda paths which would make this much simpler. Dan,
> is
> > > this still on your TODO list?
> > >
> > > ________________________________________
> > > From: M Olivia Kim [olivijuly.gmail.com]
> > > Sent: 16 July 2015 23:13
> > > To: amber.ambermd.org
> > > Subject: [AMBER] Mutation of a protein residue in TI computation with
> > > pmemd.MPI
> > >
> > > Hi,
> > >
> > > I am trying to run a FEP simulation for mutation of a protein residue
> > > using the pmemd.MPI module. I'm mutating a Phe to Ala - following the
> > > Amber14 manual, I prepared the prmtop and inpcrd files using parmed.py
> > > where the redundant bonding terms for the overlapping residues are
> > deleted.
> > > So the final pdb for my system looks like that protein consists of
> > residues
> > > 1-177 and the residues 178 is Ala (separated by TER), which is the
> mutant
> > > for the Phe in the WT.
> > >
> > > I was first trying to run the non-soft core simulation, where the
> atomic
> > > charges of the WT protein residue are changed to the charges of the
> > mutated
> > > residue along with lambda. However, since the atom numbers of the WT
> > > residue and mutant are different, I couldn't even run the minimization.
> > The
> > > part corresponding to TI in my minimization input file is:
> > >
> > > icfe = 1, clambda = 0.0,
> > > ifsc = 0,
> > > timask1 =
> > >
> >
> '.1019,1020,1021,1022,1023,1024,1025,1026,1027,1028,1029,1030,1031,1032,1033,1034,1035,1036,1037,1038',
> > > timask2 = '.2780,2781,2782,2783,2784,2785,2786,2787,2788,2789',
> > > crgmask = '.2780,2781,2782,2783,2784,2785,2786,2787,2788,2789',
> > > /
> > >
> > > And the error message I got is:
> > >
> > > TI Mask 1
> > >
> >
> .1019,1020,1021,1022,1023,1024,1025,1026,1027,1028,1029,1030,1031,1032,1033,1034,1035,1036,1037,1038;
> > > matches 20 atoms
> > > TI Mask 2 .2780,2781,2782,2783,2784,2785,2786,2787,2788,2789; matches
> > > 10 atoms
> > > ...
> > > ERROR: timask1/2 must match the same number of atoms for non-softcore
> run
> > >
> > > Is there any way I can correct my input or prmtop file so that I can
> run
> > > this simulation? Or is this type of simulation (mutation of a protein
> > > residue to another that has different number of atoms) not supported by
> > > pmemd version of TI in Amber yet? I'll look forward to any advice or
> > > comments. Thanks in advance!
> > >
> > > Best,
> > > Olivia
> > > _______________________________________________
> > > AMBER mailing list
> > > AMBER.ambermd.org
> > > http://lists.ambermd.org/mailman/listinfo/amber
> > >
> > > _______________________________________________
> > > AMBER mailing list
> > > AMBER.ambermd.org
> > > http://lists.ambermd.org/mailman/listinfo/amber
> > >
> > _______________________________________________
> > AMBER mailing list
> > AMBER.ambermd.org
> > http://lists.ambermd.org/mailman/listinfo/amber
> >
> > _______________________________________________
> > AMBER mailing list
> > AMBER.ambermd.org
> > http://lists.ambermd.org/mailman/listinfo/amber
> >
> >
> _______________________________________________
> AMBER mailing list
> AMBER.ambermd.org
> http://lists.ambermd.org/mailman/listinfo/amber
>
> _______________________________________________
> AMBER mailing list
> AMBER.ambermd.org
> http://lists.ambermd.org/mailman/listinfo/amber
>
_______________________________________________
AMBER mailing list
AMBER.ambermd.org
http://lists.ambermd.org/mailman/listinfo/amber
_______________________________________________
AMBER mailing list
AMBER.ambermd.org
http://lists.ambermd.org/mailman/listinfo/amber
Received on Sat Jul 18 2015 - 02:30:02 PDT
