Re: [AMBER] Fwd: Paramfit issues - Improvement of R^2 value.

From: David Cerutti <dscerutti.gmail.com>
Date: Mon, 30 Apr 2018 23:05:25 -0400

If you are getting that form mdgx, there is definitely something strange
about your input coordinates. mdgx is telling you that the conformation,
as input, violates some of the bonds or other high-frequency energy terms
that the force field in the topology sil.top thinks is governing the
system. It's trying to get things right but it can't. The feature was
really written for times I wanted to take coordinates straight from someone
else's quantum input deck and run them with my own MM topology, but didn't
have any guarantees their atom order was the same as mine, so I wrote mdgx
to figure it out for me (if it was even possible, which is your case it
looks like it's not--one conformation out of a whole line of them is just
really bad). Probably some QM optimization has gotten you into a situation
where the atoms have been chemically rearranged, or the paramfit
conformation generator has made something truly hideous. As long as you
start with a good, well-behaved structure and don't ask for tremendous Krst
values, the mdgx &configs module should keep things well behaved. If those
are the default values for Krst, I probably set them that way to ensure
that whatever restraints you apply will overwhelm the typical force field
tendencies and pull you out of whatever minima the system was originally
in. They should be fine for your first pass at this. I'd suggest 512-1024
for bonds: just keep restraints guiding the sampling at least 1.5-2x
stiffer than the parameters you're trying to sample, and remember it's all
harmonic, so the restraint near its own target value won't do much to
oppose a parameter that's far from equilibrium.

There are no rules governing how many samples to use with count, but I
would suggest about 20 samples for each dihedral you want to fit, and
perhaps 400 if there is an important combination of dihedrals that vary
together and could vary in any combination (see next paragraph). For bonds
and angles I'd suggest that you get about the same number of individual
samples, 20 over a range, and make sure the range covers both conformations
where the bonds are compressed as well as those where the bonds are
stretched. The fitting module later in the workflow needs to be able to
see "energy gets better, better, no improvement, worse, worse" in order to
place a harmonic minimum in the right spot.

GridSample is probably what you want to use for sampling only 2-3 degrees
of freedom. Each dihedral should be defined by the atoms you want to move,
and you can define as many as you want. However, understand that if you do
not use combine [ listing of restraints to combine, i.e. 1 2 3 ] GridSample
is going to make a linear interpolation between each restraint at the same
rate, so you won't get all combinations of all restraints, just a line of
conformations in the N-dimensional restraint space. The combine feature
says "no, I want count/Replicas restraints sampling all over this two or
three dimensional surface" and mdgx will then do that. The other way to
ensure coverage of the entire surface or hypersurface is to use
RandomSample rather than GridSample, although you won't be guaranteed a
regular sampling interval in that case. Run mdgx -CONFIGS and it'll print
everything for the &configs namelist.

Dave




On Mon, Apr 30, 2018 at 4:20 PM, Alechania Misturini <alechaniam.gmail.com>
wrote:

> Ops, forgot go ask about *count* value. There are some rule for choosing
> how many samples to use, related with the parameters for fitting? Cause
> when I was using paramfit, I tryed different values, but maybe I didn't use
> sufficient strutures or representative ones.
>
> Best Regards
>
> 2018-04-30 17:08 GMT-03:00 Alechania Misturini <alechaniam.gmail.com>:
>
>>
>>
>>
>>
>> *Thanks to answering me.I'm having some trouble for fitting with my old
>> data (conformations generated to use in paramfit):mdgx >> Rearrangement
>> required for conformation:mdgx >> Topology sil.top / Coordinates
>> coords.cdf, Frame 3179.mdgx >> 5 strained interactions, total bonded
>> energy 601.303086mdgx >> Optimization unsuccessful after 6536
>> iterations.mdgx >> 3179 conformations' energies computed.
>> mdgx >> High-energy conformations pruned.mdgx >> Energies
>> normalized, conformations weighted.mdgx >> Accumulating and checking
>> fitting matrix.mdgx >> Error: Matrix A is rank deficient.Segmentation fault
>> (core dumped)So, I'll generate structures following the mdgx tutorial,
>> since I want to change the level of theory employed. I search about the
>> Krst value used in to set the conformations, and didn't find more
>> information about it. Its a default using Krst = 64.0 for dihedrals and
>> 256.0 for angles? And how about bonds?In case of dihedrals that are
>> described with more than one set of parameters, they are defined multiple
>> times at this first input, with GridSample? Or it happens only when we do
>> the fitting?I'm defining the range for bonds with RandomPerturb, and around
>> { 0.7 3.5 }, is that alright? (0,74 A is the shortest bond length,
>> H_2)Sorry to bothering you again.*
>>
>> 2018-04-30 11:23 GMT-03:00 David Cerutti <dscerutti.gmail.com>:
>>
>>> I'd suggest you play around with it; try brst 0.0002 at first, and just
>>> leave brstcpl alone unless you are seeing zero movement in the equilibria
>>> with each run. With a few hundred data points you should literally be able
>>> to try, try, try again about as fast as you can type new input, so you can
>>> do that or script something to hit it with a spectrum of tether strengths
>>> to keep the parameters close to their original values. Just be sure that
>>> you are printing results into a different file--if you are overwriting the
>>> original parameters you'll be always tethering to the results of the last
>>> fitting attempt and you could see weird journeys around the parameter
>>> space. Of course, the hope is that the molecular mechanics force field can
>>> indeed make sense of the quantum energy surface, and if that is true the
>>> mdgx will steer you towards the same set of parameters every time. You
>>> should try to find brst and brstcpl such that the overall error relative to
>>> QM is minimal (technically, the error will be minimal with brst = 0.0, when
>>> brstcpl should be inconsequential) but the parameters are not too different
>>> from their original values (sometimes if your data is sparse you can even
>>> optimize angles with arst = 0.0 to find negative stiffness constants--no
>>> good!). Take the strongest value of brst that you can apply which still
>>> allows perhaps +/- 100 kcal/mol-A2 changes in the bond stiffnesses​ (the
>>> original values will probably be 350-550 kcal/mol-A2) with up to 0.1A
>>> changes in the equilibrium lengths. As always, the best test of your new
>>> parameters is to us them in subsequent mdgx &configs runs to create new
>>> conformations, new QM data, and verify that the MM model is tracking the QM
>>> benchmark as it guides molecules between energy minima.
>>>
>>> Dave
>>>
>>>
>>>
>>
>>
>> --
>> Alechania Misturini
>> Mestranda em Físico-Química
>> Bacharel e Licenciada em Química
>> Grupo de Estrutura Eletrônica Molecular –〈
>> GE
>> ⎜
>> EM
>> 〉
>> *Departamento de Química – CFM*
>> *Universidade Federal de Santa Catarina – UFSC*
>>
>
>
>
> --
> Alechania Misturini
> Mestranda em Físico-Química
> Bacharel e Licenciada em Química
> Grupo de Estrutura Eletrônica Molecular –〈
> GE
> ⎜
> EM
> 〉
> *Departamento de Química – CFM*
> *Universidade Federal de Santa Catarina – UFSC*
>
_______________________________________________
AMBER mailing list
AMBER.ambermd.org
http://lists.ambermd.org/mailman/listinfo/amber
Received on Mon Apr 30 2018 - 20:30:01 PDT
Custom Search