Re: [AMBER] Using the RED Server Dev. to parameterize a model transition state

From: Kamali Sripathi <ksripath.umich.edu>
Date: Sun, 15 Dec 2013 18:43:19 -0500

Dear Francois,

Thank you very much for the great information. I'll try the course you
suggested.

I have one final question for now: in a previous message, you mentioned
that I should provide a QM output (with one imaginary
frequency) as input for RED Python. Does this mean that I should leave the
imaginary frequency calculations of the form, e.g.,
*Frequencies -- xxxx.xxxx yyyy.yyyy zzzz.zzzz* in the QM output, but
delete all the real frequencies as specified here:
http://q4md-forcefieldtools.org/Tutorial/Tutorial-1.php#6?

Thank you so much again for all your help, and have a great day,

Kamali




On Fri, Dec 13, 2013 at 10:30 AM, FyD <fyd.q4md-forcefieldtools.org> wrote:

> Dear Kamali,
>
> > I also agree that I shouldn't use these results for my final
> calculations;
> > I just used the OPT_Calc = "OFF1" option as part of a "practice run" so
> > that I could see how to successfully submit a job to R.E.D Python and how
> > the outputs would look.
> >
> > In a previous message you also asked why I used the INTRA-MCC constraints
> > that I did: this was also part of the "practice". I chose atoms to
> > constrain that made sense to me, and I just wanted to see what would
> > happen.
>
> You constrain atoms (using intra-mcc) that needs to be; for instance
> to generate a molecular fragment from a whole molecule; in your case
> you study a TS so I do not think you need to generate a molecular
> fragment (except if you do want to connect it to a larger entity).
>
> > I actually don't think I need to constrain any atoms beyond the
> > default setting; is it possible to also use defaults for INTRA-MCC?
>
> No
>
> > Maybe I should start a different thread for this question, but I have
> > actually tried to submit job P8902 with OPT_Calc = OFF and providing my
> own
> > QM output. However, I've gotten errors in the R.E.D.-
> > Server-14651.master0.q4md-forcefieldtools.org file:
>
> ok I looked at the minimum you provided in the P8902 job. This is not
> the TS you are interested in... please compare your PDB input file and
> this QM output i.e. the last set of Cartesian coordinates...
>
> > *...*
> > *Reading JOB1*.log for molecule : 1*
> > *Find "Normal termination of Gaussian"*
> > *Find "Normal termination of Gaussian"*
> > *Find "Stationary point found"*
> > *Find "Standard orientation"*
> > *Element in JOB1 is not the same with PDB*
> > *Find big error in log : /home/u0263B7BZ/P8902/JOB1-gau_m1-1.log*
> > *JOB1*.log is not correct , can not continue......*
> > *Information of end*
> > *Time : Wed Dec 11 22:26:46 CET 2013*
> >
> > But to me the JOB1-gau_m1-1.log looks fine. Did I have to modify my QM
> > output (similar to what one must do for R.E.D Server Perl)?
>
> try removing the Frequency job from the Mol_red1.log file as described
> in red color at:
> http://q4md-forcefieldtools.org/Tutorial/Tutorial-1.php#6
> (see -4-The last step consists of checking and preparing the P2N and
> QM output files for the R.E.D. III.x execution)
>
> >> In these conditions _the key_ is the definition of the atom
> >> connectivities from this PDB file (atom typing & FF parameters are
> >> based on these atom connectivities). You have C-H bonds with 1.5 Angs
> >> length in your PDB input file. Are you sure an atom connectivity is
> >> created in this case (VMD does not...)?
> >
> > You're correct; those bond lengths should be fixed, and I will do that
> > before the next time I submit my structure.
> >
> >
> >> Your best bet here (because RED Python computes atom connectivities
> >> after the geometry optimization step and you skip this step...) is to
> >> provide the atom connectivities *YOU* want in the PDB input file and
> >> *prevent* their computation.
> >>
> >> See
> >>
> http://q4md-forcefieldtools.org/REDS-Development/Demo2-Files/Project.config
> >> MOLECULE1-CALCONECT = OFF
> >>
> >> then you need to carefully think which atom types you want: obviously
> >> two CT-H1 or CT-H2 bonds (CT, H1, H2 are atom types) cannot have two
> >> different equilibrium values; i.e. ~1.0 and ~1.5 Angs -> a new atom
> >> type has to be created...
> >>
> >> Moreover your phosphor atom has 5 bonds I would set up a new atom type
> >> for it...
> >
> > Thank you for the suggestion, Francois. But I am wondering if the new
> > phosphorus atom type is necessary? I ask because it is possible in nature
> > for phosphorus to form 5 bonds.
>
> look at the atom type definition in the Amber FF for the phosphorus
> atom in the phosphate group; it is related to this NA structure:
> http://en.wikipedia.org/wiki/File:DNA_chemical_structure.svg
>
> >> Finally I really hope the Cartesian coordinates you have were checked
> >> by QM: once again better providing a QM output (with one imaginary
> >> frequency) as input for RED Python...
> >
> > For job P8902, I used Gaussian09 with these options:
> >
> > #P hf/6-31G* Opt=(Tight,CalcFC) Freq SCF(Conver=8) Test
> >
> > I tried to answer this question on my own, but I got confused: why do I
> > need imaginary frequencies for my calculations, and why only just one?
>
> See http://www.gaussian.com/g_tech/g_ur/k_opt.htm
> & use the QST2 or QST3 keyword.
>
> one imaginary frequency for a TS is what I learned: I think you need
> to read about that ;-)
>
> regards, Francois
>
>
> >> Quoting FyD <fyd.q4md-forcefieldtools.org>:
> >>
> >> > Kamali,
> >> >
> >> > please check the PDB file format you use as input...
> >> > see
> >> http://q4md-forcefieldtools.org/REDS-Development/Demo2-Files/readme.txt
> >> >
> >> > regards, Francois
> >> >
> >> >
> >> > Quoting FyD <fyd.q4md-forcefieldtools.org>:
> >> >
> >> >> Dear Kamali,
> >> >>
> >> >>> I have been trying a few things out with my transition state
> structure
> >> of
> >> >>> RNA cleavage, and I had a few follow-up questions. My questions are
> >> >>> regarding my job P8824:
> >> >>
> >> >> ok
> >> >>
> >> >> in the PDB input file:
> >> >> ATOM 15 O2P RT0 1 31.511 16.369 54.925
> >> >> O
> >> >> ATOM 16 O5' RT0 6 32.667 14.580 56.630
> >> >> O
> >> >> ATOM 17 P RT0 1 30.934 15.202 55.888
> >> >> P
> >> >> ->
> >> >> ATOM 15 O2P RT0 1 31.511 16.369 54.925
> >> >> O
> >> >> ATOM 16 O5' RT0 1 32.667 14.580 56.630
> >> >> O
> >> >> ATOM 17 P RT0 1 30.934 15.202 55.888
> >> >> P
> >> >>
> >> >> & look at the residue name generated by R.E.D. Python
> >> >>
> >> >> You deactivated many keywords in the Project.config file letting the
> >> >> defaults to be executed; this is correct. I am glad to see you can do
> >> >> the job without a tutorial. We need to provide a tutorial. The
> >> >> description of the keywords in the Project.config and
> Configuration.py
> >> >> files is too limited...
> >> >>
> >> >>> Following your suggestion in the private communication (when we
> first
> >> >>> talked about this structure), for this job I just ran the PDB
> through
> >> the
> >> >>> RED Server Dev. with the option OPT_Calc = Off1 in the
> Project.config
> >> file,
> >> >>> as a practice run.
> >> >>
> >> >> OPT_Calc = Off1 means you skip the QM job and use the Cartesian
> >> >> coordinates from the PDB file to be used in MEP computation. (I am
> >> >> really not a big fan of this approach; I would prefer you provide a
> QM
> >> >> job for a Cartesian coordinate set with one imaginary frequency). ok
> >> >> for a 'dirt cheap' job ;-) (just my personal opinion)
> >> >>
> >> >>> I was wondering if there is any documentation for the
> >> >>> files generated by RED Python?
> >> >>
> >> >> Not yet - so far we work on a new version of R.E.D. Python with new
> >> >> features...
> >> >>
> >> >> Just ask for help; sorry for that...
> >> >>
> >> >>> I found two types of files in the
> >> >>> P8824/Data-R.E.D.Server/Mol_m1 folder: *sm* files (which I know are
> >> files
> >> >>> containing intra-molecular charge constraints) and the *ia*: what
> does
> >> the
> >> >>> *ia* mean? I also have two .mol2 files: Mol-ia0_m1-c1.mol2
> >> >>
> >> >> sm means files related to single molecule charge fit (without
> intra-mcc)
> >> >> ia means files related to charge fitting carried out with intra-mcc
> >> >> (because you use 'MOLECULE1-INTRA-MCC1 = 0.0 | 1 4 11 | Keep')
> >> >> mm means files related to multiple molecule fit (i.e. number of
> >> >> molecules > 1)
> >> >>
> >> >>> and Mol-sm_m1-c1.mol2, but these files appear to have two columns of
> >> >>> charges rather than one next to the rightmost **** column. Could you
> >> tell
> >> >>> me which column is the correct one for charges, or alternatively
> point
> >> me
> >> >>> to documentation that answers my questions?
> >> >>
> >> >> See
> >> >>
> >>
> http://cluster.q4md-forcefieldtools.org/~x/P8824/Data-R.E.D.Server/Mol_m1/Mol-sm_m1-c1.mol2
> >> >>
> >> >> 1 C4' 30.706000 13.015000 52.944000 CT 1 U01 0.0520 0.8780
> >> ****
> >> >>
> >> >> 0.0520 is the column of the atomic charge value
> >> >> 0.8780 is the column of the atomic polarizability
> >> >> (useful for the non-additive FF model)
> >> >>
> >> >> ---
> >> >>
> >> >> I do not understand why you use: MOLECULE1-INTRA-MCC1 = 0.0 | 1 4 11
> |
> >> Keep
> >> >> in
> >> >>
> >>
> http://cluster.q4md-forcefieldtools.org/~x/P8824/Data-R.E.D.Server/Mol_m1/Mol-ia0_m1-c1.mol2
> >> >>
> >> >> i.e. 0.6957-0.5092-0.1865=0.0000
> >> >>
> >> >> then look at the impact of the INTRA-MCC1:
> >> >>
> >>
> http://cluster.q4md-forcefieldtools.org/~x/Project/P8824/Data-R.E.D.Server/Mol_m1/Statistics_m1.txt
> >> >>
> >> >> See the table:
> >> >> Charge values for single molecule fit; without (SM) vs with (IA)
> >> INTRA-MCC:
> >> >> -> the impact seems quite important...
> >> >>
> >> >> You could also look at :
> >> >> - the RRMS values of the corresponding fits.
> >> >> - the es*.pdb files to display where the errors are.
> >> >> See http://q4md-forcefieldtools.org/RED/resp/#intro
> >> >> espdb
> >> >> -j PDB-like file with relative residual in the TempFactor field
> >> >> esqpotpdb
> >> >> -y PDB-like file with input MEP values in the TempFactor field
> >> >> (atomic units)
> >> >> esmpotpdb
> >> >> -z PDB-like file with MEP values for new charges in TempFactor field
> >> >> (atomic units)
> >> >>
> >> >> let us know if you need more information...
> >> >>
> >> >> regards, Francois
> >> >>
> >> >>
> >> >>> On Tue, Dec 10, 2013 at 8:12 AM, Kamali Sripathi
> >> >>> <ksripath.umich.edu> wrote:
> >> >>>
> >> >>>> Dear Francois,
> >> >>>>
> >> >>>> Thank you very much for the detailed information, and I'm sorry for
> >> the
> >> >>>> confusion. I would have mentioned that I had questions about a
> topic
> >> >>>> related to our private discussion, but I didn't want you to have
> >> >>>> to look up
> >> >>>> that communication, and so I tried to start a new thread :)
> >> >>>> Looking back, I
> >> >>>> see that my initial email in this thread was unfortunately less
> clear
> >> than
> >> >>>> I intended.
> >> >>>>
> >> >>>> I will try the protocol you suggested, but I think it would be
> better
> >> for
> >> >>>> me to use the RESP-A1 charge type because I'm working with RNA
> rather
> >> than
> >> >>>> carbohydrates, and using AMBER10 for my simulations.
> >> >>>>
> >> >>>> Thank you very much again, and have a great day,
> >> >>>>
> >> >>>> Kamali
> >> >>>>
> >> >>>>
> >> >>>> On Tue, Dec 10, 2013 at 2:52 AM, FyD <fyd.q4md-forcefieldtools.org
> >
> >> wrote:
> >> >>>>
> >> >>>>> Dear Kamali,
> >> >>>>>
> >> >>>>> > No, thank you very much for the offer. I have the PDF and refer
> to
> >> it
> >> >>>>> often
> >> >>>>> > because there are concepts that I'm still learning :)
> >> >>>>>
> >> >>>>> Yes a looot to learn in this work...
> >> >>>>>
> >> >>>>> > I have two final question: from examples on the RED Server
> (perhaps
> >> >>>>> this is
> >> >>>>> > not applicable to RED Perl or RED Python), it appears that one
> can
> >> also
> >> >>>>> set
> >> >>>>> > the INTRA-MCC constraints to defined values, e.g., sugar
> carbons to
> >> >>>>> > standard values in the AMBER force field instead of setting the
> >> sum of
> >> >>>>> > certain atoms to 0 or to another value. Do you suggest one way
> >> over the
> >> >>>>> > other?
> >> >>>>>
> >> >>>>> I would constrain hydrogen atoms (only these in the methylene &
> >> methyl
> >> >>>>> groups, using the 'F' flag in the intra-mcc to keep them in the FF
> >> >>>>> lib) _only_ if you use the GLYCAM force field. So the 'RESP-C2'
> >> charge
> >> >>>>> model (to be provided in the Configuration.py file to overwrite
> the
> >> >>>>> default RESP-A1):
> >> >>>>> see
> >> http://q4md-forcefieldtools.org/REDS-Development/popup/popkeyword.php
> >> >>>>>
> >> >>>>> Here once again I would use INTRA-MCC1; in general one always uses
> >> >>>>> INTRA-MCC1 except when one realizes that some atoms are not
> >> >>>>> equivalenced as they should be ;-)
> >> >>>>>
> >> >>>>> > I understand that it very much depends on one's case, but I was
> >> >>>>> > wondering if you had any thoughts. I was thinking that I might
> >> >>>>> constrain
> >> >>>>> > the charges of some of the sugar atoms that are 5' of the
> >> >>>>> transition-state
> >> >>>>> > pentacovalent phosphate to their values in the AMBER force
> fields
> >> >>>>> instead
> >> >>>>> > of summing their charges to 0.
> >> >>>>>
> >> >>>>> oh oh I got it: You continue the discussion from the private
> >> >>>>> assistance...
> >> >>>>> Not easy to follow you. In your case I would:
> >> >>>>> -1 characterize the TS using the QM program of your choice
> >> >>>>> -2 run RED Python using CHR_TYP = "RESP-C2", OPT_Calc = "Off" &
> >> >>>>> MEPCHR_Calc = "On" using the QM log file & PDB file as inputs
> >> >>>>> -3 study what you get: here you will be able to:
> >> >>>>> - check the atom connectivities generated in the mol3 file
> >> >>>>> - re-run RED Python using Re_Fit = "On" (i.e. upload/provide
> the
> >> >>>>> entire/previous RED Python job in the archive)
> >> >>>>> - provide new FF parameters in the frcmod.user (see the
> >> >>>>> frcmod.unknown file generated in the previous job) in the archive
> >> file
> >> >>>>> See
> >> >>>>>
> >> http://q4md-forcefieldtools.org/REDS-Development/RED-Server-demo1.php
> >> >>>>>
> >> >>>>>
> >> http://q4md-forcefieldtools.org/REDS-Development/popup/poptestcase.php
> >> >>>>>
> >> >>>>>
> >>
> http://q4md-forcefieldtools.org/REDS-Development/Demo2-Files/frcmod.user
> >> >>>>> - provide new FF atom types in the Project.config file for key
> >> atom
> >> >>>>> center(s)
> >> >>>>> MOLECULE1-ATMTYPE = ...
> >> >>>>> (here you can load your own FF as input to R.E.D. Python)
> >> >>>>>
> >> >>>>> When ready just request a private assistance and provide the PXXXX
> >> >>>>> R.E.D. Server job name in the body of your email.
> >> >>>>>
> >> >>>>> > I lastly wanted to double-check that the INTER-MCC and INTER-MEQ
> >> >>>>> keywords
> >> >>>>> > are only necessary if one is deriving charges for more than one
> >> >>>>> molecule,
> >> >>>>> > and so would not be applicable in my one-omolecule
> one-conformation
> >> >>>>> case?
> >> >>>>>
> >> >>>>> Yes INTRA means within a molecule; INTER means between molecules.
> >> >>>>>
> >> >>>>> regards, Francois
> >> >>>>>
> >> >>>>>
> >> >>>>> >> regards, Francois
> >> >>>>> >>
> >> >>>>> >>
> >> >>>>> >> > On Sat, Dec 7, 2013 at 4:30 PM, FyD <
> >> fyd.q4md-forcefieldtools.org>
> >> >>>>> >> wrote:
> >> >>>>> >> >
> >> >>>>> >> >> Dear Kamali,
> >> >>>>> >> >>
> >> >>>>> >> >> > I am using the RED Server Development to parameterize a
> model
> >> >>>>> compound
> >> >>>>> >> >> for
> >> >>>>> >> >> > the transition state and different protonation states of
> one
> >> or
> >> >>>>> both
> >> >>>>> >> of
> >> >>>>> >> >> the
> >> >>>>> >> >> > nonbridging oxygens in ribozyme cleavage.
> >> >>>>> >> >>
> >> >>>>> >> >> Classical geometry optimization as implemented in R.E.D.
> >> >>>>> will likely
> >> >>>>> >> >> generate a true minimum; not a TS.
> >> >>>>> >> >>
> >> >>>>> >> >> So better first characterizing the TS & loading the
> >> >>>>> corresponding QM
> >> >>>>> >> >> log file as input by setting OPT_Calc = "Off" & MEPCHR_Calc
> =
> >> "On"
> >> >>>>> in
> >> >>>>> >> >> the Configuration.py file
> >> >>>>> >> >>
> >> >>>>> >> >> > My questions are regarding a
> >> >>>>> >> >> > portion of the Project.config file that one may optionally
> >> >>>>> include to
> >> >>>>> >> >> > override default options (
> >> >>>>> >> >> >
> >> >>>>> >> >>
> >> >>>>> >>
> >> >>>>>
> >>
> http://q4md-forcefieldtools.org/REDS-Development/Demo2-Files/Project.config
> >> >>>>> >> >> ).
> >> >>>>> >> >> > There is a portion of the Project.config file that deals
> with
> >> >>>>> >> >> > intramolecular constraints (INTRA-MCC). My questions are:
> >> >>>>> >> >>
> >> >>>>> >> >> INTRA-MCC in R.E.D. Perl = INTRA-MCC1 in R.E.D. Python,
> >> >>>>> which is the
> >> >>>>> >> >> classical way to set up intra-molecular charge constraint...
> >> >>>>> >> >> See Cieplak et al. in
> >> >>>>> >> >>
> >> >>>>>
> >> http://www3.interscience.wiley.com/cgi-bin/abstract/109583237/ABSTRACT
> >> >>>>> >> >>
> >> >>>>> >> >> > - How does one decide what value to assign for
> >> intramolecular
> >> >>>>> >> >> > constraints?
> >> >>>>> >> >>
> >> >>>>> >> >> In general when one wants to remove a part of a whole
> molecule
> >> to
> >> >>>>> >> >> generate a molecular fragment the intra-mcc total charge
> >> >>>>> value takes
> >> >>>>> >> >> an integer value (& this integer value is very often zero).
> >> >>>>> >> >>
> >> >>>>> >> >> See also a discussion about the relation between the
> chemical
> >> group
> >> >>>>> >> >> involved in the constraint & the total charge of this
> >> constraint:
> >> >>>>> >> >> http://www.ncbi.nlm.nih.gov/pmc/articles/PMC2918240/
> >> >>>>> >> >>
> >> >>>>> >> >> > - How does one decide whether or not to keep the
> >> >>>>> intra-molecular
> >> >>>>> >> >> charge
> >> >>>>> >> >> > constraints through the second stage of fitting (i.e.,
> >> >>>>> INTRA-MCC1
> >> >>>>> >> vs.
> >> >>>>> >> >> > INTRA-MCC2)?
> >> >>>>> >> >>
> >> >>>>> >> >> In general you want to use the way defined by Cieplak et al.
> >> i.e.
> >> >>>>> >> >> INTRA-MCC1
> >> >>>>> >> >> See
> >> >>>>> >>
> >> http://www3.interscience.wiley.com/cgi-bin/abstract/109583237/ABSTRACT
> >> >>>>> >> >>
> >> >>>>> >> >> run both types of intra-mcc on the CH3CO group in a
> dipeptide
> >> model
> >> >>>>> to
> >> >>>>> >> >> study the differences; in short an intra-mcc1 is only
> >> >>>>> applied in the
> >> >>>>> >> >> 1st resp input, while intra-mcc2 is applied in both resp
> >> >>>>> inputs; the
> >> >>>>> >> >> later allows charge equivalencing for CH2 and/or CH3 groups
> >> when
> >> >>>>> >> >> involved in the constraint...
> >> >>>>> >> >>
> >> >>>>> >> >> > If my questions are trivial ones, please direct me to the
> >> >>>>> appropriate
> >> >>>>> >> >> > references.
> >> >>>>> >> >>
> >> >>>>> >> >> no trivial question - do not worry ;-)
> >> >>>>> >> >>
> >> >>>>> >> >> R.E.D. Python is not yet published - still coding new
> >> features...
> >> >>>>> >> >> R.E.D. Perl is published:
> >> >>>>> >> >> http://www.ncbi.nlm.nih.gov/pmc/articles/PMC2918240/
> >> >>>>> >> >>
> >> >>>>> >> >> regards, Francois
>
>
>
>
> _______________________________________________
> AMBER mailing list
> AMBER.ambermd.org
> http://lists.ambermd.org/mailman/listinfo/amber
>
_______________________________________________
AMBER mailing list
AMBER.ambermd.org
http://lists.ambermd.org/mailman/listinfo/amber
Received on Sun Dec 15 2013 - 16:00:02 PST
Custom Search