Openff-toolkit: Raise exception when generic parameters are assigned

Created on 17 Nov 2017  路  24Comments  路  Source: openforcefield/openff-toolkit

In SMIRNOFF99Frosst, generic parameters are set to fairly ridiculous values. For example, the generic bond ([*:1]~[*:2] for any two atoms connected by any bond) the length is 4.0 angstroms and the force constant is 2,000 kcal/mol/angstrom**2. For our molecule sets we have eliminated molecules where these parameters are assigned, indeed that is how we discovered what parameters were missing when developing smirnoff99Frosst to begin with. However, if these parameters are assigned to a molecule the simulation will run without any kind of warning and unless someone looks closely at the data they may never know there is a problem.

@davidlmobley and I agree that as more people begin to use our force field running quietly with these parameters is not ideal. Offline we discussed two options for how to handle this scenario:

  1. We could make the parameters slightly more reasonable and add a warning that generic parameters were assigned to the parameters.

  2. Default behavior is to have an exception raised when one of these parameters is assigned. There should be an option that allows users to run the simulation anyway where a warning would still be printed.

David and I think option 2 is the better solution, but this raises a couple of questions:

  • Currently generic parameters are assigned the number 1 in their category (i.e. b1 for Bonds and n1 for non-bonds). This numbering scheme comes from how we convert the "frocmodish" files into the smirnoff ffxml format where each parameter is numbered in order. Should we maintain this convention where we look at the parameter id? OR assume the first parameter in any section is generic?

  • Should we keep the generic parameters as is or change the numbers to be slightly more reasonable values since default behavior will be to raise an error. If we change the values what should we change them to?

@cbayly13 could we get your input on this when you have time?

enhancement question

All 24 comments

It's worth noting that the parameter ID is purely optional, so when we talk about retaining the convention where we look at the parameter ID, we're really just talking about whether WE are using it internally as a way to assess whether a generic is assigned. All our documentation so far says this is optional and has no functional importance, and I'm leery about changing that.

To explain the history here, as I remember it -- Chris really wanted generics present which were "obviously wrong" because this was something he'd used previously in the past effectively. But, as we broaden this out so that more people use it, I get increasingly uncomfortable with having these present as most people are just going to plunge ahead without noticing they are being used.

As I think about it more, @bannanc , maybe this all is an argument to extend the spec slightly to allow explicit designation of generics, e.g. adding a Generic_bond tag in addition to Bond, Generic_angle, etc. These would be explicitly designated as generics by their TAGS (not their parameter IDs, which we say have no functional importance). Then, when we go with option 2, if the user does NOT insist on proceeding with generics, they would simply get an exception about a missing parameter if a parameter is missing. But if they insist on proceeding with generics, then these would be employed.

Does that sound sensible?

As I think about it more, @bannanc , maybe this all is an argument to extend the spec slightly to allow explicit designation of generics, e.g. adding a Generic_bond tag in addition to Bond, Generic_angle, etc. These would be explicitly designated as generics by their TAGS (not their parameter IDs, which we say have no functional importance). Then, when we go with option 2, if the user does NOT insist on proceeding with generics, they would simply get an exception about a missing parameter if a parameter is missing. But if they insist on proceeding with generics, then these would be employed.

@davidlmobley This seems like a more reasonable choice to me, how much of a change to the API and forcefield code would this require?

I think it should be a minor change.

I agree with the concerns raised; we need to change the behavior. I agree that Caitlin's option 2 to raise an exception is the right way to go. I am wondering about changing the default parameters to be "slightly more reasonable": As David points out, my original reasoning is that these parameters have only one purpose and that is to tell that something really bad and unnatural is happening, never anything reasonable.

I think we should have exceptions raised, period, if the "obviously wrong" world is encountered. Don't forget that one of the use cases for the "obviously wrong" exception is a structure with broken valence, ie an illegal valence structure, which one definitely gets in large dbs like emolecules. There are no reasonable parameters for those and we want exceptions thrown in those cases. For new but not wrong structures perhaps we should put in a little bit more due diligence on _real_ generic parameters that are some kind of reasonable . I thought I put a complete set into smirnoff99Frosst, but if Caitlin's been hitting the "obviously wrong" parameters in such cases I guess I missed some (no surprise).

I like David's idea of tags; this could help us in a general way to track provenance; the lowest form of provenance is probably "generic". It could still be fit to QM though.

@cbayly13 - by "generic" I was thinking mainly parameters of the type "contains no element and will apply to anything I've never seen before". You certainly did do a good job (as far as I can tell) of covering elemental generics in a sensible way; what I'm mostly worried about now is people managing to come up with corner cases where they throw in an element we don't yet support and get a "bad" parameter assigned to it.

Caitlin could put in some more sensible non-elemental generics in these cases (e.g., not a bond length of 4 angstroms for the generic bond; maybe something like 1.5 angstroms) so that if people insist on running, it will at least run...

But it sounds like we're agreed.

@bannanc - let's talk at some point about how to do this.

@jchodera -- I wanted to get this issue on your radar for the SMIRNOFF spec updates in #86 -- we want a way to flag generic parameters that are in a forcefield to ensure coverage but are not really intended for use:

As I think about it more, @bannanc , maybe this all is an argument to extend the spec slightly to allow explicit designation of generics, e.g. adding a Generic_bond tag in addition to Bond, a Generic_angle in addition to Angle, etc. These would be explicitly designated as generics by their TAGS (not their parameter IDs, which we say have no functional importance). Then, when we go with option 2, if the user does NOT insist on proceeding with generics, they would simply get an exception about a missing parameter if a parameter is missing. But if they insist on proceeding with generics, then these would be employed.

Basically the thought is that this would allow us to build in bonds like [*:1]~[*:2] (which smirnoff99Frosst currently has, but which we don't like because it's risky) which will raise an exception UNLESS the user explicitly specifies an option like use_generics = True or something similar, in which case only a warning would be issued.

Basically the incentive is to take advantage of the fact we can have very broad generics to ensure everything can be simulated, but also clearly designate parameters which are not really intended for use.

Basically the incentive is to take advantage of the fact we can have very broad generics to ensure everything can be simulated, but also clearly designate parameters which are not really intended for use.

I'm a bit confused about this. Do we want to allow the user to simulate these or not?

We can't desire both things simultaneously, can we?

@jchodera

I'm a bit confused about this. Do we want to allow the user to simulate these or not?

We can't desire both things simultaneously, can we?

Yes, actually -- take the behavior of OpenEye's Omega as an example. There's an option for StrictAtomTypes; if True, it will refuse to generate a conformer if an exactly matching atom type is not found; if False, it will take the nearest atom type and proceed. I use both depending on which application I'm working on.

I would want people NORMALLY to not be able to use generic parameters here, but we could come up with generics which would be "good enough" for really rough simulations; e.g., to first order, treating a boron in an aromatic ring as "roughly like a carbon" will be preferable to some people over "you just can't model it until you fit the parameters", e.g. if you're doing docking or something similar where you aren't after a great deal of accuracy.

So, I think we CAN desire both things.

You can achieve this right now by simply having two versions of the same ffxml file (one including generics, one not), or, if we go with the latest proposal in https://github.com/openforcefield/openforcefield/issues/96, simply prepending an ffxml file with generics.

Does this really require extensions of the ffxml spec and the ForceField API if we can already do it with barely any effort?

We could of course add this option, but I just want to make sure that it's understood that you can essentially do this now (though we're still finishing #96).

Ah, I totally hadn't thought of that; while less user-friendly that may be a good idea since it doesn't require extending anything. @bannanc is this OK with you?

I hadn't thought about that as an option either, but I like it. We don't really want using generics to be too "user" friendly so this seems like a good solution to me.

Great! If you just create a generics.ffxml, it's just a difference between

forcefield = ForceField('smirnoff.ffxml', 'water.ffxml')

and

forcefield = ForceField('generics.ffxml', 'smirnoff.ffxml', 'water.ffxml')

presuming we have reached consensus in https://github.com/openforcefield/openforcefield/issues/96

It would take me about 5 minutes to make a generic.ffxml file, if someone can let me know where we should save it.

Right now, ForceField automatically searches the path where the data/ directory is installed, so if you put it here, you could specify

forcefield = ForceField('forcefield/generics.ffxml', 'forcefield/smirnoff.ffxml', 'forcefield/water.ffxml')

We can always change the search path here, or add more paths to search.

I realized I should document somewhere why I'm choosing the generic parameters I have chosen, I'm going to write out my thought process here and maybe I can get @cbayly13 to comment if he agrees.

My thought is we want the generics to be at least somewhat reasonable, but we also know what types of elements our force field is less likely to have information. Anything we pick is likely to be better for somethings than it is for others. Since boron is the element that keeps coming up as something we want, but don't have parameters for, I was thinking I would pick parameters that seem reasonable for those.

[*:1]     1.9    0.1

[*:1]~[*:2]    310.0    1.5

[*:1]~[*:2]~[*:3]    50.    120.

[*:1]~[*:2]~[*:3]~[*:4]    1     0.0    0.0    2

@bannanc I believe this still needs to be addressed, correct? I think this will be important to include in the 1.0 release so I was thinking of having a go at this if you're not already working on it.

Currently, in the topology branch, there are no generics, but ForceField.create_openmm_system(topology) fails silently with some unexpected behavior if no match is found for an element of the topology in the force field. A couple of examples,

  • When no LJ parameters are found for a particle, they just get assigned epsilon=0 and sigma=1.
  • When a bond is not matched, the HarmonicBondForce simply has no term for it.

In the rush to make the January workshop release, I had disabled the missing valence checks. We should re-enable those before the alpha release so exceptions can be raised if needed here.

Ah! Sorry, I had missed that. I think I'm going to reenable them, and fix the eventual test failures that come up then.

Thanks, Andrea. Sorry you're spending so much time cleaning up my messes :-/

No worries at all!! We needed to have a working prototype under time constraints so that was absolutely necessary.

OK, I was going to say that these should have been fixed by the new SMIRNOFF99Frosst that does not have generic parameters so you should get an error if you use that file if your parameters are not assigned.

A "power user" could use a generic offxml file in conjunction to get around this problem if they really wanted to run a simulation with only the generic parameters. Like John outlined above

So big picture:

  • If no parameters are assigned, there should be an exception (so those checks need re-enabling)
  • But generics are perfectly fine

I think that means THIS issue can be closed and there should be a separate issue to deal with https://github.com/openforcefield/openforcefield/issues/75#issuecomment-465211049

Thanks for your input, @bannanc, @davidlmobley .

I think that means THIS issue can be closed and there should be a separate issue

Will do!

Was this page helpful?
0 / 5 - 0 ratings