[gmx-developers] free energy calculations with restraints

[David Mobley]

Yes, I'm very enthusiastic about adding defined intermolecular restraints
using absolute atom indices. It seems like a straightforward, correct,
easy-to-understand way to deal with this simple problem.

David



On Tue, Aug 20, 2013 at 9:39 AM, Shirts, Michael (mrs5pt) <
mrs5pt at eservices.virginia.edu> wrote:

>
> The other thing I would say is that the couple-mol code is a little
> difficult to work with----there's known quirks, like handling the
> intramolecular dispersion in a way that's slightly off from free energy
> off,
> as well as some new quirks that I'll try to post about later this
> week---that it would actually be cleaner to add clearly defined
> intermolecular restraints than add complexity to the coupling feature.
>
> Best,
> ~~~~~~~~~~~~
> Michael Shirts
> Assistant Professor
> Department of Chemical Engineering
> University of Virginia
> michael.shirts at virginia.edu
> (434)-243-1821
>
>
> > From: Floris Buelens <floris_buelens at yahoo.com>
> > Reply-To: Floris Buelens <floris_buelens at yahoo.com>, Discussion list for
> > GROMACS development <gmx-developers at gromacs.org>
> > Date: Tue, 20 Aug 2013 08:35:15 -0400
> > To: Berk Hess <hess at kth.se>, Discussion list for GROMACS development
> > <gmx-developers at gromacs.org>
> > Subject: Re: [gmx-developers] free energy calculations with restraints
> >
> > Hi Berk,
> >
> > I think we might have gone out of sync there - I sent the message below
> before
> > reading yours from this morning. From the misplaced context it maybe
> sounded
> > like I was disagreeing with your '[ intermolecular-interactions ]'
> suggestion,
> > which I'm definitely not - this would clearly be preferable and probably
> > useful to a wider audience.
> >
> >
> >
> > ----- Original Message -----
> > From: Berk Hess <hess at kth.se>
> > To: Floris Buelens <floris_buelens at yahoo.com>; Discussion list for
> GROMACS
> > development <gmx-developers at gromacs.org>
> > Cc:
> > Sent: Tuesday, 20 August 2013, 12:04
> > Subject: Re: [gmx-developers] free energy calculations with restraints
> >
> > Hi,
> >
> > We could consider this.
> > I don't know how you did this now, but this would need an extra mdp
> > option, something like couple-group, to avoid one option with a moltype
> > or index group as an input.
> > Then we would also need to decide what to do with a group which is part
> > of one molecule of a certain moleculetype with multiple copies. Disallow
> > this? Or duplicate the moltype and only modify one copy (as is done now
> > for qmmm)?
> >
> > Cheers,
> >
> > Berk
> >
> > On 08/20/2013 09:24 AM, Floris Buelens wrote:
> >> Hi,
> >>
> >> I would also still argue for augmenting the current couple-moltype
> >> functionality with this modification. I can see this is a niche
> requirement,
> >> which of course does speak against inclusion. But in my view the
> replacement
> >> code is no less elegant than the current implementation (modification
> of five
> >> self-contained functions at the pre-processing stage, maybe 50-odd extra
> >> lines), we can maintain full backward compatibility, and a single simple
> >> check will remove any potential for accidental misuse.
> >>
> >>
> >>
> >> ________________________________
> >> From: David Mobley <dmobley at uci.edu>
> >> To: Discussion list for GROMACS development <gmx-developers at gromacs.org
> >
> >> Cc: Floris Buelens <floris_buelens at yahoo.com>; Michael R. Shirts
> >> <michael.shirts at virginia.edu>
> >> Sent: Monday, 19 August 2013, 20:02
> >> Subject: Re: [gmx-developers] free energy calculations with restraints
> >>
> >>
> >>
> >> Hi, Berk,
> >>
> >> I agree that this solution would be the "proper" one and would be nice.
> >> However, it's been a long time coming (this is an issue I've been
> raising for
> >> approximately 5 years, and there's no progress; normally no one even
> >> answers). Currently, for the restraints I need to use, I *have* to
> merge the
> >> ligand into my protein topology, inconvenient or not. Doing this also
> means
> >> there are some features I simply can't use (couple-moltype for example).
> >>
> >> It strikes me that Floris's solution is a good interim one, in that it
> at
> >> least solves most of the typical use cases we're facing, and has the
> >> advantage of being *already implemented*.
> >>
> >> Is there any way we could get this into the main GROMACS, at least until
> >> someone gets time to implementing the "proper" solution (which could be
> years
> >> from now)?
> >>
> >> Thanks!
> >> David
> >>
> >>
> >>
> >>
> >>
> >> On Mon, Aug 19, 2013 at 1:10 AM, Berk Hess <hess at kth.se> wrote:
> >>
> >> Hi,
> >>> I fully understood this point, maybe my reply wasn't well structured.
> >>> I was trying to argue that the proper solution would be to implement
> >>> inter-molecular restraints, instead of introducing an option which
> could me
> >>> misused. Also merging a ligand topology into your protein topology is
> very
> >>> inconvenient. This problem with the proper solution is, of course, that
> >>> someone will need to implement inter-molecular restraints/potentials.
> >>>
> >>> One issue used to be that we corrected molecules for PBC before
> calculating
> >>> bonded interactions. But in 4.6 this is usually not faster and not
> used any
> >>> more. That makes it easier to treat intra- and inter-molecular
> interactions
> >>> the same way at the mdrun level.
> >>>
> >>> Cheers,
> >>>
> >>> Berk
> >>>
> >>>
> >>> On 08/19/2013 09:55 AM, Floris Buelens wrote:
> >>>
> >>> Hi Berk,
> >>>> I think some clarification is needed. What we're talking about is only
> >>>> relevant in the context of non-bonded interactions, typically between
> a
> >>>> small molecule and a protein. To access the binding affinity of a
> ligand
> >>>> through alchemical methods, it's useful to switch off only the
> interactions
> >>>> of the ligand with its environment while maintaining intra-ligand
> >>>> potentials - this is what's provided by the 'couple-moltype' code
> path.
> >>>>
> >>>> The limitation that we're trying to circumvent comes from the fact
> that as
> >>>> you scale down interactions with the environment, the ligand is no
> longer
> >>>> held in the binding site. To counteract this, it's necessary to
> perturb in
> >>>> restraining potentials as the intramolecular nonbonded interactions
> are
> >>>> perturbed out. The most practical method makes use of a single
> distance
> >>>> restraint, two angle and three dihedral restraints, whose effect can
> later
> >>>> be accounted for analytically.
> >>>>
> >>>> The current code only allows decoupling of a whole logical 'molecule'
> as
> >>>> specified in the topology file. This precludes the use of regular
> (fully
> >>>> perturbation-aware) bonded interactions. As far as I'm aware, the
> pull code
> >>>> doesn't provide what's required. Inter-molecular bonded interactions
> would
> >>>> be a great general solution but presumably won't show up any time
> soon.
> >>>>
> >>>> The workaround here is instead to represent two physical molecules
> (e.g.
> >>>> protein and ligand) as a single logical molecule (in the topology
> file), so
> >>>> regular bonded potentials can be applied. Current Gromacs doesn't
> allow the
> >>>> decoupling ('couple-moltype') code to be used in this scenario. My
> >>>> modification allows you to target the ligand by its residue name and
> get
> >>>> this functionality back.
> >>>>
> >>>> Decoupling a single residue of a multi-residue chain is indeed
> probably not
> >>>> correct in this framework (I did call it 'weird' in my last message
> :-) ).
> >>>> An extra check would block users from trying this.
> >>>>
> >>>> I hope that clarifies the problem we're trying to solve. I agree this
> will
> >>>> be useful a relatively small number of users, but on the other hand
> it's a
> >>>> very unobtrusive, self-contained modification which can maintain full
> >>>> backwards compatibility.
> >>>>
> >>>> thanks,
> >>>>
> >>>> Floris
> >>>>
> >>>>
> >>>> ----- Original Message -----
> >>>> From: Berk Hess <hess at kth.se>
> >>>> To: Floris Buelens <floris_buelens at yahoo.com>; Discussion list for
> GROMACS
> >>>> development <gmx-developers at gromacs.org>
> >>>> Cc:
> >>>> Sent: Monday, 19 August 2013, 9:03
> >>>> Subject: Re: [gmx-developers] free energy calculations with restraints
> >>>>
> >>>> Hi,
> >>>>
> >>>> There is the rotational pull code, but that might not provide the
> exact
> >>>> functionality you want.
> >>>> Already for some time we have been discussing inter-molecular
> >>>> interactions, by specifying molecule type, molecule index and atom
> >>>> index, but there are no concrete plans for implementing this yet.
> >>>>
> >>>> An option for decoupling part of a molecule indeed sound useful. But
> in
> >>>> practice you always need to replace that part by something else, at
> >>>> least a hydrogen, and modify some potentials of connecting atoms, so I
> >>>> don't know how generally useful such an option is.
> >>>>
> >>>> Cheers,
> >>>>
> >>>> Berk
> >>>>
> >>>> On 08/19/2013 07:56 AM, Floris Buelens wrote:
> >>>>
> >>>> I suppose that would make sense. The changes are fairly minor (more
> or less
> >>>> just the function convert_moltype_couple and the four functions it
> calls)
> >>>> and mainly consist of a bit of extra juggling with nonbonded
> exclusions
> >>>> during preprocessing. The current functionality (decouple a whole
> molecule)
> >>>> could be handled as a special case of the new code. A check for bonds
> to
> >>>> the rest of the structure could stop people from trying weird stuff
> like
> >>>> decoupling residues from a chain.
> >>>>>
> >>>>>
> >>>>>
> >>>>> ________________________________
> >>>>> From: David Mobley <dmobley at gmail.com>
> >>>>> To: Floris Buelens <floris_buelens at yahoo.com>
> >>>>> Cc: Discussion list for GROMACS development <
> gmx-developers at gromacs.org>
> >>>>> Sent: Friday, 16 August 2013, 18:30
> >>>>> Subject: Re: [gmx-developers] free energy calculations with
> restraints
> >>>>>
> >>>>>
> >>>>>
> >>>>> This does sound useful, though it would be more useful if this could
> be
> >>>>> implemented into the main GROMACS rather than a separate code (since
> >>>>> otherwise it will go away as GROMACS is further developed).
> >>>>>
> >>>>> Would this be a possibility going forward, devels?
> >>>>>
> >>>>> Thanks!
> >>>>>
> >>>>>
> >>>>>
> >>>>>
> >>>>> On Fri, Aug 16, 2013 at 3:07 AM, Floris Buelens <
> floris_buelens at yahoo.com>
> >>>>> wrote:
> >>>>>
> >>>>> Hi David,
> >>>>>
> >>>>> I have the same requirement as you, but I've gone about it slightly
> >>>>> differently. I've hacked the couple-moltype code path to allow
> decoupling
> >>>>> of a specific residue (identified by name) instead of a molecule.
> This
> >>>>> allows the residue of interest to be part of another molecule block
> so you
> >>>>> can set up distance, angle and dihedral restraints using regular
> bonded
> >>>>> interactions with full perturbation support.
> >>>>>> If this is useful to you or to anyone else, let me know and I'll be
> happy
> >>>>>> to share.
> >>>>>>
> >>>>>> best,
> >>>>>>
> >>>>>> Floris
> >>>>>>
> >>>>>> ________________________________
> >>>>>> From: David Mobley <dmobley at gmail.com>
> >>>>>>
> >>>>>> To: Discussion list for GROMACS development <
> gmx-developers at gromacs.org>
> >>>>>> Sent: Friday, 14 June 2013, 21:47
> >>>>>> Subject: [gmx-developers] free energy calculations with restraints
> >>>>>>
> >>>>>>
> >>>>>>
> >>>>>>
> >>>>>> Hi, devs,
> >>>>>>
> >>>>>>
> >>>>>> I'm writing with an issue relating to the interplay of new free
> energy
> >>>>>> features with restraints.
> >>>>>>
> >>>>>> I'm very much appreciating some of the new free energy features in
> >>>>>> gromacs, such as the 'couple-moltype' option which provides a way
> to set
> >>>>>> up decoupling or annihilation of a specific molecule via free energy
> >>>>>> calculations without having to edit the topology file directly.
> This is
> >>>>>> especially great when it comes to decoupling -- charge decoupling
> was not
> >>>>>> previously possible via topology file editing, and vdW decoupling
> took
> >>>>>> substantial manipulation of the topology file.
> >>>>>>
> >>>>>> However, for binding free energies, my work employs orientational
> >>>>>> restraints between the protein and ligand. I need to be able to
> impose
> >>>>>> both dihedral and angle restraints on angles between the protein and
> >>>>>> ligand. Currently, I do this using angle-restraints and
> >>>>>> dihedral-restraints. This requires that both the protein and ligand
> be
> >>>>>> within the same logical 'molecule', which (unfortunately) means
> that I
> >>>>>> can't make use of the new free energy features above, since
> >>>>>> couple-moltype has to apply to a whole molecule, not just part of a
> >>>>>> molecule.
> >>>>>>
> >>>>>> So, my I see two possible solutions to the problem, and hence have
> these
> >>>>>> questions:
> >>>>>> 1) Can dihedral and angle restraints be applied via the pull code?
> If
> >>>>>> not, are there any plans to add that?
> >>>>>> 2) Alternatively, what about modifying the restraints code so it
> uses (or
> >>>>>> at least optionally allows) absolute atom numbering, rather than
> >>>>>> numbering within a specific molecule, thus allowing restraints
> >>>>>> (angle-restraints and dihedral-restraints) to be applied between
> >>>>>> 'molecules'?
> >>>>>>
> >>>>>> Thanks!
> >>>>>> David
> >>>>>>
> >>>>>>
> >>>>>>
> >>>>>> --
> >>>>>> David Mobley
> >>>>>> dmobley at gmail.com
> >>>>>> 949-385-2436
> >>>>>>
> >>>>>>
> >>>>>> --
> >>>>>> gmx-developers mailing list
> >>>>>> gmx-developers at gromacs.org
> >>>>>> http://lists.gromacs.org/mailman/listinfo/gmx-developers
> >>>>>> Please don't post (un)subscribe requests to the list. Use the
> >>>>>> www interface or send it to gmx-developers-request at gromacs.org.
> >>>>>>
> >>>>>>
> >>> --
> >>> gmx-developers mailing list
> >>> gmx-developers at gromacs.org
> >>> http://lists.gromacs.org/mailman/listinfo/gmx-developers
> >>> Please don't post (un)subscribe requests to the list. Use the www
> interface
> >>> or send it to gmx-developers-request at gromacs.org.
> >>>
> >>
> > --
> > gmx-developers mailing list
> > gmx-developers at gromacs.org
> > http://lists.gromacs.org/mailman/listinfo/gmx-developers
> > Please don't post (un)subscribe requests to the list. Use the
> > www interface or send it to gmx-developers-request at gromacs.org.
>
> --
> gmx-developers mailing list
> gmx-developers at gromacs.org
> http://lists.gromacs.org/mailman/listinfo/gmx-developers
> Please don't post (un)subscribe requests to the list. Use the
> www interface or send it to gmx-developers-request at gromacs.org.
>



-- 
David Mobley
dmobley at gmail.com
949-385-2436
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://maillist.sys.kth.se/pipermail/gromacs.org_gmx-developers/attachments/20130820/44912d00/attachment.html>