Previous Thread
Next Thread
Print Thread
Page 2 of 6 1 2 3 4 5 6
Re: generating parameters for pyruvate
Kenno #24080 04/22/10 03:45 PM
Joined: Jan 2010
Posts: 85
R
RD2015 Offline OP
Forum Member
OP Offline
Forum Member
R
Joined: Jan 2010
Posts: 85
It looks like the water failed to find a stable distance from H31. After each optimization cycle it kept flying further and further away.

Perhaps roh should have had a smaller initial value?

Edited:

I had something very strange occur. I played around with roh and set it to 2, 1.5, 1 and had no convergence. Then I set it to .8 and it converged to an optimized distance of 3.6287.

To check the value then I put 3.6287 value back into the initial script and it did not converge. Was the convergence I got for .8 some sort of artifact?

I also had another thought, maybe I should make the step size smaller and increase the number of steps?

Last edited by RD2015; 04/22/10 04:27 PM.
Re: generating parameters for pyruvate
RD2015 #24082 04/22/10 04:41 PM
Joined: Dec 2005
Posts: 1,535
Forum Member
Offline
Forum Member
Joined: Dec 2005
Posts: 1,535
Pyruvate is negatively charged and as we're only optimizing the distance, the water is always pointing with its oxygen atom towards the negative charge. This yields an unfavorable charge-dipole interaction. The interaction is repulsive, and as this is an entirely physical phenomenon, there's nothing to be done about it.

How do we use this information to optimize our charges? Either we just use the qualitative information ("our force field should also have a repulsive interaction here"), or we fix the water molecule at a certain distance both in the QM and the MM and we try to reproduce the (positive) interaction energy. The latter approach gives us more quantitative information, but has the important disadvantage that we introduce an arbitrary distance into the procedure. Furthermore, equilibrium distances in gas phase and in high-dielectric condensed phase are known to differ; see Figure 8 and the accompanying discussion in the CGenFF paper. This is not well-studied for repulsive interactions (it's not a trivial thing to do) so we have no idea whether we should apply and offset, and if yes, by how much.

In the end, here's what I'd recommend you to do:
  • Calculate all the QM water interaction energies
  • Check whether the two other interactions with the methyl group are repulsive as well. The optimization often seems to converge because the repulsive force becomes smaller than the convergence criterion at a certain distance, but you'll still get a positive interaction energy.
    • If all 3 interactions with the methyl group are repulsive, we might need to fix the distance. But talk to me first if this is the case.
    • If at least one of them is attractive (there are favorable dipole-dipole interactions that may overcome the unfavorable charge-dipole interaction for some interaction geometries), we have enough information to carry on with our charge optimization without having to bother about quantifying the repulsive interactions.

Re: generating parameters for pyruvate
Kenno #24084 04/22/10 06:34 PM
Joined: Jan 2010
Posts: 85
R
RD2015 Offline OP
Forum Member
OP Offline
Forum Member
R
Joined: Jan 2010
Posts: 85
I think all three might be repulsive. The H31 and H33 script crash outright. The H32 script converges, but I calculated a SCF Done: E(RHF) = -415.972756709 for "far away" (roh=50) and I calculated a SCF Done: E(RHF) = -415.968272866 for the converged value (roh = 3.20987). So I think it may be repulsive as well.

I'll continue to run the other scripts in the meantime.

Edited:

By the way, I am not sure I need to model pyruvate extremely accurately for my experiments. Am I at the point where I have a "good enough" set of parameters so that I can start minimizing and heating my protein/ligand complex (at least for now)?

Last edited by RD2015; 04/22/10 07:00 PM.
Re: generating parameters for pyruvate
RD2015 #24091 04/22/10 07:42 PM
Joined: Dec 2005
Posts: 1,535
Forum Member
Offline
Forum Member
Joined: Dec 2005
Posts: 1,535
Then we'll just forget about the methyl group. On second thought, you have only 4 degrees of freedom in your charge optimization (9 atoms - 3 aliphatic hydrogens - 1 because the carboxylate oxygens are equivalent - 1 because of the total charge constraint), so as long as you have 4 or more acceptable (attractive) water interactions, your system won't be underdefined.

Next step is to calculate the interaction energies as (energy of the complex) - (energy of the model compound) - (energy of tip3 water = -76.010530 Hartree) and then convert to kcal/mol. When you're done with that, I'll show you how to do the CHARMM water interactions.

Originally Posted By: RD2015
By the way, I am not sure I need to model pyruvate extremely accurately for my experiments. Am I at the point where I have a "good enough" set of parameters so that I can start minimizing and heating my protein/ligand complex (at least for now)?
What do you mean? You don't have any charges yet! Also, minimizing and heating usually only takes a couple of hours.

Re: generating parameters for pyruvate
Kenno #24106 04/22/10 09:15 PM
Joined: Jan 2010
Posts: 85
R
RD2015 Offline OP
Forum Member
OP Offline
Forum Member
R
Joined: Jan 2010
Posts: 85
I put the values in the .ods file inside the archive.

Attached Files
archive.tar.gz (623.94 KB, 647 downloads)
Re: generating parameters for pyruvate
RD2015 #24107 04/22/10 11:15 PM
Joined: Dec 2005
Posts: 1,535
Forum Member
Offline
Forum Member
Joined: Dec 2005
Posts: 1,535
That looks pretty good. Here's the decisions I made:
  • Between O2_0 and O2_180, I chose O2_180 because O2_0 is somewhat hindered by the electrostatic effects of the carboxylate group (as apparent from the interaction distance and energy) and because O2_180 is closer to the methyl group, a region where we have limited coverage.
  • Between O12_0 and O12_180, I chose O12_0 even though O12_180 has a more favorable interaction energy. The reasons for doing this are:
    • We already have a water molecule at the carboxylate side: O11_180 . I don't expect O12_180 to contribute a lot more information.
    • O12_0 happens to be in the region between the carboxylate carbonyl and the ketone carbonyl; it might be interesting to have some coverage there and the distance from the ketone carbonyl is still long enough to not have too much contamination.

I prepared a CHARMM input file analogous to the water_constr.inp in the tutorial. You can use this to optimize your charges.

Tip: in cases where it's impossible to get two water interactions of the same atom right, the strongest water interaction has priority. Unless the difference between the two QM interaction energies is small (we sometimes use "smaller than 2 kcal/mol" as a rule of thumb), then it's better to go for a compromise between the two water interactions.

Edit: Now may be a good time to start your MP2 dihedral potential energy scan. See the tutorial for examples; don't "fix" any coordinates.

Attached Files
pyru_water_constr.inp.txt (22.24 KB, 702 downloads)
pyru_mp2.crd (1.34 KB, 666 downloads)
Last edited by Kenno; 04/23/10 04:06 PM. Reason: Oops... forgot attachments
Re: generating parameters for pyruvate
Kenno #24110 04/23/10 05:52 PM
Joined: Jan 2010
Posts: 85
R
RD2015 Offline OP
Forum Member
OP Offline
Forum Member
R
Joined: Jan 2010
Posts: 85
I've attached the results of the run.

Is this the right idea?

i) The energy difference in 1 and 2 is less than 2kcal so I can keep the charge on O2?
ii) Increase the charge on O11 to get ene diff below 2kcal.
iii) keep the same charge on the methyl hydrogens?

Attached Files
archive.tar.gz (22.52 KB, 601 downloads)
Last edited by RD2015; 04/23/10 06:14 PM.
Re: generating parameters for pyruvate
RD2015 #24113 04/23/10 07:53 PM
Joined: Dec 2005
Posts: 1,535
Forum Member
Offline
Forum Member
Joined: Dec 2005
Posts: 1,535
You have to do some things by yourself. This is largely a matter of experimentation: adjust your charges a bit, see how the water interactions are affected, adjust them again, repeat ad nauseum... As the CHARMM run takes mere seconds, you can do several of these cycles per minute if you organize yourself a bit.

Deviation of 2 kcal/mol are not acceptable as such. In general, we aim to get all the interaction energies within 0.2 kcal/mol, as mentioned in the CGenFF paper and demonstrated in the tutorial. But there are a few exceptions. If you have two water interactions for the same HB acceptor, then in some cases it is impossible to get them all right, because any adjustment of the charges will affect both of them in the same direction. If that's the case, then one should either go for the strongest interaction, or make a compromise between the two interactions. To decide whether to make a compromise or go for the strongest, you look at the difference between the two QM interaction energies, and see if it's lower or higer than 2 kcal/mol; see my previous post. But there are cases in which you can get both of them right - I'd recommend you play with the balance of charge between the 3 carbons and see what the effect is.

Originally Posted By: RD2015
iii) keep the same charge on the methyl hydrogens?
Yes, they should be kept at 0.09. And the two carboxylate oxygens should always be given the same charge because they're chemically equivalent. And the total charge of the molecule should be kept at -1 at all time. This can be most easily accomplished by shifting charge around. For instance, you could try to subtract 0.02 charge units from each of the carboxylate oxygens (making them more negative) and move the resulting 2x0.02=0.04 charge units onto the carboxylate carbon (making it more positive).

As this leaves you only 4 degrees of freedom, this is a really easy assignment to begin with.

The initial guess water interactions look not too bad to begin with; in fact, the ones on O2 and O12 look pretty good. But do play around with the charges a bit, you might get them better. Keep an eye on the dipole moments while doing so; we want to overestimate them a bit (although they're less important than the interaction energies, especially for charged systems). And try to keep the charges physical; it's sometimes tempting to run off in an unphysical direction. If you're not sure about what's physical, try to not deviate more than 0.1 or 0.2 charge units from your Merz-Kollman initial guess.

Last edited by Kenno; 04/23/10 07:58 PM. Reason: Typos...
Re: generating parameters for pyruvate
Kenno #24167 05/03/10 08:04 PM
Joined: Dec 2005
Posts: 1,535
Forum Member
Offline
Forum Member
Joined: Dec 2005
Posts: 1,535
It's been a while... any luck with either the charge optimization or the dihedral potential energy scan?

If you can't find a way to improve the charges, well, in about 2% of cases, the Merz-Kollman charges are as good as it gets, so you might just be "lucky" wink

Edit: however, I would recommend rounding the charges to two digits after the comma. This is to preserve compatibility with some softwares that do not have a higher precision.

Last edited by Kenno; 05/03/10 08:08 PM. Reason: added note about rounding.
Re: generating parameters for pyruvate
Kenno #24259 05/11/10 02:56 PM
Joined: Jan 2010
Posts: 85
R
RD2015 Offline OP
Forum Member
OP Offline
Forum Member
R
Joined: Jan 2010
Posts: 85
Hi. Sorry I've been swamped with exams and a presentation (I just finished today). I'll post my results soon.

Last edited by RD2015; 05/11/10 02:57 PM.
Page 2 of 6 1 2 3 4 5 6

Moderated by  alex, chmgr, John Legato 

Link Copied to Clipboard
Powered by UBB.threads™ PHP Forum Software 7.7.4
(Release build 20200307)
Responsive Width:

PHP: 5.6.33-0+deb8u1 Page Time: 0.011s Queries: 37 (0.005s) Memory: 1.0005 MB (Peak: 1.1368 MB) Data Comp: Off Server Time: 2020-05-25 18:26:05 UTC
Valid HTML 5 and Valid CSS