Previous Thread
Next Thread
Print Thread
Page 1 of 2 1 2
Dihedral with k=0
#33513 03/03/14 03:11 PM
Joined: Oct 2010
Posts: 63
C
chemist Offline OP
Forum Member
OP Offline
Forum Member
C
Joined: Oct 2010
Posts: 63
Hello, I'm trying to parametrise the DOTA chelator: 1,4,7,10-tetraazacyclododecane-1,4,7,10-tetraacetic acid with the formula (CH2CH2NCH2CO2H)4. This I have attached to a standard amino acid (I've chosen Lysine).

The parameters I need to optimise all involve the DOTA tertiary amine junction to the Lysine branch. Not many cgenff parameters exist for tertiary amines.

This is a big molecule (73 atoms) so I've split the DOTA part (not Lys) into 2 fragments.

One particular dihedral is causing problem to fit. The penalty is below 50, but it must be validated:

NG2S1 CG2O1 CG321 NG301 0.4000 1 0.00 ! Kw_mol , from NG2S1 CG2O1 CG324 NG3P3, penalty= 37.9

OG2D1 CG2O1 CG321 NG301 0.0000 1 0.00 ! Kw_mol , from OG2D1 CG2O1 CG324 NG3P3, penalty= 37.9

CG2O1 CG321 NG301 CG321 2.5000 1 180.00 ! Kw_mol , from NG311 CG321 NG311 CG2R61, penalty= 146.5
CG2O1 CG321 NG301 CG321 1.5000 2 0.00 ! Kw_mol , from NG311 CG321 NG311 CG2R61, penalty= 146.5
CG2O1 CG321 NG301 CG321 0.5000 3 0.00 ! Kw_mol , from NG311 CG321 NG311 CG2R61, penalty= 146.5


The parameter is derived from OG2D1 CG2O1 CG324 NG3P3 which also has n=1, k=0.000. These non-contributing dihedrals are next to a rigid peptide bond, so I'm guessing we do not want them contributing to the potential energy function as they would disrupt planarity?

Do I fit a bond like this? Do I use the Paramchem guess?

Here I tried fitting all 3 dihedrals, and then I will use the Paramchem values for the first two, and the fitting values for the n=1,n=2, n=3 dihedral: www.peecee.dk/upload/view/433133



Also: Kenno, could I ask for the parameters you have for PAP-1 coumarin?

Thanks

Last edited by chemist; 03/03/14 10:18 PM.
Re: Dihedral with k=0
chemist #33529 03/07/14 08:15 PM
Joined: Dec 2005
Posts: 1,535
Forum Member
Offline
Forum Member
Joined: Dec 2005
Posts: 1,535
First a general note: I understand you chose the neutral amine because it's chelating a metal. If so, you'll need to worry about representing the metal chelation correctly, which can be really tricky depending on the metal, its oxidation state and its coordination. Here's a related discussion. if not, you should protonate the nitrogen.

Originally Posted By: chemist
The parameter is derived from OG2D1 CG2O1 CG324 NG3P3 which also has n=1, k=0.000. These non-contributing dihedrals are next to a rigid peptide bond, so I'm guessing we do not want them contributing to the potential energy function as they would disrupt planarity?
It's complicated. The rotation around that bond is controlled by 4 parameters: OCCN, OCCC(beta), NCCN and NCCC(beta). The final profile is the sum of the 1-4 and longer range non-bonded interactions + all the multiplicities of all the dihedrals, each of which is either added or subtracted depending on the multiplicity and geometrical considerations. Also related. We generally limit our choice of multiplicities to the ones that are geometrically and chemically sensible. Also, we want our parameters to be transferable, which leads us to constrain ourselves even further, eg. by setting some types of dihedral parameters always to a predefined value and only fitting the remaining ones, or sometimes by requiring certain fitted parameters to have the same values. Finally, when dihedral contributions with high amplitudes partially or completely cancel out, this leads to a residual angular force which is generally undesirable. All these considerations are taken into account when fitting to the QM. Then, in the next step, bulk phase simulations are done on bigger relevant model systems, and the dihedrals are tweaked to match experimental (often NMR) data, again taking into account all the above considerations. So if you look at the end product, and ask: "why does this parameter have this value", the question is not always easy to answer. You are right that setting the OG2D1 CG2O1 CG324 NG3P3 to a high value and counteracting this with another parameter would disrupt the planarity, and that likely played a role in the decision to set it to 0, but a low value should not be disastrous.

Looking forward, what we would do in the most generic case is scan and fit the N-C-C-N because that's the backbone of your molecule, as opposed to the O, which is a dead end. Prior to that, N-C-C-O is fitted to the C=O wagging motions. However, this generic reasoning is not valid here for 2 reasons. First, we know that the C=O is part of a rigid amide group and that the bond currently under consideration is a free rotor, so it only makes sense to fit the dihedrals on the amide side to the wagging motions, not the dihedral you're talking about. Second, wagging motions of carbonyl groups that are not part of rigid planar rings cannot be fit by angles and dihedrals alone, so an improper is needed, which becomes the primary agent for getting the wagging motion right.

Are you still there? Yeah, force field parameterization is a lot more complicated than meets the eye.

Edit: oh yeah, you probably want a recommendation. I think keeping this one at 0 is not a bad idea. Especially if you're using automatic fitting; it will make things easier on the poor stupid fitting algorithm. (No offense to its creators; I'm calling all fitting algorithms stupid because they don't have common sense or chemical intuition. It's really really hard to put that stuff into a box.)

Also, you attached an energy plot, but such plots don't say much if you don't also show the parameters. If you have a bad plot and the parameters are weird, then it is likely something has gone wrong with the fitting. If you get a good plot and the parameters are weird, that's still not good; they may be overfitted and blow up in your face when transfered. Conversely, if you get a mediocre plot but the parameters are sensible, it might just be that the force field is too approximate to capture physical reality more accurately than that.

Originally Posted By: chemist
Also: Kenno, could I ask for the parameters you have for PAP-1 coumarin?
I finally posted them.

Last edited by Kenno; 03/07/14 10:48 PM.
Re: Dihedral with k=0
chemist #33616 03/25/14 03:07 PM
Joined: Oct 2010
Posts: 63
C
chemist Offline OP
Forum Member
OP Offline
Forum Member
C
Joined: Oct 2010
Posts: 63
Hello Kenno,

Thanks for this very detailed and clear answer.

I still have some questions about the parametrisation of the DOTA chelator in the presence of Gd(III) ion:

Quote:


If you really must, you can put some dummy parameters on the metal ion and use constraints to get the coordination geometry right, but even then, it will be tricky to get the electrostatics right (and you'd have to be careful to avoid destabilizing the MD simulation).




1. I wish to parametrise DOTA with the metal.

Does that mean I:
I) Carry out QM calculations (water interaction, Hessian, and dihedral scans) on the full system.
II) Perform the MM part with a dummy parameter for Gd?


There are CHARMM parameters for the LJ parameters of Gd(III). I can use those for the dummy atom.


2. If I use a dummy for the MM fitting part:

I set the charge of the dummy to a value which will balance out to an total integer charge. I have RESP values of the partial charges, as well as GD-DOTA bonded parameters published by Henriques et al., Int. Journal of Quantum Chemstry, 1999.

3. For the CHARMM parametrisation:
You advocate a splitting approach, but when there's a metal bound, how do I split the chelator? I can't see how this would be done.


4. For the actual MD simulation:

I have seen in the literature that a harmonic restraint approach is used on the dummy bond, angle and dihedral using e.g. the NAMD Colvars feature.

In my case, I can access 2 conformations according to NMR. Do I fit / restrain relative to the minimum energy conformer?






--------


Alternatively:

You write I should optimise the chelator with protonated N. Fine, but at the end of the procedure, I will have a topology file for protonated DOTA. How could I modify this to account for the metal?
I can think of 2 approaches:
I) Create a second topology file without N protonated. The charge imbalance is accounted for by the presence of a dummy metal parameter with charge balancing the 4 protons?


II) Parametrise the whole DOTA chelator partial charges to add up to zero when excluding the protons.
Then define a DOTA topology file
And define a GD topology file containing GD with integer +3 charge. Use bond patches to link DOTA residue with GD residue.
This ignores polarisation of charges.






Are you aware of any resources for metal parametrisation in CHARMM?



Thanks

Re: Dihedral with k=0
chemist #33620 03/25/14 04:41 PM
Joined: Dec 2005
Posts: 1,535
Forum Member
Offline
Forum Member
Joined: Dec 2005
Posts: 1,535
Originally Posted By: chemist
II) Perform the MM part with a dummy parameter for Gd?
No, I talked about the dummy parameters assuming you do not want to parameterize them metal. Which I thought was a safe assumption; you're really getting yourself into a hairy situation if you want to do that.

Originally Posted By: chemist
There are CHARMM parameters for the LJ parameters of Gd(III). I can use those for the dummy atom.
Not in the official release, there aren't. I know there are a lot of unofficial parameter sets for various metals out there, but their quality is very spotty, so I'd need to have a look at what kind of parameterization effort is behind it.


Originally Posted By: chemist
I set the charge of the dummy to a value which will balance out to an total integer charge. I have RESP values of the partial charges, as well as GD-DOTA bonded parameters published by Henriques et al., Int. Journal of Quantum Chemstry, 1999.
If you go for the pragmatic route and don't parameterize the metal (restraining its coordination instead), you can do pretty much anything. If you do want to parameterize the metal, then the charges you got are useless. Among our Frequently Asked Questions, this is a particularly frequent one.

Originally Posted By: chemist
You advocate a splitting approach, but when there's a metal bound, how do I split the chelator? I can't see how this would be done.
I'd need to have some idea of the 3D structure of the complex in order to answer that question (and to give more definite advice on whether there's any point at all in doing the full parameterization, or you should just go for the pragmatic route).

Originally Posted By: chemist
In my case, I can access 2 conformations according to NMR. Do I fit / restrain relative to the minimum energy conformer?
Same as above, but if they are both relevant, and differ to the extent they can't readily interconvert, then the best thing to do would probably be two independent simulations for the two species.

Originally Posted By: chemist
You write I should optimise the chelator with protonated N.
No, I wrote you should do that if not [chelating a metal]. Sorry if I formulated that paragraph ambiguously.

Originally Posted By: chemist
Are you aware of any resources for metal parametrisation in CHARMM?
For free metals in solution, Benoît Roux et al.'s work is authoritative. For metals that are bound in a complex, well, there have been various efforts with various (usually very high) degrees of approximation, but there's very little consistency, and nothing stands out as a "gold standard". Nothing except QM/MM, that is.

Re: Dihedral with k=0
chemist #33621 03/25/14 04:46 PM
Joined: Oct 2010
Posts: 63
C
chemist Offline OP
Forum Member
OP Offline
Forum Member
C
Joined: Oct 2010
Posts: 63
Hello Kenno,

My top priority is the ligand conformation, so let's say I put the metal parameters as lowest priority, I would use a dummy.


Could you explain how exactly you would parametrise that system? If I don't parametrise the metal, do I put a dummy during the parametrisation?

Can I give the dummy a formal charge?

For the quantum calculations, I assign the dummy a type X.
For the MM fitting, what charge sum will I aim for the chelator? The opposite value of the formal charge of the dummy metal?

Thanks

Re: Dihedral with k=0
chemist #33623 03/25/14 04:52 PM
Joined: Oct 2010
Posts: 63
C
chemist Offline OP
Forum Member
OP Offline
Forum Member
C
Joined: Oct 2010
Posts: 63
I try to clarify: What's the chemistry of the system I would parametrise?

since it is DOTA-GD, if I don't parametrise the metal, then would I be parametrising DOTA-DZ, where DZ is my dummy.



http://upload.wikimedia.org/wikipedia/commons/thumb/f/f8/Gadoteric_acid.png/220px-Gadoteric_acid.png

As you see, splitting this is extremely tricky if a metal or a dummy is bound.

Re: Dihedral with k=0
chemist #33624 03/25/14 04:54 PM
Joined: Dec 2005
Posts: 1,535
Forum Member
Offline
Forum Member
Joined: Dec 2005
Posts: 1,535
  • If your ligand is tightly bound to the metal, you can't really consider it independently from the metal, especially in the QM calculations.
  • I wrote "dummy parameters" not "dummy atom". Using a dummy atom for your metal in the QM will have the strangest results.
  • I realize what I wrote so far may seem contradictory; I really feel I cannot give sensible advice without more information.

Re: Dihedral with k=0
chemist #33626 03/25/14 05:02 PM
Joined: Dec 2005
Posts: 1,535
Forum Member
Offline
Forum Member
Joined: Dec 2005
Posts: 1,535
OK, our posts crossed each other. Going by the 2D structure you posted, you'll definitely have to follow the pragmatic route, but some details still need to be filled in. I have to sleep over it for a night, and at any rate, my self-imposed CHARMM forum time quotum is exceeded for today. Remind me if I don't reply in the next 36h or so. In the meanwhile, it would be really helpful if you could give me some 3D information on the 2 conformations.

Re: Dihedral with k=0
chemist #33648 03/26/14 12:07 PM
Joined: Oct 2010
Posts: 63
C
chemist Offline OP
Forum Member
OP Offline
Forum Member
C
Joined: Oct 2010
Posts: 63
Hello Kenno,

Firstly, I must clarify I'm working on a branched version of DOTA, with a lysine branch.

For me, this parametrisation is extremely urgent, so I prefer to work on the fastest method available. We put as lowest priority the accurate representation of the metal, which is poorly treated in Class I ff anyway.


DOTA-GD comes in 2 geometries: Square antiprismatic and Inverse square antiprismatic. These differ in the orientation (clockwise or counterclockwise) of the carboxylate branches.

so from what I understand you advise:

I) Run QM calculations on DOTA-GD, using the relevant basis sets to converge the metal correctly.

# MP2(FULL)/GEN 6D 7F Opt=(Redundant) SCF=Tight pseudo=read Geom=PrintInputOrient



II) Perform MM fit on the chelator. Use dummy parameters on metal for fitting?

III) Run MD with restrained dummy parameters on metal.




I upload the square antiprismatic (Lys-branched) DOTA-GD chelator.

http://peecee.dk/upload/view/434778
http://peecee.dk/upload/view/434779
http://peecee.dk/upload/view/434780

Attached Files
DOTA.png (85.11 KB, 388 downloads)
Re: Dihedral with k=0
chemist #33649 03/26/14 12:10 PM
Joined: Oct 2010
Posts: 63
C
chemist Offline OP
Forum Member
OP Offline
Forum Member
C
Joined: Oct 2010
Posts: 63
Or alternative approach:


I) Parametrising chelator in the absence of metal (not in PDB/PSF)

II) Adding restrained dummy parameters to mimic the metal during the MD

Page 1 of 2 1 2

Moderated by  alex, lennart, rmv 

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.015s Queries: 36 (0.008s) Memory: 0.9978 MB (Peak: 1.1331 MB) Data Comp: Off Server Time: 2020-10-25 02:13:24 UTC
Valid HTML 5 and Valid CSS