GROMACS version: 2019.3
GROMACS modification: Yes/No
Here post your question
Dear all,
I am trying to do MD simultion on protein in complex with Fe2+ and/or Mn2+. Does anyone know where I can find the force field parameters for these ions and how should I include them in the simulation?
With many thanks!
Best regards
Mn2+ parameters for AMBER: http://research.bmh.manchester.ac.uk/bryce/amber
Mn2+ parameters for CHARMM: dx.doi.org/10.1021/jp309150r
Fe2+ should already be present in both, though is generally designed to only be used in heme (at least for CHARMM), so your mileage may vary for other species.The paper above for CHARMM has new Fe2+ parameters but I have never used them so I cannot vouch for their quality.
Add the relevant parameters to ffnonbonded.itp and if you need pdb2gmx to create a topology with them, also add relevant information to the .rtp and .atp files.
Hello Jalemkul,
Thank you for your reply. I am trying to build the protein-Fe complex with charmm36-mar2019.ff. In the charmm36-mar2019.ff directory, there are the atomtypes.atp, ffnonbonded.itp and merged.rtp files. These files all contains lines for FE. For example, in the ffnonbonded.itp file, there is a line reading:
FE 26 55.847000 0.000 A 0.115816833358 -0.00000
In the atomtypes.atp file, there is a line reading:
FE 55.84700 ; heme iron 56
In the merged.rtp file, there are lines reading:
[ HEME ]
[ atoms ]
FE FE 0.240 0
However, when I use pdb2gmx to build the topology, I get the following error:
Fatal error:
Residue ‘FE’ not found in residue topology database
So Could you further explain how to deal with this situation?
Best regards
There is no standalone residue entry for iron in the force field. There are parameters for the heme Fe2+, and you see it appear in the HEME residue, but there is no FE residue.
What is it that you’re actually dealing with? Free ions? Part of some other cofactor? We need more detail to be able to help.
Hello Jalemkul,
I am just dealing with a Fe-Protein complex, and there is no other co-factor in the complex.
The coordinate for this complex can be found at:
https://drive.google.com/file/d/1jNx8Sfl7v8PY-dGF4mRx5OCKbkr3AduI/view?usp=sharing
Best regards
Follow http://manual.gromacs.org/current/how-to/topology.html#adding-a-new-residue
Be advised that the CHARMM Fe2+ parameters are only developed for heme and have not been validated for His and Asp/Glu interactions like you have in your complex. You should do some validation against QM interaction energies and dimer geometries before trying to just apply the existing parameters. I suspect you will find you need NBFIX terms.
I added the following lines to the merged.rtp file:
[ FE ]
[ atoms ]
FE FE 2.000 0
Now the pdb2gmx has worked through. However, could you also give me some instruction on:
(1) how to validate for Fe2+ interactions with His and Asp/Glu?
(2) how to validate against QM interaction energies?
Best regards
I would recommend you thoroughly read the CHARMM force field papers, many of which describe these aspects in great detail.
Do you mean the paper for the metal ion charmm ff you listed above all Just the Charmm ff paper?
And I wonder what is the most reliable method to calculate the binding free energy between the Fe2+ and the protein, the SMD method or the alchemical method?
Best regards.
Any review of the CHARMM force field will tell you all about how the force field is parametrized, including the proper QM methods.
Decoupling charged species is tricky. My first try would be with SMD, but given the partial covalent nature of the actual interactions of these species with proteins, obtaining something physically realistic may be extremely challenging.
Hello everyone,
I am trying to do MD for the above mentioned protein-Fe2+ complex. I first did minimization, which finished briefly with the following output:
Steepest Descents converged to machine precision in 40 steps,
but did not reach the requested Fmax < 10.
Potential Energy = -1.1980662e+06
Maximum force = 3.3884829e+08 on atom 2806
Norm of force = 1.8291110e+06
However, when I did the NVT step, I got the following error:
step 5: One or more water molecules can not be settled.
Check for bad contacts and/or reduce the timestep if appropriate.
Back Off! I just backed up step5b.pdb to ./#step5b.pdb.1#
Back Off! I just backed up step5c.pdb to ./#step5c.pdb.1#
Wrote pdb files with previous and current coordinates
Segmentation fault (core dumped)
Is there a way to tackle this problem?
Best regards
There’s your problem. The system is not physically sound.
Yeah, the normal magnitude of the system maximum force is 10^2, isn’t it?
The atom 2806 is the OE1 atom of Glu189, one of the contact atom of the Fe2+. So is the problem arose by the fact that the OE1 atom is too close to the FE2+ in the original coordinates? How to deal with this problem?
No, it means the interaction is physically unstable. Have you done the QM interaction energy between acetate and Fe2+ and compared it to the force field value?
Hello,
I have no idea how to do QM interaction energy between acetate and Fe2+. Could you suggest any tutorial for this?
This won’t deal with metals but you can learn more about the logic here: http://ffparam.umaryland.edu/
If you want to study transition metals, you will have to do some force field parametrization of your own. Existing parameters are generally far from robust.
HI Jason,
been facing the same issue and have reached a repeated dead end of segmentation fault. would be helpful if you shared if you were able to solve your problem of running MD with an Fe2+ atom
I have the same problem, my protein contains FE-S cluster, and I don’t know how to add force field parameters
I have found this paper about QM interaction:
https://onlinelibrary.wiley.com/doi/abs/10.1002/jcc.23287