Issues while simulating new PSII molecules in GROMACS

User discussions
1

GROMACS version: 2021.4
GROMACS modification: Yes/No
Here post your question

I am a new gromacs user with previous experience in md simulations and I have to simulate a part of the photosystem II complex. My system has a protein core bonded with carotenoids and chlorophylls embedded in a POPC lipid bilayer.
After obtaining pdb structure files from pdb data bank, I obtained the parameters for the new molecules in my system, namely the carotenoid and chlorophylls, from academic resources, in the form of .itp files. Reading around in the manual and forums, I got to know that these files are usually included in topol.top files, as was done in the protein-ligand complex tutorial by Dr. Lemkul. The issue I am facing is that the pdb crystal structure files lack atoms like hydrogen, resulting in a mismatch in the number of atoms from the itp files. Due to this mismatch, I face errors like mismatch in number of coordinates in my gro and top files during solvation and other such steps.
I have tried looking for solutions, here is what I have come across so far and why it isn’t working for me:

  • For adding missing hydrogen atoms to my pdb files, pdb2gmx can be used, but in my case since these molecules are new, I get a fatal error “Residue XXX not found in residue topology database”
  • Making some edits in the GMXLIB directory files, like editing the rtp, residuetypes.dat files and so on. But in this case my new molecules aren’t amino acids, DNA, RNA etc, instead carotenoids and chlorophylls so how do I deal with such molecules?
  • Parametrizing these by myself . I haven’t ever done that and the gromacs manual says that its too much work even for an expert. Also, since they are such big molecules and extremely important for my study, I prefer using previously verified parameters instead of coming up with new ones, to make my research reliable. Having said so, I am willing to learn more about it and get to doing the parametrization but even for that I’d need some tips and guidance.

All in all, Can anyone please tell me a way in which I can simulate my system involving new molecules, preferably using my previously obtained .itp parameter? If any of the above mentioned avenue can help me out, Could anyone please guide me as to how to apply the solution to my system in particular?

I would really be grateful for any sort of help. Thank you so much for reading through. I am tagging a few people who I think might be able to help as my previous posts haven’t been able to reach much people.
@jalemkul @MagnusL @milosz.wieczor @alevilla
Thank you so much for your time and help.

2

There are a few tools in Gromologist that can help. For example, there’s an automated rtp/hdb maker from .itp used in pathway 2 here:

so if you have working .itp files already, and need just to pass things through pdb2gmx, you can use it to make and add .rtp/.hdb files to your force field library.

One important question in that case would be, are your extra molecules bound through covalent or only non-covalent bonds? That will affect the workflow quite a bit.

3

Thank you so much for your response. The extra molecules are bound through non covalent bonds in my case. I took a look at the link you shared so I think procedure 1 would suffice for my case. So what we’re basically doing is to create a pdb file for the complete system including the new molecules which then would be able to be passed through the pdb2gmx tool and can be used in further steps? Also, Sorry if this is too basic to ask, but I have not used gromologist before so how can I go about using the steps available on the link? Do i need to install the package and it will be up and running or do I need to do other things to set it up?
One more thing that I think I should mention is that while the protein topology is written using the gromos54a7 ff, the cofactors have been parametrized using the gromos 53a6 ff. So like could this cause issues while combining the files?
Thanks again for your help, I’m really grateful!

4

Procedure 1 only works if your ligands have hydrogen atoms assigned, because eventually the structure and topology need to have matching/consistent atom numbering. If you want the entire thing to pass through pdb2gmx, go for option 2, because this way you’re also creating the respective residue templates.

Gromologist is installable with pip, so there should be no problem.

I’m not an expert on versions of the Gromos FF, adjacent FF versions are usually similar enough that most ligand parameterizations should be compatible (that’s my guess) but 53a6 and 54a7 are 8 years apart. You’d have to see what specifically was updated between these two versions, and you can also try out Amber/CHARMM FFs that tend to be more broadly tested/validated.

5

Thank you, I’ll try using your suggestions and see how it works out :)