Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Very large (creating entire rings) proposals #53

Closed
pgrinaway opened this issue Dec 8, 2015 · 18 comments
Closed

Very large (creating entire rings) proposals #53

pgrinaway opened this issue Dec 8, 2015 · 18 comments
Assignees

Comments

@pgrinaway
Copy link
Member

So, I just tried a very large proposal of transforming erlotinib into gefitinib. A couple things:

  1. While it didn't die in dynamics (1000 steps of verlet at 300K, no NCMC), the potential did get much larger (~1e12 kJ/mol)
  2. The aromatic ring is really messed up, presumably due to the need for impropers.
  3. The non-aromatic ring is somewhat less messed up, though could definitely use a bit of help.

For viewing pleasure, the results are available here. The configuration after dynamics looks pretty messed up, and may be on its way to NaNdyland.

@pgrinaway pgrinaway mentioned this issue Dec 8, 2015
@jchodera
Copy link
Member

jchodera commented Dec 8, 2015

+1 for naming this file "gefitini". This will definitely be a drink on Friday's menu.

@jchodera
Copy link
Member

jchodera commented Dec 8, 2015

Ring creation will definitely need some TLC. We may need to do more than just consider impropers.

@pgrinaway
Copy link
Member Author

+1 for naming this file "gefitini".

I wish I had done that on purpose...

Ring creation will definitely need some TLC. We may need to do more than just consider impropers.

Yeah, you're probably right about that.

@jchodera
Copy link
Member

jchodera commented Dec 8, 2015

The big question here: If we add impropers, are there still multiple ways the bonds around the ring could go so that the ring doesn't close properly, or is there only one maximum in the torsion p(phi) after the first bond of the ring is placed?

@jchodera
Copy link
Member

jchodera commented Dec 8, 2015

The best model system to start to test on is probably benzene to biphenyl, since it involves the creation of a phenyl ring onto another phenyl ring.

@jchodera
Copy link
Member

jchodera commented Dec 9, 2015

Without adding in improper torsions, I did a quick test of benzene to biphenyl using test_run_geometry_engine() from perses/tests/test_geometry_engine.py, and found that it clearly generated three distinct behaviors:
biphenyl
For a six-membered planar ring, our proposal scheme has several ways it could go wrong. I believe that creating a properly fused six-membered planar ring requires the random torsion selection make two essentially equiprobable decisions a particular way, resulting in a success rate of maybe 25%.

@jchodera
Copy link
Member

jchodera commented Dec 9, 2015

I think we might benefit from looking into ensemble chain-growth or Rosenbluth-type methods for growing a good proposal configuration.

Here's a review of pruned-enriched Rosenbluth methods: http://arxiv.org/abs/1107.1214

The difficulty would be once we start using ensemble-growth and/or pruning-type methods, it may be difficult to compute the Hastings proposal ratio for the Metropolis-Hastings acceptance criteria. Some of these methods track weights or an estimate of the partition function for the generated ensemble, so it is possible this could provide what we need.

@jchodera
Copy link
Member

jchodera commented Dec 9, 2015

Another possibility is to use parallel Monte Carlo to generate the geometry proposal before relaxing it with NCMC. If we only used the valence potential energy for this part, we that would remove the contribution of valence energy terms from the overall NCMC + RJMC + NCMC log acceptance probability, and still allow us to use NCMC to anneal the other contributions from intermolecular interactions.

@jchodera
Copy link
Member

jchodera commented Dec 9, 2015

Another idea: currently, there are no potential terms to keep rings in a ring like shake because they are not needed when rings are covalently closed. We could instead introduce a guide field potential that functions solely to prevent rings from being created in a malformed geometry like two of the three examples above. It should really only kick in for the malformed cases, and should still be finite everywhere, but sufficiently large to greatly increase the yield of properly closed rings.

It might require some brainstorming and analysis of MCSS pairings to figure out the simplest way to implement this, but it would be very simple from an implementation and algorithmic point of view.

@pgrinaway
Copy link
Member Author

Thanks for the suggestions--a lot of interesting stuff! Going through the PERM review now, and gathering some thoughts on this...

@jchodera
Copy link
Member

Lots and lots of great stuff on the related configurational-bias Monte Carlo (CBMC) here:
http://towhee.sourceforge.net/algorithm/cbmc.html

@jchodera
Copy link
Member

Some slides with an overview of methods here: http://www.acmm.nl/molsim/molsim2011/Presentations/Lecture8.pdf

@jchodera
Copy link
Member

From the Towhee CBMC page:

Fixed Endpoint CBMC
Cyclic molecules are substantially more difficult to grow using CBMC because their conformational space is severely limited by the constraints of having cyclic portions of the molecule. Attempting to grow a cyclic molecule using standard CBMC methods, and just hoping it closes itself up properly, has an acceptance rate that is nearly zero. It is generally believed that an additional biasing is required during the growth procedure in order to nudge the growth into positions that will result in reasonable ring closures. There are a variety of biasing procedures in the literature. One that is notable, but not currently implemented into Towhee, is the self-adapting fixed-endpoint (SAFE) CBMC algorithm of Wick and Siepmann 2000. The problem with SAFE-CBMC is it can use a large amount of memory in order to keep track of all of the adapting fixed-endpoint biasing functions. The version implemented into Towhee uses some analytical biasing functions based upon a crude, but consistent, transformation of the distance between growth atoms and target ring atoms, into a bias function based loosely upon dihedral, bending, and vibrational energies. Considerable research is still needed in this area to determine optimal biasing strategies. The algorithm implemented into Towhee was first used, but not satisfactorily described, in Martin and Thompson 2004. A more detailed description of the biasing functions was later published in Martin and Frischknecht 2006.

@pgrinaway
Copy link
Member Author

This is resolved.

@jchodera
Copy link
Member

CBMC stuff may still be useful. Maybe split this off into separate issue?

@pgrinaway
Copy link
Member Author

Oh, I agree completely that CBMC stuff is still useful. I'm personally very interested in it and would love to see where we could take it. That said, I'd rather not have a very general theory question as an issue that could hang around forever--I'd personally prefer to keep the issues to near-term goals and software problems. I'm not wedded to this view--if you think that's a bad plan and that it'd be best to have these discussions here feel free to reopen.

@pgrinaway
Copy link
Member Author

Or open a new issue/whatever you think is best.

@jchodera
Copy link
Member

This plan is fine with me.

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Projects
None yet
Development

No branches or pull requests

2 participants