SAM co-factor with CHARMM force field for GROMACS silmualtion

Hello All,

Following the GROMACS tutotial on Protein-ligand complex, I am trying to simulate a protein having SAM as one of its ligands using the CHARMM force field. All my attempts to get the topology file for SAM failed. I would like to know if anyone has succeeded in writing the topology (itp file for SAM with CHARMM FF) of this well-known ligand to help me do so as well.

Best,
Bourgeois

For CHARMM compatible topologies, I think your choices are CHARMM-GUI, https://cgenff.com/ or Brooks Lab - Match - University of Michigan (possibly via STaGE, GROMACS / STaGE · GitLab). Have you tried all these alternatives?

I have tried CHARMM-GUI and CGenFF but could not figure it out.

Hi Gabbi,

I would like to obtain CHARMM ff for SAM molecule but couldn’t achieve it via CHARMM-GUI, and CGENFF. Then I tried to form it manually by collecting info from toppar_prot and toppar_na, however I still don’t have the full bonded terms list. The reference (Senn et al., 2005, 127, 13643-13655) list newly generated terms but some are missing or I couldn’t find. Have you obtained the parameter file?

Hello Gulshah,

Getting the CHARMM force field for SAM can look tricky at first, but it’s actually more straightforward than it seems. SAM is already defined in charmm36-jul2022.ff.

If you’re simulating SAM together with a protein, I recommend the following steps:

  1. Prepare your PDB file

    • Clean your protein PDB so that it contains the protein plus the SAM ligand.

    • If your original PDB has other ligands besides SAM, extract them into separate PDB files and treat them individually. Keep only the protein and SAM together in the same file.

  2. Generate topology with GROMACS

    • Run:

      gmx pdb2gmx -f protein_SAM.pdb -o processed.gro -p topol.top -ff charmm36-jul2022
      
      
    • Alternatively, if you just want SAM by itself, you can run:

      gmx pdb2gmx -f SAM.pdb -o SAM.gro -p topol.top -ff charmm36-jul2022
      
      
  3. Check the output

    • For the combined system, you’ll usually get two .itp files:

      • topol_Protein_chain_A.itp (protein)

      • topol_Other_chain_A2.itp (SAM)

    • Both are automatically included in the topol.top file.

This way, you don’t need to manually collect bonded terms–everything is handled within CHARMM36.

Best regards,
Gabbi

Thank you very much Gabbi, it worked!

Dear Gabbi,

Hope you are doing well. When I run the following command to generate the topology “gmx pdb2gmx -f protein_SAM.pdb -o processed.gro -p topol.top -ff charmm36-jul2022"I encounter the following fatal error

“Atom S in residue SAM 246 was not found in rtp entry SAM with 50 atoms while sorting atoms.”

It seems that the sulfur atom present in the SAM residue of my PDB file is not being recognized by the corresponding SAM entry in the CHARMM36-jul2022 residue topology (rtp) file. Could you please advise on how best to resolve this issue?

Any guidance would be greatly appreciated

Regards,

Archana

The atom names in your PDB file are simply inconsistent with the residue definition in the .rtp file. The S atom should be named SD. If needed, change any incorrect atom names in the PDB file to match the residue definition.

Thank you so much for your response. When I compared the .rtp file with the PDB file, I noticed that the atom names are different in the two files. For example, in the PDB file, the atoms are simply named as O, N, S, and H, whereas in the .rtp file they are explicitly defined as HN1, HN2, CA, SD, etc. Because of this, it is difficult to identify the corresponding atoms and rename them correctly.

I downloaded the sam.sdf file from PubChem and attempted molecular dynamics simulations using PDB files that were either directly converted from the SDF file or generated after docking. Could you kindly let me know which database you used to download the PDB file that has atom names consistent with the .rtp file?

Regards

Archana

There is no guarantee that a PDB file will agree with the force field definitions, outside of “core” residues like canonical amino acids and nucleobases. SDF files are even worse :) You need to correct the atom names in the PDB file so they agree with the residue definition, simple as that. It’s a common task.

Hello,

Thank you for your help. I performed MD simulations of a methyltransferase with SAM retained in the protein structure and a purine-like small-molecule ligand, using the SAM parameters and coordinates provided in the CHARMM 2022 package. The simulation runs stably, but the ligand diffuses out of the binding pocket after ~25 ns.

I would appreciate guidance on:

  1. Recommended ways to verify that SAM coordinates and parameters (from aminoacid.rtp) are correctly interpreted during system setup and dynamics.

  2. Whether ligand diffusion in this context typically reflects force-field/parameter issues, or is more commonly due to factors such as protonation states, missing waters, or incomplete stabilizing interactions.

  3. Best practices for equilibration or validation (e.g., temporary restraints) to ensure proper system relaxation without biasing ligand binding.

Thank you for your help.