/
2021-06-16 Meeting notes

2021-06-16 Meeting notes

Date

Jun 16, 2021

Participants

  • @Jeffrey Wagner

  • @Connor Davel

  • @Owen Madin

Discussion topics

Item

Notes

Item

Notes

Detecting proton migration

  • CD โ€“ Worked on other things

  • CD + OM โ€“ We should prioritize making a database of this

Non-geomeTRIC fixes:



  • CD found some sqm keywords that may help -- maxcyc and SHAKE, and one other

    • CD โ€“ Iโ€™ve tested this out.

GeomeTRIC background and starting points

  • How to start a restriction scheme?

Database of failing cases/molecules

  • Proton vs. heavy atom migrations

  • How to make a database of these?

    • Initial test case

    • Zwitterions from protein structures

    • Non-protein cases

    • Non-zwitterion case

  • CD โ€“ Diversity in functional groups?

  • JW โ€“ Possibly add N1NC(=N)C(=N)1)

  • CD โ€“ Which molecule datasets do we have available?

    • Minidrugbank.sdf openff-toolkit/openff/toolkit/data/molecules at main ยท openforcefield/openff-toolkit

      • Are perception differences here a major issue?

        • JW โ€“ No, in most cases both ways to interpret an ambiguous molecule are valid. So weโ€™ll just want to make sure to regenerate the conformers to place hydrogens consistently/accoding to the graph

        • CD โ€“ I already do generate new conformers. So Iโ€™d load an SDF, generate new conformers (disregarding the roiginal), output those conformers to NEW sdfs, and then run those with antechamber.

        • JW โ€“ Perfect.

    • LOW VALUE, POSSIBLY HARD: DrugBank_X.mol2, AlhEthOH files, and FreeSolv

      • Antechamber can read mol2 directly. But there are three โ€œsub specificationsโ€ of mol2. So you can try to chunk out the mol2 files from the link above and see if antechamber can read any of them.

    • burn-in-sdf

  • ย 

  • ย 

ย 

Whatโ€™s the โ€œgold standardโ€?

  • Is OpenEye AM1 always the goal to match? Or should we look at something like RESP, or physical property accuracy?

Reproducing Case

Two conformers of the following molecule were tested as a proof-of-concept for an AM1-restriction scheme:

Conformer 0 (no proton rearrangement)

Conformer 1 (proton rearrangement from the highlighted atom to the oxygen)

Conformer 0 (no proton rearrangement)

Conformer 1 (proton rearrangement from the highlighted atom to the oxygen)

ย 

Before AM1 minimization

ย 

Before AM1 minimization
After AM1 minimization

ย 

After AM1 minimization

ย 

Conformer 0 does not undergo any proton rearrangement whereas conformer 1 sees proton migration from a nitrogen to an oxygen. Since we are concerned with the affect of this proton migration on partial charge assignment, the resulting partial charges were plotted using both amber and openeye. The black dots represent openeyeโ€™s charge assignments, which where very consistent between conformer 0 and 1 and show up as a single point. In the amber case, notice differences between the charge assignment of conformer 0 and 1.

ย 

For reference, atom 15 is the highlighted nitrogen, atoms 28 and 29 are the carboxylate oxygens, and atoms 16 and 17 are the nitrogen protons (one of which migrates). Rearrangement of one of the conformers results in a large charge discrepancy between the two conformers. This is most apparent on the nitrogen and oxygens, the atoms to/from which the proton migrates.

Solutions

The following section quickly explores 3 solutions that could be easily implemented without much change in how sqm is run/called. Each solution will be compared for internal consistency (do both conformers give the ~same partial charges?) and for consistency with other methods (in this first case, openeye). While interesting and good to know, these initial solution were ultimately ineffective at producing โ€œgoodโ€ partial charges. The solutions are:

  1. Restricting any and all geometry changes with amberโ€™s โ€œ-ek maxcyc=0โ€ flag, which restricts the number of iterations in which the geometry changes to 0.

  2. Restricting atom-proton bonds with amberโ€™s โ€œqmshake=1โ€ flag, which holds all atom-proton bonds fixed and so does not allow any proton migrations to occur.

  3. First minimizing the geometry with MD (forcefield openff-1.3.0.offxml) and then running AM1 with โ€œ-ek maxcyc=0โ€. This intends to first find a global energy minimum before running AM1, hopefully eliminating the problems with method 1.

To help visualize, the resulting molecules of each simulation are tabulated below. Only Conformer 1 (the second conformer) is shown.

ย 

Conformer 1 before simulation

Conformer 1 after simulation

ย 

Conformer 1 before simulation

Conformer 1 after simulation

Method 1

โ€œmaxcyc=0โ€

(no geometry changes allowed)

ย 

no change in geometry

ย 

Method 2

โ€œqmshake=1โ€

ย 

folding of the molecule

ย 

Method 3

MD energy minimization followed by AM1 partial charge calculation with โ€œmaxcyc=0โ€

ย 

very small change in geometry

ย 

Each of these solutions โ€œsolvedโ€ the problem in their own way. That is to say, protons 16 and 17 are still attached to the nitrogen. Two things to note are: the folding of the molecule in Method 2 (it seems that the protons are still somehow very strong attracted to the carboxylate oxygens during simulation!) and the very slight (almost unnoticeable) change in molecule geometry in method 3. It is clear the MD simulations cannot very accurately mimic the geometry changes that occur in qm simulations.

In a similar way as before, the partial charges were compared to openeye charges.

Method 1 (maxcyc=0)

Note how the partial charge assignments of amber are very consistent between conf0 and conf1 (the red/magenta points are nearly on top of each other); however, the difference between openeye and ambertools results is concerning (especially in atoms 15, 16, and 17)

Method 2 (qmshake=1)

These results are similar to method 1. The amber points indicate that the charge assignments in ambertools are now fairly consistent between conf0 and conf1; however, these assigned charges are far from the assigned charges of openeye in atoms 15, 16, and 17). It is also clear that the folding of the molecule in conf1 introduced some subtle differences in the assigned charges compared to conf0 (that is, the red/magenta points are slightly displaced from eachother)

ย 

Method 3 (MD geometry minimization followed by AM1 with maxcyc=0)

The results of Method 3 are most similar to those of method 2. First minimizing the geometry with openff MD did not appear to bring ambertools AM1 results into closer alignment with openeye.

Two things are clear at this point:

  1. A new metric of success is needed. Openeye is a great first metric to test our results against. Still, in the future, it may be useful to use other qm simulations to compare the success of our methods. Another option could be physical property tests.

  2. A restriction scheme is the next step. It would be useful if we could apply restraints to the starting internal geometry of the molecule, or at least of certain parts of the molecule. This is likely where this project is heading.

Action items

Decisions