CHARMM Development Project Forums User Discussion & Questions Parameter Set Discussion H5 Hydrogen on 1,2,3-triazole in CGENFF

Hi,

Using the General CHARMM parameters (CGENFF2B6) I haved cobbled up a molecule which contains a 1,2,3-trizole. I used the topology/parameters for 1,2,3-triazole (TRZ3) in CGENFF but there is a hydrogen on the N5 nitrogen which I can't find in other published work on this triazole. Can someone please guide me as to why the hydrogen is there in the MM parameters and absent from the previous QM work?

Thank you

You must be referring to the tautomeric equilibrium. If so, I'd like to point out that CGenFF is parameterized towards MD simulations of drug-like molecules in water. We don't care that gas-phase experiments and calculations show a strong preference for 2H protonation; as far as I know, solution-phase experiments (as well as some simulations that take into account the solvent) indicate that both tautomers coexist in significant quantities in water. Also relevant is the fact that these heterocycles are meant as scaffolds in larger drug-like molecules. Substitutions of any of the nitrogens are common in drug-like molecules, and substitutions of the carbons but not the nitrogens often make 1H protonation more favorable.

That said, thank you for bringing it to my attention that there's only one out of of two significant tautomers in the topology file. The person who parameterized 1,2,3-triazole probably just picked the tautomer that was depicted in the maybridge heterocyclic ring numbering chart; I will add the other one to our to-do list.

As a side note, I find using the CGenFF program (available through the ParamChem interface) a great deal easier and less error-prone than "cobbling up a molecule" by hand. That's right, we're eating our own dog food. See the sticky at the top of the "Parameter Set Discussion" forum. Yum! It probably won't even do a disastrous job on the 2H tautomer.

Kenno,

Thanks a lot for the quick and detailed response. I have already run some MDs on the system and was a bit worried whether I need to redo it or not. Now with your explanation, I think I know what I did wrong. By the way, I'm not a chemist.

As you correctly pointed out, it's not a single triazole. It's part of a larger molecule resulting from a click reaction of an azide with an alkyne. So using the numbering scheme you point to from Maybridge, N1 would be the substituted nitrogen, that simple fact helps me a lot. I think I made a mistake by taking N3 as the substituted Nitrogen which leads to a triazole ring with a hydrogen placed on N5 (3H in the 1,2,3-triazole tautomerization language?). I think I need to go back and redo the charge assignments and (sigh) the dihedral parameter assignments by analogy. The latter was real torture when I did it first.

As with ParamChem, thank you for the hint. If it can do what it says, it will be of IMMENSE help to me. I did a quick test supplying it with my already working PDB and it gave me the following error :
===
Processing uploaded File resulted in a few errors and/or warnings: attype warning: amide base not supported;
skipped molecule.
===
I schemed the faqs but didn't find anything related. I was supplying a PDB and NOT mol2. I'll try that later and see how it works. This is a dipeptide with a PRO at c-term (it's a skethcy in-progress PDB, needs tweeking itself, but just a yes/no test). The PRO is really an azidoproline clicked to an alkyne with a single benzene ring side chain.

Thanks again for your time.

To make the discussion extra confusing, the numbering in the topology file is incorrect. The maybridge chart has the correct IUPAC numbering, which I used in my previous post and will continue using throughout this discussion.

If N1 is the substituted nitrogen, and there are no other substitutions on the ring, then you should not have any problem; as long as the carbons are not substituted, N1 and N3 are chemically equivalent.

And that would be the reason why the CGenFF program rejected it. Since you asked, I took the liberty of taking a peek at your molecule on the server, and the protonation state is invalid. Specifically, your molecule contains another nitrogen (not in the 1,2,3-triazole ring) that has two single bonds to carbon atoms and no other bonds. Therefore, the interface assumes you're trying to work on a deprotonated negatively charged nitrogen, an "amide base" in organic chemistry speak. As these "superbases" cannot exist in watery solution, they are (as yet) not supported by CGenFF. The error message is spot on!

More generally spoken, it is extremely important that you submit the molecule to the CGenFF program in the exact protonation state you're expecting in your simulation conditions (this is usually either at physiological pH or in interaction with an active site that may influence the protonation state). At physiological pH, I would expect your 1,2,3-triazole ring to be unprotonated (the 1-substitution already counts for a proton in this context), and the secondary amine group would be protonated and positive. In other words, starting from the structure you submitted, you'd have to remove the proton from your 1,2,3-triazole ring and add two protons to your secondary amine.

Edit: this is not an easy molecule. Even a chemist with CGenFF experience would have to be very careful in order to assign correct atom types and initial guess parameters. Usage of the CGenFF program is highly recommended in this situation. It usually does a pretty good job. That said, we recently identified an issue in the current beta version of the program that may result in bad parameters with low penalties involving the substituted carbon in your 1,2,3-triazole ring, so take the parameters that involve this carbon with an extra grain of sand. We hope to be able to fix this issue in a few weeks.

Yeah I started by using the Maybridge numbering and then in the end mixed up with the param numbering (N5)

The original molecule I pathced up has a ferrocenyl ring system in place of the phenyl (talk about easy molecule ). In getting the uploaded system (just a phenyl), I just did a dirty replacement of the cyclopentadienyl ring with the phenyl ring and used the charge assigments of the biphenyl (BFL) as a first guess.

Thank you for going that far to check the actual file. I don't know what I was thinking when I uploaded it since when I look at it right now, it looks ridiculous : As I said it should have a PRO but for whatever reason the PDB was truncated. You are right in all your guesses about the protonation sates of those nitrogens. I fixed the PDB but to fix the main problem (the hydrogen) I want to first fix the charges in my topology file and r-assign the atom types (based on CGENFF) and then remove the hydrogen. After that I will submit it to the server. Right now it gives me the following error :
===
readmol2 warning: non-unique atoms were renamed.
attype warning: carbon radical, carbocation or carbanion not supported;
skipped molecule.
===
I don't know what that means but I will retry when I fix the main topology and regenerate the PDB. Are the error codes documented anywhere I can look at?

I'll update the thread when I fix the hydrogens and resubmit the molecule.

Thank you

It's the same error message as in your previous question. Do you mean we should put up a chart with 2D structures of all ~70 unsupported chemical groups (amide bases, carbon radical, carbocation, carbanions,...) so that the user won't have to look them up on google or wikipedia?

OK, I admit that in this case, there's a slight twist to it. From the ParamChem upload page:What happened is that you uploaded a pdb file with a distorted geometry, which caused OpenBabel to produce a flawed connectivity, which caused the CGenFF program to assume that you're trying to generate parameters for a carbon radical, carbocation or carbanion. You might have known that if you would have followed the above recommendation. Anyway, I updated the ParamChem FAQ to clarify this problem to other users; you are far from alone.

Generally spoken, I'm not the greatest fan of that OpenBabel functionality in ParamChem. People do these kind of things all the time no matter how many warnings we put up, and then they conclude that our program is not working properly. In my personal opinion, it would be better if everyone would be forced to generate a proper connectivity in a format that contains bond order information (such as mol2, mol and sdf). Unfortunately, this strongly goes against the will of the ever growing "bomlev -5 crowd". I sometimes feel like the whole field is going to bomlev -5!

Worked!

OK I'm kind of excited as I didn't get any errors (Getting excited about molecular mechanics on a Friday afternoon? :-/). What the results mean is a whole different issue.

First I forgot to state in my previous post that I DID submit a mol2 file in addition to a PDB but got some error about wild card atom types or something. I simply loaded the PDB in VMD and saved it as mol2. I guess it was crappy from the get go.

Since that and a couple of other dirty trials in VMD/text editing the PDB/Par files didn't work; I just moved back to Gaussian, made a clean azidoproline with the clicked functionality and saved it as a mol2 file. You would think that should give you good results right?

Nope!

It was complaining about unknown bond types. The bond type saved by Gaussian in the mol2 file was "Ar" for some bonds whereas the mol2 documentation calls it "ar". I just changed it to "ar" in the file and voila! Worked!

So just a difference in letter cases? I liked that! easy to fix.

Now by just looking at the results (disregarding the analogy scores of 53 and such!) I'm pretty impressed! Making this file manually, just mining out the analogies and getting all the values and connectivities right was a huge torture. Even if the analogies are pure trash, I love it! So Kenno thank you both for setting up this server (and whoever else is involved) and also for your help and guidance.

I guess analysis of the results is a whole different issue.

( sampler :
RESI Molecule -1.000 ! param penalty= 71.500 ; charge penalty= 53.371
)

That would explain it. VMD's support for bond orders is somewhat flaky:

• When loading a pdb of a small molecule into VMD, it marks all bonds it identifies as single (one could argue that assigning bond orders in not trivial and inherently imprecise, but so is assigning the existence/absence of a bond).
• The VMD mol2 plugin seems to be able to load aromatic and double bonds form a mol2 file, but upon saving, all aromatic bonds (but not the double bonds) are converted to single.
• Modifying bond orders is cumbersome (topo addbond a b -bondorder x), and further frustrated by the fact that bonds with order 1.5 are displayed as double bonds (but saved as single bonds ).

In the end, it's possible to prepare a single resonance structure (ie. a structure without fractional bonds) with VMD, which is acceptable for the CGenFF program because aromatic bonds are internally converted to a single resonance structure anyway. However, doing so is so cumbersome that it's easier to just save a mol2 file with incorrect bond orders and correct them manually with a text editor.

Yes indeed. I have to confess I'm using GaussView too for all my CGenFF work, because no other molecular editor I tried lets you build a chemically correct organic structure so quickly and smoothly. There's only 2 minor issues with GaussView:

• it violates the mol2 standard for aromatic bonds (it clearly says "ar" in there, not "Ar")
• the molecule name gets written into a comment line instead of the designated @<TRIPOS>MOLECULE record

So I open my GaussView-generated mol2 file in vi, change the name under @<TRIPOS>MOLECULE , do :%s/Ar/ar/ and I'm good to go. Still faster than using a different editor. Nevertheless, I'm contemplating making the mol2 reader case-insensitive (or at the least have it accept "Ar" for aromatic) to work around GaussView's violation of the standard.


I'm happy you like it, there's quite a lot of work under the hood there. It will be submitted to a journal soon. Additionally, I like your cautious stance towards the parameters - no amount of sensible analogy can substitute for actual validation. As I mentioned in one of my previous posts, the first thing to do is to have a very close look at the parameters involving the substituted carbon in your 1,2,3-triazole ring. I see two angle parameters with very low penalties and one bond parameter that are visibly terrible, and the dihedrals around that bond are probably trash as well. As I said, I'm working on fixing this particular shortcoming but it's going to be weeks.